Cleanup Multiplicative Parameters handling

src/MOI_wrapper.jl

@@ -74,10 +74,6 @@ function _cache_multiplicative_params!(
     for term in f.pv
         push!(model.multiplicative_parameters_pv, term.variable_1.value)
     end
-    for term in f.pp
-        push!(model.multiplicative_parameters_pp, term.variable_1.value)
-        push!(model.multiplicative_parameters_pp, term.variable_2.value)
-    end
     return
 end
 
@@ -91,15 +87,22 @@ function _cache_multiplicative_params!(
             term.scalar_term.variable_1.value,
         )
     end
-    for term in f.pp
-        push!(
-            model.multiplicative_parameters_pp,
-            term.scalar_term.variable_1.value,
-        )
-        push!(
-            model.multiplicative_parameters_pp,
-            term.scalar_term.variable_2.value,
-        )
+    return
+end
+
+function _cache_multiplicative_params!(
+    model::Optimizer{T},
+    f::ParametricCubicFunction{T},
+) where {T}
+    for term in f.pv
+        push!(model.multiplicative_parameters_pv, term.variable_1.value)
+    end
+    for term in f.pvv
+        push!(model.multiplicative_parameters_pv, term.index_1.value)
+    end
+    for term in f.ppv
+        push!(model.multiplicative_parameters_pv, term.index_1.value)
+        push!(model.multiplicative_parameters_pv, term.index_2.value)
     end
     return
 end
@@ -137,7 +140,6 @@ function MOI.is_empty(model::Optimizer)
            isempty(model.vector_affine_constraint_cache) &&
            #
            isempty(model.multiplicative_parameters_pv) &&
-           isempty(model.multiplicative_parameters_pp) &&
            isempty(model.dual_value_of_parameters) &&
            model.number_of_parameters_in_model == 0 &&
            isempty(model.parameters_in_conflict) &&
@@ -197,8 +199,6 @@ function MOI.empty!(model::Optimizer{T}) where {T}
     empty!(model.vector_affine_constraint_cache)
     # multiplicative_parameters_pv
     empty!(model.multiplicative_parameters_pv)
-    # multiplicative_parameters_pp
-    empty!(model.multiplicative_parameters_pp)
     # dual_value_of_parameters
     empty!(model.dual_value_of_parameters)
     # [SKIP] evaluate_duals

src/ParametricOptInterface.jl

@@ -186,7 +186,6 @@ mutable struct Optimizer{T,OT<:MOI.ModelLike} <: MOI.AbstractOptimizer
         ParametricVectorAffineFunction{T},
     }
     multiplicative_parameters_pv::Set{Int64}
-    multiplicative_parameters_pp::Set{Int64}
     dual_value_of_parameters::Vector{T}
     evaluate_duals::Bool
     number_of_parameters_in_model::Int64
@@ -254,8 +253,6 @@ mutable struct Optimizer{T,OT<:MOI.ModelLike} <: MOI.AbstractOptimizer
             DoubleDicts.DoubleDict{ParametricVectorAffineFunction{T}}(),
             # multiplicative_parameters_pv
             Set{Int64}(),
-            # multiplicative_parameters_pp
-            Set{Int64}(),
             # dual_value_of_parameters
             T[],
             # evaluate_duals

src/cubic_objective.jl

