Implement "parametric" Expressions

[sstroemer]

This is a pretty complete (but rather basic) implementation of parametric Expressions. It currently uses the JuMP.ParameterRef approach. Tests are passing, a new example has been added. I probably missed something somewhere, but it should not break too many non-parametric things.

I would consider this highly experimental for now, so I'd like to "just get that in" and see if we can get some actual application feedback from the MGA topic or from some rolling window optimization (or MPC) to improve on it. Felt like there will be some learnings that we probably can't easily "think through" without just playing around with it.

---

@daschw What are your thoughts on the naming for the public functions modify!(expression, new_value) and query(expression)?

I do not really want to use "parameter" anywhere (since we already have model parameters, etc., and I fear that people will start mixing them up). My first thought was calling these "unknowns", with set_unknown!(expression, new_value) and get_unknown(expression), but I thought "modify" is more expressive. Wasn't super happy with "query" then, but did not want to go for "access" or "get" (because we are already using those for other stuff) and "query" was the only reasonable expressive thing I came up with...

Motivation

This is based on example 07, with 8760 snapshots. Comparing "generate and solve each time" versus "pre-generate once" with "modify and re-solve each time".

Non-parametric:

Benchmark: 3 samples with 5 evaluations
        824.651 ms (13720309.20 allocs: 598.877 MiB, 36.25% gc time)
        841.027 ms (13720331.60 allocs: 598.878 MiB, 37.79% gc time)
        892.024 ms (13720400.40 allocs: 598.882 MiB, 42.72% gc time)

Parametric:

Benchmark: 4 samples with 5 evaluations
        514.215 ms (12063866 allocs: 386.921 MiB, 7.79% gc time)
        571.527 ms (12063866 allocs: 386.921 MiB, 16.12% gc time)
        575.096 ms (12063866 allocs: 386.921 MiB, 16.26% gc time)
        579.452 ms (12063866 allocs: 386.921 MiB, 16.84% gc time)

Approximately 120 ms are spent in HiGHS solving the model.
Therefore, roughly 730 ms versus 440 ms for the IESopt part, resulting in a 40% speedup.
Note: Around 150 ms are spent in result extraction, so disabling this ups that to a 50% speedup (580 ms vs. 290 ms).

Rough idea

Doing this

config:
  optimization:
    problem_type: PARAMETRIC+LP

allows us to prepare a proper objective function (QuadExpr) assuming that some Parameters might be added to it.

The following then allows to parameterize some scalar settings, while providing default values:

my_decision:
  type: Decision
  lb: $(0)
  ub: $(10)
  cost: $(0)

One can use cost: $() to "not" provide a default (since we still need one internally "no" always means 0).

Similarly for temporal expression the following works:

demand:
  type: Profile
  carrier: electricity
  node_from: node
  value: $(col@file)

Again, value: $(t) acts as the "no" default syntax for a temporal expression.

(Re-)optimizing

The approach then looks like this:

using IESopt

# Construct and solve the default config.
model = generate!("config.iesopt.yaml")
optimize!(model)

# Change the cost of the Decision and re-solve.
modify!(get_component(model, "my_decision").cost, 100)
optimize!(model)

# Extract the (current) value of the demand.
ts = query(get_component(model, "demand").value)

# Set the demand to "static 15" (a temporal parametric expression does not need
# a vector-valued setter value), and re-solve.
modify!(get_component(model, "demand").value, 15.0)
optimize!(model)

# Set the demand to the initial profile time series, scaled with "+ 25%", and re-solve.
modify!(get_component(model, "demand").value, 1.25 .* ts)
optimize!(model)

The "nice" thing behind that is, that we can directly work with the actual Expression objects, which can be accessed from Python in exactly the same way. No need to "register" parameters in the model, look them up in any way, etc.

---

The model also prints pretty nicely now (shortened the output to show some lines as example):

