Extensions
Here we provide guidance to various ways InfiniteOpt
can be extended.
Overview
Extendibility is one of the core ideas of InfiniteOpt
so that it can serve as a convenient tool for those developing and implementing advanced techniques for infinite dimensional optimization problems. Thus, InfiniteOpt
is developed in a modular manner to readily accommodate user-defined functionality and/or to serve as useful base in writing a JuMP
extension. Admittedly, this modularity is not perfect and comments/suggestions are welcomed to help us improve this.
Infinite Domains
Infinite domains are used to characterize the behavior of infinite parameters and used to govern the behavior of supports in InfiniteOpt
. Here we walk through how user-defined domains can be added to various degrees of functionality. A template is provided in ./test/extensions/infinite_domain.jl
. The extension steps employed are:
- Define the new
struct
infinite domain type (only thing required as bare minimum) - Extend
InfiniteOpt.round_domain
(enables safe use of significant digit rounding) - Extend
InfiniteOpt.supports_in_domain
(enables error checking of supports) - Extend
InfiniteOpt.generate_support_values
(enables support generation vianum_supports
keyword arguments) - If a lower bound and upper bound can be reported, extend
JuMP
lower bound and upper bound methods (enables automatic bound detection inintegral
) - Extend
InfiniteOpt.MeasureToolbox.generate_expect_data
(enables the use ofexpect
)
As an example, let's create a univariate disjoint interval domain as an infinite domain type. This corresponds to the domain $[lb_1, ub_1] \cup [lb_2, ub_2]$ where $ub_1 \leq lb_2$. First, we need to create the DataType
with inheritance from InfiniteScalarDomain
:
using InfiniteOpt
struct DisjointDomain <: InfiniteOpt.InfiniteScalarDomain
lb1::Float64
ub1::Float64
lb2::Float64
ub2::Float64
# constructor
function DisjointDomain(lb1::Number, ub1::Number, lb2::Number, ub2::Number)
if lb1 > ub1 || lb2 > ub2 || ub1 > lb2
error("Invalid bounds")
end
return new(convert(Float64, lb1), convert(Float64, ub1),
convert(Float64, lb2), convert(Float64, ub2))
end
end
Notice that we also define the constructor function to error check and convert as needed (this is recommended, but not required). For basic functionality this is all we have to do to add a domain in InfiniteOpt
.
We can now define infinite parameters using this domain via @infinite_parameter
both anonymously and explicitly:
julia> model = InfiniteModel();
julia> t = @infinite_parameter(model, domain = DisjointDomain(0, 1, 3, 4), base_name = "t")
t
julia> @infinite_parameter(model, t in DisjointDomain(0, 1, 3, 4))
t
Once defined (without further extension), these parameters can be used as normal with the following limitations:
- Supports must be specified manually (
num_supports
is not enabled) - Supports will not be checked if they are in the domain of the infinite domain
- Domain bounds cannot be queried.
- The
DiscreteMeasureData
orFunctionalDiscreteMeasureData
must be provided explicitly to evaluate measures
However, all of these limitations except for the last one can be eliminated by extending a few functions as outlined above. To address the last one, we need to extend generate_integral_data
. See [Measure Evaluation Techniques
] for details.
To enable support domain checking which is useful to avoid strange bugs, we will extend InfiniteOpt.round_domain
which rounds the domain to use proper significant digits and InfiniteOpt.supports_in_domain
which returns a Bool
whether a vector of supports is in the domain:
function InfiniteOpt.round_domain(
domain::DisjointDomain,
sig_digits::Int
)
lb1 = round(domain.lb1, sigdigits = sig_digits)
ub1 = round(domain.ub1, sigdigits = sig_digits)
lb2 = round(domain.lb2, sigdigits = sig_digits)
ub2 = round(domain.ub2, sigdigits = sig_digits)
return DisjointDomain(lb1, ub1, lb2, ub2)
end
function InfiniteOpt.supports_in_domain(
supports::Union{Number, Vector{<:Number}},
domain::DisjointDomain
)
return all((domain.lb1 .<= supports .<= domain.ub1) .| (domain.lb2 .<= supports .<= domain.ub2))
end
Now the checks are enabled, so the following would yield an error because the support is not in the domain:
julia> @infinite_parameter(model, domain = DisjointDomain(0, 1, 3, 4), supports = 2)
ERROR: At none:1: `@infinite_parameter(model, domain = DisjointDomain(0, 1, 3, 4), supports = 2)`: Supports violate the domain bounds.
To enable automatic support generation via the num_supports
keyword and with functions such as fill_in_supports!
, we will extend InfiniteOpt.generate_support_values
:
struct DisjointGrid <: InfiniteOpt.PublicLabel end
function InfiniteOpt.generate_support_values(
domain::DisjointDomain;
num_supports::Int = InfiniteOpt.DefaultNumSupports,
sig_digits::Int = InfiniteOpt.DefaultSigDigits
)
length_ratio = (domain.ub1 - domain.lb1) / (domain.ub1 - domain.lb1 + domain.ub2 - domain.lb2)
num_supports1 = Int64(ceil(length_ratio * num_supports))
num_supports2 = num_supports - num_supports1
supports1 = collect(range(domain.lb1, stop = domain.ub1, length = num_supports1))
supports2 = collect(range(domain.lb2, stop = domain.ub2, length = num_supports2))
return round.([supports1; supports2], sigdigits = sig_digits), DisjointGrid
end
Now automatic support generation is enabled, for example:
julia> par = @infinite_parameter(model, domain = DisjointDomain(0, 2, 3, 4), num_supports = 10)
noname
julia> supports(par)
10-element Vector{Float64}:
0.0
0.333333333333
0.666666666667
1.0
1.33333333333
1.66666666667
2.0
3.0
3.5
4.0
We can extend the appropriate JuMP
upper and lower bound functions if desired which are:
JuMP.has_lower_bound
JuMP.lower_bound
JuMP.set_lower_bound
JuMP.has_upper_bound
JuMP.upper_bound
JuMP.set_upper_bound
However, if we want has_lower_bound = false
and has_upper_bound = false
then no extension is needed. For our current example we won't do this since lower and upper bounds aren't exactly clear for a disjoint interval. Please refer to the template in ./InfiniteOpt/test/extensions/infinite_domain.jl
to see how this is done.
Finally, we can optionally enable the use of expect
taken with respect to infinite parameters with this new domain type by extending InfiniteOpt.MeasureToolbox.generate_expect_data
:
function InfiniteOpt.MeasureToolbox.generate_expect_data(domain::DisjointDomain,
pref::GeneralVariableRef,
num_supports::Int;
kwargs...
)
for (k, _) in kwargs
error("Keyword argument `$k` not supported for expectations over ",
"disjoint domains.")
end
coeff_func = (supps) -> ones(size(supps)[end]) ./ size(supps)[end]
return InfiniteOpt.FunctionalDiscreteMeasureData(pref, coeff_func, 0, All)
end
The above implementation simply sums over all the supports associated with pref
and divides by the total number. Now we can use expect
:
julia> @variable(model, y, Infinite(t))
y(t)
julia> expect(y, t)
𝔼{t}[y(t)]
Derivative Evaluation Methods
Derivative evaluation methods are used to dictate how we form the auxiliary derivative evaluation equations (derivative constraints) when we evaluate derivatives in InfiniteOpt. Users may wish to implement their own methods beyond the finite difference and orthogonal collocation ones we natively provide. Thus, we provide an API to do just this. A complete template is provided in ./test/extensions/derivative_method.jl
to help streamline this process. The extension steps are:
- Define the new method
struct
that inherits from the correctAbstractDerivativeMethod
subtype - Extend
InfiniteOpt.allows_high_order_derivatives
- Extend
InfiniteOpt.generative_support_info
if the method is aGenerativeDerivativeMethod
- Extend
InfiniteOpt.derivative_expr_data
- Extend
InfiniteOpt.make_indexed_derivative_expr
.
To exemplify this process let's implement 1st order explicit Euler which is already implemented via FiniteDifference(Forward())
, but let's make our own anyway for the sake of example. For a first order derivative $\frac{d y(t)}{dt}$ explicit Euler is expressed:
\[y(t_{n+1}) = y(t_n) + (t_{n+1} - t_{n})\frac{d y(t_n)}{dt}, \ \forall n = 0, 1, \dots, k-1\]
Let's get started with step 1 and define our new method struct:
using InfiniteOpt
struct ExplicitEuler <: NonGenerativeDerivativeMethod end
Notice that our method ExplicitEuler
inherits from NonGenerativeDerivativeMethod
since explicit Euler uses the existing support scheme without adding any additional supports. If our desired method needed to add additional supports (e.g., orthogonal collocation over finite elements) then we would need to have used GenerativeDerivativeMethod
.
Now we need to decide if this method will directly support higher order derivatives. In this case, let's say it won't and define:
InfiniteOpt.allows_high_order_derivatives(::ExplicitEuler) = false
Conversely, we could set the output to true
if we wanted to directly support higher order derivatives. In which case, we would need to take the order into account in steps 4 and 5.
Since, this is a NonGenerativeDerivativeMethod
we skip step 3. This is however exemplified in the extension template.
For step 4, we extend InfiniteOpt.derivative_expr_data
. This function generates all the needed data to make the expressions necessary to build the derivative evaluation equations (derivative constraints). We assume these relations to be of the form $h = 0$ where $h$ is a vector of expressions. Thus, mathematically $h$ should be of the form:
\[\begin{aligned} &&& y(t_{2}) - y(t_{1}) - (t_{2} - t_{1})\frac{d y(t_1)}{dt} \\ &&& \vdots \\ &&& y(t_{n+1}) - y(t_n) - (t_{n+1} - t_{n})\frac{d y(t_n)}{dt} \\ \end{aligned}\]
The required data must include the support indices used for each derivative variable and then any other constants needed. In this case, we will need the indices $\{1, \dots, n\}$ and no additional data (additional data is exemplified in the extension template). With this in mind let's now extend InfiniteOpt.derivative_expr_data
:
function InfiniteOpt.derivative_expr_data(
dref::GeneralVariableRef,
order::Int,
supps::Vector{Float64},
method::ExplicitEuler
)
# generate the support indices to be used for each call of `make_indexed_derivative_expr`
idxs = 1:length(supps)-1
# return the indexes and the other iterators
return (idxs, ) # output must be a tuple
end
Finally, we just need to extend InfiniteOpt.make_indexed_derivative_expr
. This will be used to create derivative expressions for each index determined (and additional datum) produced by derivative_expr_data
.
function InfiniteOpt.make_indexed_derivative_expr(
dref::GeneralVariableRef,
vref::GeneralVariableRef,
pref::GeneralVariableRef,
order::Int,
idx,
supps::Vector{Float64}, # ordered
write_model::Union{InfiniteModel, AbstractTransformationBackend},
::ExplicitEuler,
# put extra data args here (none in this case)
)
# generate the derivative expression h corresponding to the equation of
# the form h = 0
d = InfiniteOpt.make_reduced_expr(dref, pref, supps[idx], write_model)
v1 = InfiniteOpt.make_reduced_expr(vref, pref, supps[idx], write_model)
v2 = InfiniteOpt.make_reduced_expr(vref, pref, supps[idx + 1], write_model)
return JuMP.@expression(write_model, -(supps[idx+1] - supps[idx]) * d + v2 - v1)
end
We used InfiniteOpt.make_reduced_expr
as a convenient helper function to generate the semi-infinite variables/expressions we need to generate at each support point.
If your new derivative method is not compatible can not be broken up into the derivative_expr_data
-make_indexed_derivative_expr
workflow, then you can instead extend InfiniteOpt.evaluate_derivative
. This is discouraged where possible since it may make your method incompatible with backends that depend on the preferred workflow.
Now that we have completed all the necessary steps, let's try it out!
julia> model = InfiniteModel();
julia> @infinite_parameter(model, t in [0, 10], num_supports = 3,
derivative_method = ExplicitEuler());
julia> @variable(model, y, Infinite(t));
julia> dy = deriv(y, t);
julia> evaluate(dy)
julia> derivative_constraints(dy)
2-element Vector{InfOptConstraintRef}:
-5 d/dt[y(t)](0) + y(5) - y(0) = 0
-5 d/dt[y(t)](5) + y(10) - y(5) = 0
We implemented explicit Euler and it works! Now go and extend away!
Measure Evaluation Techniques
Measure evaluation methods are used to dictate how to evaluate measures. Users may wish to apply evaluation methods other than Monte Carlo sampling and/or Gaussian quadrature methods. To create multiple measures using the same new evaluation methods, users may want to embed the new evaluation method under the integral
function that does not require explicit construction of AbstractMeasureData
.
Creating a DiscreteMeasureData Object
The basic way to do that is to write a function that creates DiscreteMeasureData
object and pass the object to measure
. This considers a measure approximation of the form:
\[\sum_{i \in I} \alpha_i f(\tau_i) w(\tau_i)\]
where $\alpha_i$ are coefficients, $f(\cdot)$ is the expression being measured, $w(\cdot)$ is a weighting function, and $i \in I$ indexes the support points. Let's consider defining a function that enables the definition of a uniform grid for a univariate infinite parameter in IntervalDomain
. This example approximation uses a uniformly spaced supports $\tau_i$ with $\alpha_i = \frac{ub - lb}{|I|}$:
function uniform_grid(param, num_supports)
lb = lower_bound(param)
ub = upper_bound(param)
supps = collect(LinRange(lb, ub, num_supports))
coeffs = ones(num_supports) / num_supports * (ub - lb)
return DiscreteMeasureData(param, coeffs, supps, lower_bound = lb, upper_bound = ub)
end
It is necessary to pass the infinite parameter reference since the construction of measure data object needs parameter information. Now let's apply the new uniform_grid
function to infinite parameters in intervals. We consider a time parameter t
and 2D spatial parameter x
, and two variables f(t)
and g(x)
parameterized by t
and x
, respectively:
julia> m = InfiniteModel();
julia> @infinite_parameter(m, t in [0, 5]);
julia> @variable(m, y, Infinite(t));
Now we can use uniform_grid
to construct a DiscreteMeasureData
and create a measure using the measure data, as shown below:
julia> tdata = uniform_grid(t, 6);
julia> y_meas = measure(y, tdata)
measure{t ∈ [0, 5]}[y(t)]
julia> expand(y_meas)
0.8333333333333333 y(0) + 0.8333333333333333 y(1) + 0.8333333333333333 y(2) + 0.8333333333333333 y(3) + 0.8333333333333333 y(4) + 0.8333333333333333 y(5)
Integral Evaluation Methods
For integrals, we can implement a new approximation method via the extension of InfiniteOpt.MeasureToolbox.generate_integral_data
. This will allow users to use their custom measure evaluation methods in the integral
function that does not explicitly require a measure data object. A template for how such an extension is accomplished is provided in ./test/extensions/measure_eval.jl
. In general, such an extension can be created as follows:
- Define a new empty
struct
(e.g.my_new_fn
) that dispatches your function - Extend
InfiniteOpt.MeasureToolbox.generate_integral_data
, wheremethod
is of the typemy_new_fn
, anddomain
needs to be a subtype ofAbstractInfiniteDomain
that you wish to apply the new evaluation method to.
Note that this procedure can be used to generate new measure evaluation methods not only for existing infinite domains, but also for user-defined infinite domains.
For example, an extension of InfiniteOpt.MeasureToolbox.generate_integral_data
that implements uniform grid for univariate and multivariate parameters in IntervalDomain
can be created as follows:
struct UnifGrid <: InfiniteOpt.MeasureToolbox.AbstractUnivariateMethod end
function InfiniteOpt.MeasureToolbox.generate_integral_data(
pref::InfiniteOpt.GeneralVariableRef,
lower_bound::Real,
upper_bound::Real,
method::UnifGrid;
num_supports::Int = InfiniteOpt.DefaultNumSupports,
weight_func::Function = InfiniteOpt.default_weight
)
increment = (upper_bound - lower_bound) / (num_supports - 1)
supports = [lower_bound + (i - 1) * increment for i in 1:num_supports]
coeffs = ones(num_supports) / num_supports * (upper_bound - lower_bound)
return InfiniteOpt.DiscreteMeasureData(
pref, coeffs, supports,
weight_function = weight_func,
lower_bound = lower_bound,
upper_bound = upper_bound)
end
Also notice that users are free to pass keyword arguments for their new functions in addition to the required positional arguments. This might be needed in case if the new evaluation method requires additional information not captured in the default positional arguments. For example, the multivariate parameter version above needs to know if the multivariate parameter is independent in order to throw a warning when needed.
We create measure for y
using the uniform_grid
method:
julia> y_int = integral(y, t, num_supports = 6, eval_method = UnifGrid())
∫{t ∈ [0, 5]}[y(t)]
julia> expand(y_int)
0.8333333333333333 y(0) + 0.8333333333333333 y(1) + 0.8333333333333333 y(2) + 0.8333333333333333 y(3) + 0.8333333333333333 y(4) + 0.8333333333333333 y(5)
Here we go! We can freely use UnifGrid
for infinite parameters residing in IntervalDomain
s now.
Measure Data
Measures are used to evaluate over infinite domains. Users may wish to employ measure abstractions that cannot be readily represented with coefficients and discretized supports, and thus may wish to extend InfiniteOpt
's measure framework to accommodate other paradigms. This can be accomplished by implementing a user-defined measure data structure that inherits from AbstractMeasureData
. A template for how such an extension is accomplished is provided in ./test/extensions/measure_data.jl
. The extension steps employed are:
- Define the new data struct inheriting from
AbstractMeasureData
(required) - Extend
InfiniteOpt.parameter_refs
(required) - Extend
InfiniteOpt.expand_measure
(required) - Extend
InfiniteOpt.supports
(required if parameter supports are employed in any way) - Extend
InfiniteOpt.add_supports_to_parameters
(required if parameter supports are employed in measure evaluation) - Extend
InfiniteOpt.coefficients
(useful getter method if applicable) - Extend
InfiniteOpt.weight_function
(useful getter method if applicable) - Extend
InfiniteOpt.support_label
(needed to enable deletion if supports are added.) - Extend
InfiniteOpt.generative_support_info
(Needed if the measure will cause the creation of generative supports) - Make simple measure constructor wrapper of
measure
to ease definition.
To illustrate how this process can be done, let's consider extending InfiniteOpt
to include measure support for assessing the variance of random expressions. The variance of an expression $f(x, \xi)$ where $x \in \mathbb{R}^n$ are finite variables and $\xi \in \mathbb{R}^m$ are random infinite parameters is defined:
\[\mathbb{V}[f(x, \xi)] = \mathbb{E}\left[(f(x, \xi) - \mathbb{E}[f(x, \xi)])^2 \right].\]
Note, we could just accomplish this by nested use of expect
, but we implement this example to illustrate the mechanics of extension.
First, let's define our new struct
inheriting from AbstractMeasureData
:
using InfiniteOpt, Distributions
struct DiscreteVarianceData <: AbstractMeasureData
parameter_refs::Union{GeneralVariableRef, Vector{GeneralVariableRef}}
supports::Vector
label::DataType
# constructor
function DiscreteVarianceData(
parameter_refs::Union{GeneralVariableRef, AbstractArray{<:GeneralVariableRef}},
supports::Vector,
label::DataType = InfiniteOpt.generate_unique_label()
)
# convert input as necessary to proper array format
if parameter_refs isa AbstractArray
parameter_refs = convert(Vector, parameter_refs)
supports = [convert(Vector, arr) for arr in supports]
end
return new(parameter_refs, supports, label)
end
end
We have defined our data type, so let's extend the measure data query methods to enable its definition. These include:
function InfiniteOpt.parameter_refs(data::DiscreteVarianceData)
return data.parameter_refs
end
function InfiniteOpt.supports(data::DiscreteVarianceData)
return data.supports
end
function InfiniteOpt.support_label(data::DiscreteVarianceData)
return data.label
end
We also need to extend InfiniteOpt.add_supports_to_parameters
since support points will be used for measure evaluation later:
function InfiniteOpt.add_supports_to_parameters(data::DiscreteVarianceData)
pref = parameter_refs(data)
supps = supports(data)
label = support_label(data)
add_supports(pref, supps, label = label)
return
end
Note that extending supports
is not needed for abstractions that don't involve discretization of the infinite parameter(s), such as the case for certain outer approximation techniques. Our extension is now sufficiently constructed to allow us to define out the new variance measure via measure
. For example:
# Setup the infinite model
model = InfiniteModel()
@infinite_parameter(model, xi ~ Normal(), num_supports = 2) # few for simplicity
@variable(model, y, Infinite(xi))
@variable(model, z)
# Define out new variance measure
data = DiscreteVarianceData(xi, supports(xi))
mref = measure(2y + z, data, name = "Var")
# output
Var{xi}[2 y(xi) + z]
Thus, we can define measure references that employ this our new data type.
We can define variance measures now, but now let's extend expand_measure
so that they can be expanded into finite expressions:
function InfiniteOpt.expand_measure(
expr::JuMP.AbstractJuMPScalar,
data::DiscreteVarianceData,
write_model::Union{InfiniteModel, AbstractTransformationBackend}
)
# define the expectation data
expect_data = DiscreteMeasureData(
data.parameter_refs,
1 / length(data.supports) * ones(length(data.supports)),
data.supports, is_expect = true, label = data.label)
# define the mean
mean = measure(expr, expect_data)
# return the expansion of the variance using the data mean
return expand_measure((copy(expr) - mean)^2, expect_data, write_model)
end
Notice that we reformulated our abstraction in terms of measures with DiscreteMeasureData
so that we could leverage the existing expand_measure
library. Now, new the measure type can be expanded and moreover infinite models using this new type can be optimized. Let's try expanding the measure we already defined:
julia> expand(mref)
y(-0.556026876146)² + 0 z*y(-0.556026876146) - 2 y(-0.44438335711)*y(-0.556026876146) + 0 z² + 0 z*y(-0.44438335711) + y(-0.44438335711)²
Finally, as per recommendation let's make a wrapper method to make defining variance measures more convenient:
function variance(
expr::Union{JuMP.GenericAffExpr, GeneralVariableRef},
params::Union{GeneralVariableRef, AbstractArray{GeneralVariableRef}};
name::String = "Var",
num_supports::Int = 10,
use_existing::Bool = false
)
# get the supports
if use_existing
supps = supports.(params)
else
supps = generate_support_values(infinite_domain(first(params)),
num_supports = num_supports)
end
# make the data
data = DiscreteVarianceData(params, supps)
# built the measure
return measure(expr, data, name = name)
end
Notice in this case that we only permit linear expressions for expr
since it will be squared by our new measure and we currently only support quadratic expressions. (This could be overcome by defining a place holder variable for expr
.
Now let's use our constructor to repeat the above measure example:
julia> expand(variance(2y + z, xi, use_existing = true))
y(-0.556026876146)² + 0 z*y(-0.556026876146) - 2 y(-0.44438335711)*y(-0.556026876146) + 0 z² + 0 z*y(-0.44438335711) + y(-0.44438335711)²
We have done it! Now go and extend away!
Generative Support Information
As discussed in the Generative Supports section, generative supports help enable measure and/or derivative evaluation techniques that require the creation of generative supports (e.g., orthogonal collocation). Natively, we provide UniformGenerativeInfo
to help accomplish this which works for creating generative supports uniformly over finite elements as is the case for orthogonal collocation (note this includes scaling them as need to the size of each finite element). However, more complex generative support schemes can be enabled by defining a new concrete AbstractGenerativeInfo
subtype. This section will detail how this can be accomplished in InfiniteOpt
. A template for implementing this is provided in ./test/extensions/generative_info.jl
.
A new generative support information type can be created via the following:
- Define a concrete subtype of
AbstractGenerativeInfo
(required) - Make a unique support label that inherits
InternalLabel
(recommended) - Extend
InfiniteOpt.support_label
(required) - Extend
InfiniteOpt.make_generative_supports
(required).
For the sake of example, let's suppose we want to make a method that generates a certain amount of random supports for each finite element. First, let's define our struct RandomGenerativeInfo
:
using InfiniteOpt, Random
struct RandomGenerativeInfo <: InfiniteOpt.AbstractGenerativeInfo
amount::Int # amount of random supports per finite element
end
With that done, let's define a unique support label RandomInternal
for these types of supports and extend support_label
:
struct RandomInternal <: InternalLabel end
function InfiniteOpt.support_label(info::RandomGenerativeInfo)
return RandomInternal
end
Finally, let's extend make_generative_supports
to create a vector of the generative supports based on a RandomGenerativeInfo
and the existing model supports which are passed in the function as input:
function InfiniteOpt.make_generative_supports(info::RandomGenerativeInfo, pref, supps)
num_existing = length(supps)
num_existing <= 1 && error("`$pref` doesn't have enough supports.")
num_internal = info.attr
gen_supps = Float64[]
for i = 1:num_existing-1
lb = supps[i]
ub = supps[i+1]
append!(gen_supps, rand(num_internal) * (ub - lb) .+ lb)
end
return gen_supps
end
Our extension is done and now RandomGenerativeInfo
can be incorporated by a GenerativeDerivativeMethod
we create or an AbstractMeasureData
object of our choice like FunctionalDiscreteMeasureData
.
Transformation Backends
InfiniteOpt
provides a convenient interface and abstraction for modeling infinite-dimensional optimization problems. By default, InfiniteModel
s are reformulated into a solvable JuMP.Model
via TranscriptionOpt.TranscriptionBackend
which discretizes the model in accordance with the infinite parameter supports. However, users may wish to employ some other transformation method to produce the transformation backend. This section will explain how this can be done in InfiniteOpt
. A template for implementing this extension is provided in ./test/extensions/backend.jl
. Our default sub-module InfiniteOpt.TranscriptionOpt
also serves as a good example.
A new transformation approach and its corresponding transformation backend can be extended using the following steps:
- Define a
mutable struct
for variable/constraint mappings and other needed info (required) - Define an
AbstractTransformationBackend
(required) - Extend
Base.empty!
for the backend (required) - Extend
build_transformation_backend!
(required) - If appropriate and NOT a
JuMPBackend
, extend the following:transformation_model
transformation_data
JuMP.set_attribute
(including the suggested attributes)JuMP.get_attribute
(including the suggested attributes)JuMP.optimize!
JuMP.set_optimizer
JuMP.bridge_constraints
JuMP.add_bridge
JuMP.print_active_bridges
JuMP.print_active_bridges
JuMP.compute_conflict!
JuMP.copy_conflict
JuMP.mode
JuMP.backend
JuMP.unsafe_backend
- Extend the following, if possible (also enables result queries for
JuMPBackend
s): - Extend the following, if appropriate:
- As appropriate and if NOT a
JuMPBackend
, extend the following:- The remaining result related attributes listed in
JuMP.get_attribute
JuMP.lp_sensitivity_report
- The remaining result related attributes listed in
- If Step 6 was skipped and/or the backend is NOT a
JuMPBackend
then extend the following:InfiniteOpt.map_value
(enablesJuMP.value
)InfiniteOpt.map_infinite_parameter_value
(enablesJuMP.value
for infinite parameters)InfiniteOpt.map_optimizer_index
(enablesJuMP.optimizer_index
)InfiniteOpt.map_reduced_cost
(enablesJuMP.reduced_cost
)InfiniteOpt.map_shadow_price
(enablesJuMP.shadow_cost
)InfiniteOpt.map_dual
(enablesJuMP.dual
)
- Extend
InfiniteOpt.add_point_variable
andInfiniteOpt.add_semi_infinite_variable
to useexpand_measure
without modifying the infinite model.
This may seem like a lot a work, but the majority of the above steps can be skipped for JuMPBackend
s as exemplified below. A complete extension, showing all the above is provided in the extension template file.
For the sake of example, let's suppose we want to define a reformulation method for InfiniteModel
s that are 2-stage stochastic programs (i.e., only DistributionDomain
s are used, infinite variables are random 2nd stage variables, and finite variables are 1st stage variables). In particular, let's make a simple method that replaces the infinite parameters with their mean values, giving us the deterministic mean-valued problem.
First, let's define the (potentially mutable) struct
that will be used to store our variable and constraint mappings. This case it is quite simple since our deterministic model will have a 1-to-1 mapping:
using InfiniteOpt, Distributions
struct DeterministicData
# variable and constraint mapping
infvar_to_detvar::Dict{GeneralVariableRef, VariableRef}
infconstr_to_detconstr::Dict{InfOptConstraintRef, ConstraintRef}
# constructor
function DeterministicData()
return new(Dict{GeneralVariableRef, VariableRef}(),
Dict{InfOptConstraintRef, ConstraintRef}())
end
end
Now let's define the transformation backend based on JuMPBackend
that will use a tag Deterministic
:
struct Deterministic <: AbstractJuMPTag end
const DeterministicBackend = JuMPBackend{Deterministic, Float64, DeterministicData}
# Constructor
function DeterministicBackend(; kwargs...)
return JuMPBackend{Deterministic}(Model(; kwargs...), DeterministicData())
end
function DeterministicBackend(optimizer_constructor; kwargs...)
backend = DeterministicBackend(; kwargs...)
set_optimizer(backend.model, optimizer_constructor)
return backend
end
With the constructor we can now specify that a given InfiniteModel
uses a DeterministicBackend
instead of a TranscriptionBackend
or via set_transformation_backend
:
using Ipopt
# Make model using Ipopt and DeterministicModels
dbackend = DeterministicBackend(optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0))
model = InfiniteModel(dbackend)
# Or equivalently
model = InfiniteModel()
set_transformation_backend(model, DeterministicBackend())
set_optimizer(model, optimizer_with_attributes(Ipopt.Optimizer, "print_level" => 0))
Now model
uses a DeterministicBackend
as its transformation backend! With that we can build our InfiniteModel
as normal, for example:
@infinite_parameter(model, ξ ~ Uniform())
@variable(model, y[1:2] >= 0, Infinite(ξ))
@variable(model, z)
@objective(model, Min, z + expect(y[1] + y[2], ξ))
@constraint(model, 2y[1] - z <= 42)
@constraint(model, y[2]^2 + ξ == 2)
@constraint(model, sin(z) >= -1)
print(model)
# output
Min z + 𝔼{ξ}[y[1](ξ) + y[2](ξ)]
Subject to
y[1](ξ) ≥ 0.0, ∀ ξ ~ Uniform
y[2](ξ) ≥ 0.0, ∀ ξ ~ Uniform
2 y[1](ξ) - z ≤ 42.0, ∀ ξ ~ Uniform
y[2](ξ)² + ξ = 2.0, ∀ ξ ~ Uniform
sin(z) - -1 ≥ 0.0
We have defined our InfiniteModel
, but now we need to specify how to reformulate it into a DeterministicBackend
. This is accomplished by extending build_transformation_backend!
which will enable the use of optimize!
. A necessary preliminary step though, is to define Base.empty!
for DeterministicData
:
function Base.empty!(data::DeterministicData)
empty!(data.infvar_to_detvar)
empty!(data.infconstr_to_detconstr)
return data
end
This enables the backend to be cleared out before it is rebuilt which is necessary to allow for modifications to the model. Now, let's define an internal function _make_expression
that will use dispatch to convert an InfiniteOpt
expression into a JuMP
expression using the mappings stored in backend
's DeterministicData
:
## Make dispatch methods for converting InfiniteOpt expressions
# GeneralVariableRef
function _make_expression(backend::DeterministicBackend, expr::GeneralVariableRef)
return _make_expression(backend, expr, index(expr))
end
# IndependentParameterRef
function _make_expression(
backend::DeterministicBackend,
expr::GeneralVariableRef,
::IndependentParameterIndex
)
return mean(infinite_domain(expr).distribution) # assuming univariate
end
# FiniteParameterRef
function _make_expression(
backend::DeterministicBackend,
expr::GeneralVariableRef,
::FiniteParameterIndex
)
return parameter_value(expr)
end
# DependentParameterRef
function _make_expression(
backend::DeterministicBackend,
expr::GeneralVariableRef,
::DependentParameterIndex
)
return mean(infinite_domain(expr).distribution) # assuming valid distribution
end
# DecisionVariableRef
function _make_expression(
backend::DeterministicBackend,
expr::GeneralVariableRef,
::Union{InfiniteVariableIndex, FiniteVariableIndex}
)
return backend.data.infvar_to_detvar[expr]
end
# MeasureRef --> assume is expectation
function _make_expression(
backend::DeterministicBackend,
expr::GeneralVariableRef,
::MeasureIndex
)
return _make_expression(backend, measure_function(expr))
end
# AffExpr/QuadExpr/NonlinearExpr
function _make_expression(backend::DeterministicBackend, expr::Union{GenericAffExpr, GenericQuadExpr, GenericNonlinearExpr})
return map_expression(v -> _make_expression(backend, v), expr)
end
For simplicity in example, above we assume that only DistributionDomain
s are used, there are not any PointVariableRef
s, and all MeasureRef
s correspond to expectations. Naturally, a full extension should include checks to enforce that such assumptions hold. Notice that map_expression
is useful for converting expressions.
Now let's extend build_transformation_backend!
for DeterministicBackend
s. This should build out the backend in-place and thus we should also be sure to have it clear out any previous builds with Base.empty!
:
function InfiniteOpt.build_transformation_backend!(
model::InfiniteModel,
backend::DeterministicBackend
)
# TODO check that `model` is a stochastic model
# empty the model for a build/rebuild
empty!(backend)
backend.model.operator_counter = 0 # clears out any previous user defined operators
# add user-defined nonlinear operators if there are any
add_operators_to_jump(backend.model, model)
# add variables
for vref in all_variables(model)
if index(vref) isa InfiniteVariableIndex
start = NaN # simple hack for sake of example
else
start = start_value(vref)
start = isnothing(start) ? NaN : start
end
lb = has_lower_bound(vref) ? lower_bound(vref) : NaN
ub = has_upper_bound(vref) ? upper_bound(vref) : NaN
if is_fixed(vref)
lb = fix_value(vref)
end
info = VariableInfo(!isnan(lb), lb, !isnan(ub), ub, is_fixed(vref), lb,
!isnan(start), start, is_binary(vref), is_integer(vref))
new_vref = add_variable(backend.model, ScalarVariable(info), name(vref))
backend.data.infvar_to_detvar[vref] = new_vref
end
# add the objective
obj_func = _make_expression(backend, objective_function(model))
set_objective(backend.model, objective_sense(model), obj_func)
# add the constraints
for cref in all_constraints(model, Union{GenericAffExpr, GenericQuadExpr, GenericNonlinearExpr})
constr = constraint_object(cref)
new_func = _make_expression(backend, constr.func)
new_constr = build_constraint(error, new_func, constr.set)
new_cref = add_constraint(backend.model, new_constr, name(cref))
backend.data.infconstr_to_detconstr[cref] = new_cref
end
return
end
Note that Step 5 can be skipped since we are using the JuMPBackend
API which inherits all the needed methods. Now we can build our backend automatically and enable the use of optimize!
:
optimize!(model)
print(transformation_model(model))
# output
Min z + y[1] + y[2]
Subject to
sin(z) - -1.0 ≥ 0
2 y[1] - z ≤ 42
y[2]² = 1.5
y[1] ≥ 0
y[2] ≥ 0
Note that better variable naming could be used with the reformulated infinite variables. Moreover, in general extensions of build_transformation_backend!
should account for the possibility that InfiniteModel
contains constraints with DomainRestrictions
as accessed via domain_restrictions
.
Now that we have optimized our InfiniteModel
via the use the of a DeterministicBackend
, we probably want to access the results. All queries are enabled via Step 6 where we extend:
to return the variable(s)/expression(s)/constraint(s) in the backend. These will use the DeterministicData
and should error if no mapping can be found.
function InfiniteOpt.transformation_variable(
vref::GeneralVariableRef,
backend::DeterministicBackend
)
map_dict = backend.data.infvar_to_detvar
haskey(map_dict, vref) || error("Variable $vref not used in the transformation backend.")
return map_dict[vref]
end
function InfiniteOpt.transformation_expression(
expr::JuMP.AbstractJuMPScalar,
backend::DeterministicBackend
)
return _make_expression(backend, expr)
end
function InfiniteOpt.transformation_constraint(
cref::InfOptConstraintRef,
backend::DeterministicBackend
)
map_dict = backend.data.infconstr_to_detconstr
haskey(map_dict, cref) || error("Constraint $cref not used in the transformation backend.")
return map_dict[cref]
end
With these extensions we can now access all the result queries (skipping Steps 8-9). For example:
julia> termination_status(model)
LOCALLY_SOLVED::TerminationStatusCode = 4
julia> result_count(model)
1
julia> value.(y)
2-element Vector{Float64}:
-9.164638781941642e-9
1.224744871391589
julia> optimizer_index(z)
MathOptInterface.VariableIndex(3)
We also skip steps 7 and 10 since these are not applicable to this particular example.
That's it!
Wrapper Packages
InfiniteOpt
provides a convenient modular interface for defining infinite dimensional optimization problems, implementing many tedious JuMP
extensions such as facilitating mixed variable expressions. Thus, InfiniteOpt
can serve as a base package for specific types of infinite dimensional problems and/or solution techniques. These extension packages can implement any of the extensions shown above and likely will want to introduce wrapper functions and macros to use package specific terminology (e.g., using random variables instead of infinite variables).
Please reach out to us via the discussion forum to discuss your plans before starting this on your own.