@@ -45,10 +45,13 @@ function MOI.set(
     # 5. Clear old caches
     _empty_objective_function_caches!(model)
 
-    # 6. Store new cache
+    # 6. Cache multiplicative parameters
+    _cache_multiplicative_params!(model, cubic_func)
+
+    # 7. Store new cache
     model.cubic_objective_cache = cubic_func
 
-    # 7. Store original for retrieval if option is enabled
+    # 8. Store original for retrieval if option is enabled
     if model.save_original_objective_and_constraints
         MOI.set(
             model.original_objective_cache,

src/duals.jl

@@ -16,6 +16,9 @@ function _compute_dual_of_parameters!(model::Optimizer{T}) where {T}
     if model.quadratic_objective_cache !== nothing
         _update_duals_from_objective!(model, model.quadratic_objective_cache)
     end
+    if model.cubic_objective_cache !== nothing
+        _update_duals_from_objective!(model, model.cubic_objective_cache)
+    end
     return
 end
 
@@ -116,6 +119,69 @@ function _update_duals_from_objective!(model::Optimizer{T}, pf) where {T}
     return
 end
 
+function _update_duals_from_objective!(
+    model::Optimizer{T},
+    pf::ParametricQuadraticFunction{T},
+) where {T}
+    is_min = MOI.get(model.optimizer, MOI.ObjectiveSense()) == MOI.MIN_SENSE
+    sign = ifelse(is_min, one(T), -one(T))
+    # p terms: ∂(c·p_i)/∂p_i = c
+    for term in pf.p
+        model.dual_value_of_parameters[p_val(term.variable)] +=
+            sign * term.coefficient
+    end
+    # pp terms: ∂(c·p_i·p_j)/∂p_i = c·p_j
+    for term in pf.pp
+        mult = sign * term.coefficient
+        if term.variable_1 == term.variable_2
+            mult /= 2
+        end
+        model.dual_value_of_parameters[p_val(term.variable_1)] +=
+            mult * model.parameters[p_idx(term.variable_2)]
+        model.dual_value_of_parameters[p_val(term.variable_2)] +=
+            mult * model.parameters[p_idx(term.variable_1)]
+    end
+    return
+end
+
+function _update_duals_from_objective!(
+    model::Optimizer{T},
+    pf::ParametricCubicFunction{T},
+) where {T}
+    is_min = MOI.get(model.optimizer, MOI.ObjectiveSense()) == MOI.MIN_SENSE
+    sign = ifelse(is_min, one(T), -one(T))
+    # p terms: ∂(c·p_i)/∂p_i = c
+    for term in pf.p
+        model.dual_value_of_parameters[p_val(term.variable)] +=
+            sign * term.coefficient
+    end
+    # pp terms: ∂(c·p_i·p_j)/∂p_i = c·p_j (diagonal: c/2·2·p_i = c·p_i)
+    for term in pf.pp
+        mult = sign * term.coefficient
+        if term.variable_1 == term.variable_2
+            mult /= 2
+        end
+        model.dual_value_of_parameters[p_val(term.variable_1)] +=
+            mult * model.parameters[p_idx(term.variable_2)]
+        model.dual_value_of_parameters[p_val(term.variable_2)] +=
+            mult * model.parameters[p_idx(term.variable_1)]
+    end
+    # ppp terms: ∂(c·p_i·p_j·p_k)/∂p_i = c·p_j·p_k
+    for term in pf.ppp
+        coef = sign * term.coefficient
+        p1_val = model.parameters[p_idx(term.index_1)]
+        p2_val = model.parameters[p_idx(term.index_2)]
+        p3_val = model.parameters[p_idx(term.index_3)]
+        model.dual_value_of_parameters[p_val(term.index_1)] +=
+            coef * p2_val * p3_val
+        model.dual_value_of_parameters[p_val(term.index_2)] +=
+            coef * p1_val * p3_val
+        model.dual_value_of_parameters[p_val(term.index_3)] +=
+            coef * p1_val * p2_val
+    end
+    return
+end
+
 function MOI.get(
     model::Optimizer{T},
     attr::MOI.ConstraintDual,

test/test_JuMP.jl

@@ -576,6 +576,80 @@ function test_jump_dual_objective_max()
     return
 end
 
+function test_jump_dual_objective_pv()
+    # p*x is multiplicative → dual query should error
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer))
+    set_silent(model)
+    @variable(model, x >= 0)
+    @variable(model, p in Parameter(-2.0))
+    @objective(model, Min, p * x + x^2)
+    optimize!(model)
+    @test_throws ErrorException dual(ParameterRef(p))
+    return
+end
+
+function test_jump_dual_objective_pp_same()
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer))
+    set_silent(model)
+    @variable(model, x >= 1)
+    @variable(model, p in Parameter(3.0))
+    @objective(model, Min, x + p^2)
+    optimize!(model)
+    # ∂(p²)/∂p = 2p = 6
+    @test dual(ParameterRef(p)) ≈ 6.0 atol = 1e-4
+    return
+end
+
+function test_jump_dual_objective_pp_different()
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer))
+    set_silent(model)
+    @variable(model, x >= 1)
+    @variable(model, p1 in Parameter(2.0))
+    @variable(model, p2 in Parameter(3.0))
+    @objective(model, Min, x + p1 * p2)
+    optimize!(model)
+    # ∂(p1*p2)/∂p1 = p2 = 3, ∂(p1*p2)/∂p2 = p1 = 2
+    @test dual(ParameterRef(p1)) ≈ 3.0 atol = 1e-4
+    @test dual(ParameterRef(p2)) ≈ 2.0 atol = 1e-4
+    return
+end
+
+function test_jump_dual_objective_mixed_terms()
+    # min x + 2p + p^2  s.t.  x >= p,  p = 3
+    # Optimal: x* = 3, obj = 3 + 6 + 9 = 18
+    # ∂f/∂p from obj = 2 + 2p = 8
+    # Constraint x - p >= 0: dual λ = 1, ∂g/∂p = -1, contribution = 1
+    # Total = 8 + 1 = 9
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer))
+    set_silent(model)
+    @variable(model, x)
+    @variable(model, p in Parameter(3.0))
+    @constraint(model, x >= p)
+    @objective(model, Min, x + 2 * p + p^2)
+    optimize!(model)
+    @test objective_value(model) ≈ 18.0 atol = 1e-4
+    @test dual(ParameterRef(p)) ≈ 9.0 atol = 1e-4
+    return
+end
+
+function test_jump_dual_objective_parameter_update()
+    # min x + p^2  s.t.  x >= 1,  p = 2 then p = 3
+    # ∂(p^2)/∂p = 2p
+    # p=2: dual = 4, p=3: dual = 6
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer))
+    set_silent(model)
+    @variable(model, x)
+    @variable(model, p in Parameter(2.0))
+    @constraint(model, x >= 1)
+    @objective(model, Min, x + p^2)
+    optimize!(model)
+    @test dual(ParameterRef(p)) ≈ 4.0 atol = 1e-4
+    set_parameter_value(p, 3.0)
+    optimize!(model)
+    @test dual(ParameterRef(p)) ≈ 6.0 atol = 1e-4
+    return
+end
+
 function test_jump_dual_multiple_parameters_1()
     model = Model(() -> POI.Optimizer(HiGHS.Optimizer))
     set_silent(model)

