In our recent paper, we presented a new approach for generating a guaranteed affine underestimator for a black-box convex function on a box domain, by tractable derivative-free sampling. Felix Bonhoff's master's thesis (RWTH Aachen, 2023) presented a variant of this approach that requires fewer samples.

This repository provides a Julia implementation of our new sampling approach, tested in Julia v1.9. If you make use of this implementation in your own work, please cite the accompanying article:

Yingkai Song, Huiyi Cao, C. Mehta, and Kamil A. Khan, Bounding convex relaxations of process models from below by tractable black-box sampling, Computers & Chemical Engineering, 153:107413 (2021). doi:10.1016/j.compchemeng.2021.107413

In the Julia REPL, enter the following command:

import Pkg; Pkg.add(url="https://github.com/kamilkhanlab/ConvexSampling.jl")

Now the package's module can be invoked with using ConvexSampling.

Consider the following convex quadratic function f:

using LinearAlgebra

A = [65.0 56.0; 56.0 65.0]
b = [6.0, 2.0]
c = 23.0

f(x) = dot(x, A, x) + dot(b, x) + c

xL = [-5.0, -3.0]
xU = [5.0, 3.0]

on the box domain: all(xL .<= x .<= xU). Suppose we wish to construct affine underestimators and/or lower bounds of f on its box domain.

1. Once the package is installed, we can load its module:

   using ConvexSampling

2. Now, f can be sampled at either 5 domain points using the approach of Song et al.:

   data = sample_convex_function(f, xL, xU; stencil=:compass)

   or at 4 domain points using the approach of Bonhoff:

   data = sample_convex_function(f, xL, xU; stencil=:simplex)

   In general, for a function of n variables, stencil=:compass will sample this function 2n+1 times, while stencil=:simplex will sample it n+2 times, but may ultimately yield a looser relaxation.

   Our approach handles univariate and multivariate functions differently, so our implementation has Julia dispatch over typeof(xL):

   • If xL is a Float64, then f must have the signature f(x::Float64) -> Float64, and our univariate formulas will be applied. Any univariate function must be entered this way.
   • If xL is a Vector{Float64} with $\geq 2$ components, then f must have the signature f(x::Vector{Float64}) -> Float64, and our multivariate formulas will be applied.

   Either way, xU must have the same type and number of components as xL.

3. Using the sampled information, we can construct a guaranteed affine underestimator of f on its box domain:

   fAffine = construct_underestimator(data)

   The constructed function fAffine underestimates f on its box domain, so fAffine(x) <= f(x) whenever all(xL .<= x .<= xU). We can instead obtain this underestimator as its constant coefficients:

   w0, b, c = evaluate_underestimator_coeffs(data)

   in which case fAffine(x) == c + dot(b, x - w0) for all x.

   We can also evaluate a guaranteed lower bound of f on its box domain:

   fL = evaluate_lower_bound(data)

   Then, we will have f(x) >= fL for each x in the domain.

4. The function f may be plotted with its sampling-based underestimator fAffine and lower bound fL:

   graph = plot_underestimator(data, f)
   @show graph

   

   Note that if the plot_sampling_underestimator function is entered directly in the REPL, the @show command is not required.

Suppose we have a convex function $f$ of $n$ variables, defined on a box domain $X = [\mathbf{x}^L, \mathbf{x}^U]$. Our new underestimating approach samples $f$ at $(2n+1)$ domain points: the midpoint of $X,$ and perturbations of this midpoint in each positive/negative coordinate direction. These sampled values are then tractably assembled using new finite difference formulas to yield guaranteed affine underestimators and guaranteed lower bounds for $f$ on $X.$ These underestimators are guaranteed by convex analysis theory; roughly, the sampled information is sufficient to infer a compact polyhedral set that encloses all subgradients of $f$ at the midpoint of $X.$ Using this information, we can deduce the "worst-case" convex functions that are consistent with the sampled data.

As in our paper, this implementation also allows for absolute error in evaluating $f$, and for off-center sampling stencils. When $n=1$, we additionally exploit the fact that each domain point is collinear with all three sampled points.

A newer procedure is implemented in Bonhoff's master's thesis (RWTH Aachen, 2023), requiring $(n+2)$ samples instead of $(2n+1)$ samples.

• Yingkai Song, Huiyi Cao, C. Mehta, and Kamil A. Khan, Bounding convex relaxations of process models from below by tractable black-box sampling, Computers & Chemical Engineering, 153:107413, 2021, DOI: 10.1016/j.compchemeng.2021.107413
• Felix Bonhoff, master's thesis, RWTH Aachen, 2023