Min 100 co2_emissions.aux_value[1] + 50 create_gas.aux_value[1]
Subject to
 node1.nodalbalance : -demand1.par.value[1] + plant_gas.conversion[1] - conn.flow[1] = 0
 ...
 build.value_lb : -build.par.lb + build.value ≥ 0
 ...
 conn.flow_lb[1] : conn.par.capacity + conn.flow[1] ≥ 0
 ...
 plant_wind.conversion_ub[1] : -10 plant_wind.par.availability_factor[1] + plant_wind.conversion[1] ≤ 0
 ...
 demand2.par.value[1] ∈ MathOptInterface.Parameter{Float64}(4.88)
 ...
 plant_wind.par.availability_factor[1] ∈ MathOptInterface.Parameter{Float64}(0.88)
 conn.par.capacity ∈ MathOptInterface.Parameter{Float64}(5.0)

Note that the .par. is just cosmetic, to be in line with the other things. There is no actual build.par.lb entry in the "opt container dict", because it's just some_component.lb (the actual field) that already contains this. For stuff like Profiles that might be mistaken for the actual "value" (the variable) in the print, which is called foo.var.value. I felt that foo.par.value is better than just foo.value.

What's missing?

Updates to the documentation and the Python wrapper to expose the functions.

Quick notes on issues to keep in mind

• set = JuMP.Parameter(0) is not properly supported by HiGHS (and others?) if using a direct model. Using fixed variables should in theory still be possible, so we should investigate further.
• It may be, that there is a more efficient way than using an optimize hook to "fix" the quadratic objective function pre-solver (the main issue is that HiGHS does not detect that it's actually an LP, but maybe we can report that?).
• There's an old bridge-based proof of concept hidden somewhere in a JuMP (or MOI?) PR from roughly two years ago; that could help with the currently still existing problem of quadratic constraints (which occur for various Unit related configurations).
• The implementation is pretty general for any Expression, but I've only tested it with the "common" configurations and only partially for Units. It can be expected that unit commitment or similar more complex formulations may be bricked by that (note: I've tested capacity and availability_factor as the most common ones and they work).
• Conversion expressions are currently not supported at all (but they might need a slightly different approach anyways).
• A fully parametric model might be really easy to "docify" - we should keep track of that opportunity.

---

Auto-generated summary

This pull request introduces significant changes to support parametric expressions and optimization in the IESopt package. The most important changes include the addition of parametric capabilities, modifications to the objective function handling, and updates to various components to support the new features.

Parametric Capabilities:

• Added support for parametric expressions in the configuration file assets/examples/51_parametric_expressions.iesopt.yaml.
• Introduced the _is_parametric function to check if a model is parametric.
• Implemented the _optimize_hook_parametric function to handle quadratic objectives for parametric models.

Objective Function Handling:

• Modified the _build_model! function to include a hook for parametric models and handle the current objective accordingly.
• Updated the _optimize! function to include a TODO for substituting quadratic objectives with affine expressions in parametric models.

Component Updates:

• Changed the Decision component to use expressions for lb, ub, and cost fields.
• Added a new function _decision_con_value_bounds! to handle value bounds constraints for decisions.
• Updated the Profile component to handle parametric costs and values. [1] [2]

Expression Handling:

• Added functions to modify and query parametric expressions, and updated the _convert_to_expression function to handle parametric cases. [1] [2]
• Updated the _prepare and access functions to support parametric expressions. [1] [2]

Miscellaneous:

• Removed the string_names_on_creation configuration check from the _build_model! function.
• Updated various functions to handle parametric expressions and costs, including _connection_obj_cost! and _decision_obj_fixed!. [1] [2]

[sstroemer]

Just leaving a small ping @hkoller: this could be really useful for re-solving when running a lot of LPs.

[daschw]

Sorry for my very late response. I could not really think of anything better than modify! and query so far. It would be nice to have some set_...! and get_... variant but I agree that parameter would be confusing to our users.