test/test_MathOptInterface.jl

@@ -1451,8 +1451,16 @@ function test_qp_objective_parameter_times_parameter()
         0.0,
         atol = ATOL,
     )
-    @test MOI.get(optimizer, MOI.ConstraintDual(), cy) == 0.0
-    @test MOI.get(optimizer, MOI.ConstraintDual(), cz) == 0.0
+    @test isapprox(
+        MOI.get(optimizer, MOI.ConstraintDual(), cy),
+        1.0,
+        atol = ATOL,
+    )
+    @test isapprox(
+        MOI.get(optimizer, MOI.ConstraintDual(), cz),
+        1.0,
+        atol = ATOL,
+    )
     MOI.set(optimizer, MOI.ConstraintSet(), cy, MOI.Parameter(2.0))
     MOI.optimize!(optimizer)
     @test isapprox(MOI.get(optimizer, MOI.ObjectiveValue()), 2.0, atol = ATOL)

test/test_cubic.jl

@@ -1623,6 +1623,79 @@ function test_jump_cubic_direct_model_ppp()
     return
 end
 
+# ============================================================================
+# Contribution of objective parameters in duals
+# ============================================================================
+
+function test_cubic_dual_ppp_terms()
+    # min x + p^3  s.t.  x >= 1,  p = 2
+    # Optimal: x = 1, obj = 1 + 8 = 9
+    # ∂(p^3)/∂p = 3p^2 = 12
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
+    set_silent(model)
+    @variable(model, x)
+    @variable(model, p in MOI.Parameter(2.0))
+    @constraint(model, x >= 1)
+    @objective(model, Min, x + p^3)
+    optimize!(model)
+    @test termination_status(model) in (OPTIMAL, LOCALLY_SOLVED)
+    @test dual(ParameterRef(p)) ≈ 12.0 atol = ATOL
+    return
+end
+
+function test_cubic_dual_ppv_terms()
+    # min p1*p2*x + x^2  s.t.  x >= 0,  p1 = 1, p2 = -2
+    # ppv is multiplicative → dual query should error
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
+    set_silent(model)
+    @variable(model, x >= 0)
+    @variable(model, p1 in MOI.Parameter(1.0))
+    @variable(model, p2 in MOI.Parameter(-2.0))
+    @objective(model, Min, p1 * p2 * x + x^2)
+    optimize!(model)
+    @test termination_status(model) in (OPTIMAL, LOCALLY_SOLVED)
+    @test value(x) ≈ 1.0 atol = ATOL
+    @test_throws ErrorException dual(ParameterRef(p1))
+    @test_throws ErrorException dual(ParameterRef(p2))
+    return
+end
+
+function test_cubic_dual_pvv_terms()
+    # min p*x^2 - 3x  s.t.  0 <= x <= 10,  p = 1
+    # pvv is multiplicative → dual query should error
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
+    set_silent(model)
+    @variable(model, 0 <= x <= 10)
+    @variable(model, p in MOI.Parameter(1.0))
+    @objective(model, Min, p * x^2 - 3 * x)
+    optimize!(model)
+    @test termination_status(model) in (OPTIMAL, LOCALLY_SOLVED)
+    @test value(x) ≈ 1.5 atol = ATOL
+    @test_throws ErrorException dual(ParameterRef(p))
+    return
+end
+
+function test_cubic_dual_ppp_three_distinct()
+    # min x + p1*p2*p3  s.t.  x >= 1,  p1=2, p2=3, p3=4
+    # ∂/∂p1 = p2*p3 = 12
+    # ∂/∂p2 = p1*p3 = 8
+    # ∂/∂p3 = p1*p2 = 6
+    model = Model(() -> POI.Optimizer(HiGHS.Optimizer()))
+    set_silent(model)
+    @variable(model, x)
+    @variable(model, p1 in MOI.Parameter(2.0))
+    @variable(model, p2 in MOI.Parameter(3.0))
+    @variable(model, p3 in MOI.Parameter(4.0))
+    @constraint(model, x >= 1)
+    @objective(model, Min, x + p1 * p2 * p3)
+    optimize!(model)
+    @test termination_status(model) in (OPTIMAL, LOCALLY_SOLVED)
+    @test dual(ParameterRef(p1)) ≈ 12.0 atol = ATOL
+    @test dual(ParameterRef(p2)) ≈ 8.0 atol = ATOL
+    @test dual(ParameterRef(p3)) ≈ 6.0 atol = ATOL
+    return
+end
+
 end  # module
 
 TestCubic.runtests()