martinbiel / StochasticPrograms.jl

Julia package for formulating and analyzing stochastic recourse models.
MIT License
75 stars 25 forks source link

How to stop the optimization of a stochastic program #34

Closed enri07 closed 2 years ago

enri07 commented 3 years ago

Hi,

I'm having some troubles trying to interrupt the optimization of a model that I have built by using StochasticPrograms.

After a lot of test I decided to use progressive-hedging algorithms to solve the model, but it is still difficult to find the optimal solution. Thus, I decided to put an optimality gap to stop the execution, but it seem that it doesn't work. This is what I put in the code:

op_model_user = instantiate(model_user, possible_scenarios, optimizer = ProgressiveHedging.Optimizer)
set_optimizer_attribute(op_model_user, SubProblemOptimizer(), CPLEX.Optimizer)
set_optimizer_attribute(op_model_user, PrimalTolerance(), 5e-1)
MOI.set(op_model_user, RawSubProblemOptimizerParameter("CPX_PARAM_EPGAP"), 5e-1)

optimize!(op_model_user)

These are some lines of the optimization after 30 minutes from the start:

CPXPARAM_MIP_Tolerances_MIPGap                   0.10000000000000001
1 of 1 MIP starts provided solutions.
MIP start 'm1' defined initial solution with objective 1.4695e+10.
Tried aggregator 2 times.
MIQCP Presolve eliminated 11393 rows and 11395 columns.
Aggregator did 1078 substitutions.
Reduced MIQCP has 29012 rows, 45200 columns, and 93954 nonzeros.
Reduced MIQCP has 0 binaries, 8 generals, 0 SOSs, and 0 indicators.
Reduced MIQCP has 1 quadratic constraints.
Presolve time = 0.10 sec. (78.89 ticks)

Root node processing (before b&c):
  Real time             =    4.34 sec. (2335.55 ticks)
Parallel b&c, 16 threads:
  Real time             =    0.00 sec. (0.00 ticks)
  Sync time (average)   =    0.00 sec.
  Wait time (average)   =    0.00 sec.
                          ------------
Total (root+branch&cut) =    4.34 sec. (2335.55 ticks)
Version identifier: 20.1.0.0 | 2020-11-10 | 9bedb6d68
CPXPARAM_MIP_Tolerances_MIPGap                   0.10000000000000001
1 of 1 MIP starts provided solutions.
MIP start 'm1' defined initial solution with objective 1.5373e+09.
Tried aggregator 2 times.
MIQCP Presolve eliminated 11393 rows and 11392 columns.
Aggregator did 1078 substitutions.
Reduced MIQCP has 29013 rows, 45203 columns, and 93954 nonzeros.
Reduced MIQCP has 0 binaries, 8 generals, 0 SOSs, and 0 indicators.
Reduced MIQCP has 1 quadratic constraints.
Presolve time = 0.10 sec. (78.90 ticks)

Root node processing (before b&c):
  Real time             =    4.09 sec. (2202.25 ticks)
Parallel b&c, 16 threads:
  Real time             =    0.00 sec. (0.00 ticks)
  Sync time (average)   =    0.00 sec.
  Wait time (average)   =    0.00 sec.
                          ------------
Total (root+branch&cut) =    4.09 sec. (2202.25 ticks)

Moreover, sometimes the status of the solution appears and it seems that the Primal Gap is very good:

Schermata

Why execution is still running? Is there something I can do to interrupt it? If I use Ctrl+C it stops but there aren't solutions saved, so I can't have any information. I also tried to set a time limit, but also in this case I'm not able to find which is the optimal solution found.

Thank you very much for your time.

martinbiel commented 3 years ago

Your DualTolerance is probably still set to the default value of 1e-5 so since your dual gap is still large the procedure continues. You could set the dual tolerance to some really large value, or evenInf, if you want to inspect the solution you get with the nice primal convergence, although you should not get an implementable solution with agreement over all scenarios since there is no dual progress.

enri07 commented 3 years ago

Thank you very much for your reply, now it works.

Do you know if is it possible to use the macro @expression for the deterministic equivalent? Because if I run the same program with the LShaped optimizer or with the ProgressiveHedging optimizer it works, but if I try to convert the model into the Determinist Equivalent problem it shows me the following error:

ERROR: LoadError: MethodError: no method matching name(::DecisionAffExpr{Float64})
Closest candidates are:
  name(::DecisionRef) at /home/user/.julia/packages/StochasticPrograms/SjnV4/src/types/decisions/variable_interface.jl:256
  name(::DecisionVariable) at /home/user/.julia/packages/StochasticPrograms/SjnV4/src/types/decision_variable.jl:259
  name(::VariableRef) at /home/user/.julia/packages/JuMP/Xrr7O/src/variables.jl:307
  ...
Stacktrace:
 [1] generate!(::StochasticProgram{2,Tuple{StochasticPrograms.Stage{NamedTuple{(),Tuple{}}},StochasticPrograms.Stage{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}},DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}}, ::DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}, ::Int64) at /home/user/.julia/packages/StochasticPrograms/SjnV4/src/methods/deterministic_equivalent/generation.jl:106
 [2] generate!(::StochasticProgram{2,Tuple{StochasticPrograms.Stage{NamedTuple{(),Tuple{}}},StochasticPrograms.Stage{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}},DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}}, ::DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}) at /home/user/.julia/packages/StochasticPrograms/SjnV4/src/methods/deterministic_equivalent/generation.jl:12
 [3] generate!(::StochasticProgram{2,Tuple{StochasticPrograms.Stage{NamedTuple{(),Tuple{}}},StochasticPrograms.Stage{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}},DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}}) at /home/user/.julia/packages/StochasticPrograms/SjnV4/src/methods/generation.jl:145
 [4] instantiate(::StochasticModel{2,Tuple{StageParameters{NamedTuple{(),Tuple{}}},StageParameters{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}}}, ::Array{Scenario_Load_Renewable,1}; instantiation::UnspecifiedInstantiation, optimizer::Type{T} where T, defer::Bool, direct_model::Bool, kw::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/user/.julia/packages/StochasticPrograms/SjnV4/src/methods/api.jl:37
 [5] top-level scope at /home/user/File_1/Single_user_stoch_nogen.jl:256
 [6] include(::String) at ./client.jl:439
 [7] top-level scope at REPL[1]:1
in expression starting at /home/user/File_1/Single_user_stoch_nogen.jl:256

Again, thank you very much for your replies and your time!

martinbiel commented 3 years ago

No, that should be a fixable oversight. Could you give a MWE of the @expression lines that do not work for the deterministic equivalent?

martinbiel commented 3 years ago

Nevermind, I found it. It should work on master now.

enri07 commented 3 years ago

Thank you very much!

enri07 commented 2 years ago

I'm sorry to re-open this issue, but I found out a new error message which I think it's a consequence of the previous change in the code. Now, if I declare an @expression using the deterministic equivalent in master mode it shows the following error message:

ERROR: LoadError: An object of name P_tot_P_overestimate is already attached to this model. If this
    is intended, consider using the anonymous construction syntax, e.g.,
    `x = @variable(model, [1:N], ...)` where the name of the object does
    not appear inside the macro.

    Alternatively, use `unregister(model, :P_tot_P_overestimate)` to first unregister
    the existing name from the model. Note that this will not delete the
    object; it will just remove the reference at `model[:P_tot_P_overestimate]`.

Stacktrace:
 [1] error(::String) at ./error.jl:33
 [2] _error_if_cannot_register at /home/user/.julia/packages/JuMP/Xrr7O/src/macros.jl:100 [inlined]
 [3] macro expansion at /home/user/.julia/packages/JuMP/Xrr7O/src/macros.jl:134 [inlined]
 [4] #169 at /home/user/Modello_1/Single_user_stoch_nogen.jl:131 [inlined]
 [5] (::var"#169#338")(::Model, ::NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}, ::Scenario_Load_Renewable) at ./none:0
 [6] generate!(::StochasticProgram{2,Tuple{StochasticPrograms.Stage{NamedTuple{(),Tuple{}}},StochasticPrograms.Stage{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}},DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}}, ::DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}, ::Int64) at /home/user/.julia/packages/StochasticPrograms/B9Ra4/src/methods/deterministic_equivalent/generation.jl:63
 [7] generate!(::StochasticProgram{2,Tuple{StochasticPrograms.Stage{NamedTuple{(),Tuple{}}},StochasticPrograms.Stage{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}},DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}}, ::DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}) at /home/user/.julia/packages/StochasticPrograms/B9Ra4/src/methods/deterministic_equivalent/generation.jl:33
 [8] generate!(::StochasticProgram{2,Tuple{StochasticPrograms.Stage{NamedTuple{(),Tuple{}}},StochasticPrograms.Stage{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}},DeterministicEquivalent{2,1,Tuple{Array{Scenario_Load_Renewable,1}}}}) at /home/user/.julia/packages/StochasticPrograms/B9Ra4/src/methods/generation.jl:167
 [9] instantiate(::StochasticModel{2,Tuple{StageParameters{NamedTuple{(),Tuple{}}},StageParameters{NamedTuple{(:price_MG_P, :price_MG_N_V, :price_MG_N_F, :price_F, :price_IMB, :c_UP),Tuple{Dict{Int64,Float64},Dict{Int64,Float64},Dict{Int64,Float64},Float64,Dict{Int64,Float64},Dict{String,Float64}}}}}}, ::Array{Scenario_Load_Renewable,1}; instantiation::UnspecifiedInstantiation, optimizer::Type{T} where T, defer::Bool, direct_model::Bool, kw::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}) at /home/user/.julia/packages/StochasticPrograms/B9Ra4/src/methods/api.jl:59
 [10] top-level scope at /home/user/Modello_1/Single_user_stoch_nogen.jl:265
 [11] include(::String) at ./client.jl:439
 [12] top-level scope at REPL[1]:1
