Implement "parametric" Expressions

assets/examples/51_parametric_expressions.iesopt.yaml

@@ -0,0 +1,92 @@
+config:
+  general:
+    version:
+      core: 2.7.0
+  optimization:
+    problem_type: PARAMETRIC+LP
+    snapshots:
+      count: 168
+    solver:
+      name: highs
+  results:
+    enabled: false
+  files:
+    data: example_data.csv
+
+carriers:
+  electricity: {}
+  gas: {}
+  co2: {}
+
+components:
+  node1:
+    type: Node
+    carrier: electricity
+
+  node2:
+    type: Node
+    carrier: electricity
+    has_state: true
+    state_lb: $(0)
+    state_ub: $(10)
+
+  conn:
+    type: Connection
+    capacity: $(5)
+    node_from: node1
+    node_to: node2
+
+  # We could also use: `$(10) out:electricity` to model a parametric capacity.
+  plant_wind:
+    type: Unit
+    outputs: {electricity: node2}
+    conversion: ~ -> 1 electricity
+    capacity: 10 out:electricity
+    availability_factor: $(ex07_plant_wind_availability_factor@data)
+
+  plant_gas:
+    type: Unit
+    inputs: {gas: gas_grid}
+    outputs: {electricity: node1, co2: total_co2}
+    conversion: 1 gas -> 0.4 electricity + 0.2 co2
+    capacity: build:value out:electricity
+
+  build:
+    type: Decision
+    lb: $(0)
+    ub: $(10)
+    cost: $(0)
+
+  demand1:
+    type: Profile
+    carrier: electricity
+    node_from: node1
+    value: $(ex07_demand1_value@data)
+
+  demand2:
+    type: Profile
+    carrier: electricity
+    node_from: node2
+    value: $(ex07_demand2_value@data)
+
+  gas_grid:
+    type: Node
+    carrier: gas
+
+  total_co2:
+    type: Node
+    carrier: co2
+
+  create_gas:
+    type: Profile
+    carrier: gas
+    node_to: gas_grid
+    mode: create
+    cost: $(50)
+
+  co2_emissions:
+    type: Profile
+    carrier: co2
+    node_from: total_co2
+    mode: destroy
+    cost: $(100.)

src/IESopt.jl

@@ -20,12 +20,6 @@ include("docify/docify.jl")
 function _build_model!(model::JuMP.Model)
     @info "[build] Begin creating JuMP formulation from components"
 