in expression starting at /home/user/Modello_1/Single_user_stoch_nogen.jl:265

The problem is that I'm sure that I've never declared P_tot_P_overestimate before. Do you think this can be fixed? Thank you for your time.

martinbiel commented 2 years ago

P_tot_P_overestimate is an @expression?

enri07 commented 2 years ago

Yes it is an @expression declared in the @stage2.

martinbiel commented 2 years ago

Could you ctrl-v the full @expression line?

enri07 commented 2 years ago

Of course:

@expression(model, P_tot_P_overestimate[u in user_set],
            max(0,
            sum(Float64[users_data[u].max_capacity[c] for c in users_data[u].asset_names if users_data[u].asset_type[c] == CONVERTER]) # Maximum capacity of the converters
            + maximum(sum(Float64[users_data[u].max_capacity[r]*Ren[u][r][t] for r = users_data[u].asset_names if users_data[u].asset_type[r] == RENEWABLE]) for t in time_set) # Maximum dispatch of renewable assets
            + sum(Float64[users_data[u].max_capacity[g]*users_data[u].max_dispatch[g] for g in users_data[u].asset_names if users_data[u].asset_type[g] == GENERATOR]) #Maximum dispatch of the fuel-fired generators
            - minimum(Load[u][t] for t in time_set)  # Minimum demand
            )
        )

user_set goes from 1 to n, where n is usually 4.

martinbiel commented 2 years ago

I am not sure why that does not work. Did a small test with an indexed @expression and it worked fine. If you e-mail me your model as a zip file I can have a look at it when I have time to search for the issue.

enri07 commented 2 years ago

Hi Martin,

I'm sorry for the delay but I was trying to make a minimal example of the above error.

I've a zip file ready to sent, where can I find your e-mail? Thank you very much for your help.

odow commented 2 years ago

You can attach files to this GitHub thread

image
enri07 commented 2 years ago

Thank you. For example, this zip file reproduce the error that I was trying to describe. Error.zip

martinbiel commented 2 years ago

Thank you. I had missed an edge case with @expressions that result in constant values. It should work on master now :slightly_smiling_face:

enri07 commented 2 years ago

Thank you, I am very happy if I can give an help to the improvement of the package, even if minimally! Now it works for the example attached before, but I would like to notice you that if we insert a dependence in the expression it gives the same problem. I have just modified the previous example with a simple dependence on the expression P_tot_P_overestimate to show that to you. Error1.zip

martinbiel commented 2 years ago

Yes, I should of course cover arrays of constant expressions as well. Should be fixed now!

enri07 commented 2 years ago

Perfect, everything works now!