-    if @config(model, general.performance.string_names, Bool) != model.set_string_names_on_creation
-        new_val = @config(model, general.performance.string_names, Bool)
-        @debug "Overwriting `string_names_on_creation` to `$(new_val)` based on config"
-        JuMP.set_string_names_on_creation(model, new_val)
-    end
-
     # This specifies the order in which components are built. This ensures that model parts that are used later on, are
     # already initialized (e.g. constructing a constraint may use expressions and variables).
     build_order = [
@@ -133,16 +127,52 @@ function _build_model!(model::JuMP.Model)
     end
 
     if !_is_multiobjective(model)
-        current_objective = @config(model, optimization.objective.current)
-        isnothing(current_objective) && @critical "[build] Missing an active objective"
-        @objective(model, Min, internal(model).model.objectives[current_objective].expr)
+        if !_is_parametric(model) || _is_qp(model)
+            current_objective = @config(model, optimization.objective.current)
+            isnothing(current_objective) && @critical "[build] Missing an active objective"
+            @objective(model, Min, internal(model).model.objectives[current_objective].expr)
+        else
+            current_objective = @config(model, optimization.objective.current)
+            isnothing(current_objective) && @critical "[build] Missing an active objective"
+
+            # Set an optimizer hook that helps resolving quadratic objectives for parametric models.
+            JuMP.set_optimize_hook(model, _optimize_hook_parametric)
+        end
     else
         @objective(
             model,
             Min,
             [internal(model).model.objectives[obj].expr for obj in @config(model, optimization.multiobjective.terms)]
         )
     end
+
+    return nothing
+end
+
+function _optimize_hook_parametric(model::JuMP.Model)
+    obj = internal(model).model.objectives[@config(model, optimization.objective.current)].expr
+
+    if obj isa JuMP.QuadExpr
+        @debug "[optimize > hook] Resolving quadratic objective"
+        obj_aff_aux = zero(JuMP.AffExpr)
+
+        # Add the affine part of the objective.
+        JuMP.add_to_expression!(obj_aff_aux, obj.aff)
+
+        # Add the quadratic part of the objective.
+        for (term, coeff) in obj.terms
+            if JuMP.is_parameter(term.a)
+                JuMP.add_to_expression!(obj_aff_aux, term.b, JuMP.parameter_value(term.a) * coeff)
+            else
+                JuMP.add_to_expression!(obj_aff_aux, term.a, JuMP.parameter_value(term.b) * coeff)
+            end
+        end
+
+        # Overwrite the objective.
+        @objective(model, Min, obj_aff_aux)
+    end
+
+    return JuMP.optimize!(model; ignore_optimize_hook=true)
 end
 
 function _prepare_model!(model::JuMP.Model)

src/config/config.jl

@@ -80,3 +80,5 @@ _has_representative_snapshots(model::JuMP.Model) = false  # TODO
 _is_multiobjective(model::JuMP.Model) = (:mo in @config(model, optimization.problem_type))::Bool
 _is_lp(model::JuMP.Model) = (:lp in @config(model, optimization.problem_type))::Bool
 _is_milp(model::JuMP.Model) = (:milp in @config(model, optimization.problem_type))::Bool
+_is_qp(model::JuMP.Model) = (:qp in @config(model, optimization.problem_type))::Bool
+_is_parametric(model::JuMP.Model) = (:parametric in @config(model, optimization.problem_type))::Bool

src/core/connection.jl

@@ -201,7 +201,7 @@ function _prepare!(connection::Connection)
 
         if _isempty(connection.capacity)
             # Only calculate capacity if it is not given by the user
-            connection.capacity = _convert_to_expression(model, connection.pf_V * connection.pf_I)
+            connection.capacity = _convert_to_expression(model, connection.pf_V * connection.pf_I, "")
         end
         # todo: convert B1
     end

src/core/connection/obj_cost.jl

@@ -23,12 +23,12 @@ function _connection_obj_cost!(connection::Connection)
 
     model = connection.model
 
-    connection.obj.cost = JuMP.AffExpr(0.0)
+    connection.obj.cost = _isparametric(connection.cost) ? zero(JuMP.QuadExpr) : zero(JuMP.AffExpr)
     for t in get_T(connection.model)
         JuMP.add_to_expression!(
             connection.obj.cost,
             connection.var.flow[t],
-            _weight(model, t) * access(connection.cost, t, Float64),
+            _weight(model, t) * access(connection.cost, t),
         )
     end
 

src/core/decision.jl

@@ -20,17 +20,17 @@ component's settings, as well as have associated costs.
     raw"""```{"mandatory": "no", "values": "numeric", "unit": "-", "default": "`0`"}```
     Minimum size of the decision value (considered for each "unit" if count allows multiple "units").
     """
-    lb::_OptionalScalarInput = 0
+    lb::Expression = @_default_expression(0.0)
 
     raw"""```{"mandatory": "no", "values": "numeric", "unit": "-", "default": "``+\\infty``"}```
     Maximum size of the decision value (considered for each "unit" if count allows multiple "units").
     """
-    ub::_OptionalScalarInput = nothing
+    ub::Expression = @_default_expression(nothing)
 
     raw"""```{"mandatory": "no", "values": "numeric", "unit": "monetary (per value)", "default": "`0`"}```
     Cost that the decision value induces, given as ``cost \cdot value``.
     """
-    cost::_OptionalScalarInput = nothing
+    cost::Expression = @_default_expression(nothing)
 
     raw"""```{"mandatory": "no", "values": "numeric", "unit": "-", "default": "-"}```
     If `mode: fixed`, this value is used as the fixed value of the decision. This can be useful if this `Decision` was
@@ -165,6 +165,7 @@ include("decision/con_fixed.jl")
 include("decision/con_sos_value.jl")
 include("decision/con_sos1.jl")
 include("decision/con_sos2.jl")
+include("decision/con_value_bounds.jl")
 include("decision/obj_fixed.jl")
 include("decision/obj_sos.jl")
 include("decision/obj_value.jl")
@@ -183,7 +184,9 @@ function _construct_constraints!(decision::Decision)
     _decision_con_fixed!(decision)
     _decision_con_sos_value!(decision)
     _decision_con_sos1!(decision)
-    return _decision_con_sos2!(decision)
+    _decision_con_sos2!(decision)
+    _decision_con_value_bounds!(decision)
+    return nothing
 end
 
 function _construct_objective!(decision::Decision)

src/core/decision/con_value_bounds.jl

@@ -0,0 +1,24 @@
+@doc raw"""
+    _decision_con_value_bounds!(decision::Decision)
+
+to be added
+"""
+function _decision_con_value_bounds!(decision::Decision)
+    if !_isempty(decision.lb)
+        decision.con.value_lb = @constraint(
+            decision.model,
+            decision.var.value >= access(decision.lb),
+            base_name = make_base_name(decision, "value_lb")
+        )
+    end
+
+    if !_isempty(decision.ub)
+        decision.con.value_ub = @constraint(
+            decision.model,
+            decision.var.value <= access(decision.ub),
+            base_name = make_base_name(decision, "value_ub")
+        )
+    end
+
+    return nothing
+end

src/core/decision/obj_fixed.jl

@@ -11,7 +11,7 @@ function _decision_obj_fixed!(decision::Decision)
 
     model = decision.model
 
-    decision.obj.fixed = JuMP.AffExpr(0.0)
+    decision.obj.fixed = _isparametric(profile.cost) ? zero(JuMP.QuadExpr) : zero(JuMP.AffExpr)
     if decision.mode === :sos1
         for i in eachindex(decision.sos)
             if haskey(decision.sos[i], "fixed_cost")

src/core/decision/obj_value.jl

@@ -13,7 +13,7 @@ function _decision_obj_value!(decision::Decision)
     end
 
     model = decision.model
-    decision.obj.value = decision.var.value * decision.cost
+    decision.obj.value = decision.var.value * access(decision.cost)
     push!(internal(model).model.objectives["total_cost"].terms, decision.obj.value)
 
     return nothing

src/core/decision/var_value.jl

@@ -17,17 +17,6 @@ function _decision_var_value!(decision::Decision)
 
     if decision.mode === :fixed
         JuMP.fix(decision.var.value, decision.fixed_value)
-    else
-        if !isnothing(decision.ub) && (decision.lb == decision.ub)
-            JuMP.fix(decision.var.value, decision.lb)
-        else
-            if !isnothing(decision.lb)
-                JuMP.set_lower_bound(decision.var.value, decision.lb)
-            end
-            if !isnothing(decision.ub)
-                JuMP.set_upper_bound(decision.var.value, decision.ub)
-            end
-        end
     end
 
     return nothing