From 4344d13bdfa2d7f9a5f1a1548e010ff5eebffe69 Mon Sep 17 00:00:00 2001 From: Vinicius Justen Pinto Date: Tue, 10 Dec 2024 16:05:25 -0300 Subject: [PATCH 1/3] fix cvar ub and add cvar tests --- Project.toml | 1 + src/cut_strategies/cuts_base.jl | 2 +- src/cut_strategies/multi_cut.jl | 18 +++++- src/policy.jl | 13 +++++ .../benders_training_iterations.jl | 1 - src/progress_logs/deterministic_equivalent.jl | 1 - src/training_strategies/benders_job_queue.jl | 8 +-- src/training_strategies/benders_serial.jl | 6 +- test/script_newsvendor_benders_job_queue.jl | 55 ++++++++++++++++--- test/test_newsvendor_benders.jl | 43 ++++++++++++--- 10 files changed, 118 insertions(+), 30 deletions(-) diff --git a/Project.toml b/Project.toml index 718393a..2a8c5a9 100644 --- a/Project.toml +++ b/Project.toml @@ -5,6 +5,7 @@ version = "0.2.0" [deps] EnumX = "4e289a0a-7415-4d19-859d-a7e5c4648b56" +HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b" JobQueueMPI = "32d208e1-246e-420c-b6ff-18b71b410923" JuMP = "4076af6c-e467-56ae-b986-b466b2749572" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" diff --git a/src/cut_strategies/cuts_base.jl b/src/cut_strategies/cuts_base.jl index 8bc9d8a..3eb56c2 100644 --- a/src/cut_strategies/cuts_base.jl +++ b/src/cut_strategies/cuts_base.jl @@ -37,7 +37,7 @@ function get_future_cost(model::JuMP.Model, policy_training_options)::Float64 if policy_training_options.cut_strategy == CutStrategy.SingleCut return LightBenders.get_single_cut_future_cost(model) else - return LightBenders.get_multi_cut_future_cost(model) + return LightBenders.get_multi_cut_future_cost(model, policy_training_options) end end diff --git a/src/cut_strategies/multi_cut.jl b/src/cut_strategies/multi_cut.jl index 2cac30f..e44c784 100644 --- a/src/cut_strategies/multi_cut.jl +++ b/src/cut_strategies/multi_cut.jl @@ -115,10 +115,22 @@ function add_all_cuts!(model::JuMP.Model, pool::CutPoolMultiCut, policy_training return nothing end -function get_multi_cut_future_cost(model::JuMP.Model)::Float64 +function get_multi_cut_future_cost(model::JuMP.Model, policy_training_options)::Float64 if !haskey(model, :epi_multi_cut) return 0.0 end - alphas = model[:epi_multi_cut]::Vector{JuMP.VariableRef} - return mean(JuMP.value.(alphas)) + alphas = JuMP.value.(model[:epi_multi_cut]) + if policy_training_options.risk_measure isa RiskNeutral + return mean(alphas) + elseif policy_training_options.risk_measure isa CVaR + discount_rate_multiplier = (1.0 - policy_training_options.discount_rate) + z_explicit_cvar = JuMP.value(model[:z_explicit_cvar]) + delta_explicit_cvar = JuMP.value.(model[:delta_explicit_cvar]) + fcf = z_explicit_cvar * discount_rate_multiplier * (policy_training_options.risk_measure.lambda) + for scen in 1:policy_training_options.num_scenarios + fcf += alphas[scen] * discount_rate_multiplier * (1 - policy_training_options.risk_measure.lambda) / policy_training_options.num_scenarios + fcf += delta_explicit_cvar[scen] * discount_rate_multiplier * (policy_training_options.risk_measure.lambda) / ((1 - policy_training_options.risk_measure.alpha) * policy_training_options.num_scenarios) + end + return fcf + end end diff --git a/src/policy.jl b/src/policy.jl index 9c9e850..5d2a532 100644 --- a/src/policy.jl +++ b/src/policy.jl @@ -12,3 +12,16 @@ end function upper_bound(policy::Policy) return last_upper_bound(policy.progress) end + +function second_stage_upper_bound_contribution(policy_training_options::PolicyTrainingOptions, objectives::Vector{Float64}) + num_scenarios = policy_training_options.num_scenarios + expected_value = sum(objectives[s] / num_scenarios for s in 1:num_scenarios) + if policy_training_options.risk_measure isa RiskNeutral + return expected_value + elseif policy_training_options.risk_measure isa CVaR + alpha = policy_training_options.risk_measure.alpha + lambda = policy_training_options.risk_measure.lambda + weights = build_cvar_weights(objectives, alpha, lambda) + return dot(weights, objectives) + end +end diff --git a/src/progress_logs/benders_training_iterations.jl b/src/progress_logs/benders_training_iterations.jl index 1de9bcd..4250dab 100644 --- a/src/progress_logs/benders_training_iterations.jl +++ b/src/progress_logs/benders_training_iterations.jl @@ -11,7 +11,6 @@ function BendersTrainingIterationsLog(policy_training_options::PolicyTrainingOpt println("Benders Training") println(" ") println("Training options:") - println("Number of stages: 2") println("Number of scenarios: ", policy_training_options.num_scenarios) println("Cut strategy: ", policy_training_options.cut_strategy) println("Risk measure: ", policy_training_options.risk_measure) diff --git a/src/progress_logs/deterministic_equivalent.jl b/src/progress_logs/deterministic_equivalent.jl index 4d065e3..f206a4a 100644 --- a/src/progress_logs/deterministic_equivalent.jl +++ b/src/progress_logs/deterministic_equivalent.jl @@ -2,7 +2,6 @@ function DeterministicEquivalentLog(num_scenarios::Int) println(" ") println("Deterministic Equivalent") println(" ") - println("Number of stages: 2") println("Number of scenarios: ", num_scenarios) return nothing end diff --git a/src/training_strategies/benders_job_queue.jl b/src/training_strategies/benders_job_queue.jl index d7ffb92..00d1a80 100644 --- a/src/training_strategies/benders_job_queue.jl +++ b/src/training_strategies/benders_job_queue.jl @@ -16,7 +16,6 @@ mutable struct SecondStageMessage end mutable struct SecondStageAnswer - UB::Float64 coefs rhs obj @@ -93,7 +92,6 @@ function job_queue_benders_train(; if !isnothing(job_answer) message = JQM.get_message(job_answer) if message isa SecondStageAnswer - progress.UB[progress.current_iteration] += message.UB store_cut!( local_pools, message.coefs, @@ -111,6 +109,9 @@ function job_queue_benders_train(; # Cuts here can be following the single cut strategy or # the multi cut strategy store_cut!(pool, local_pools, state, policy_training_options, t) + progress.UB[progress.current_iteration] += second_stage_upper_bound_contribution( + policy_training_options, local_pools.obj + ) report_current_bounds(progress) convergence_result = convergence_test(progress, policy_training_options.stopping_rule) @@ -224,10 +225,7 @@ function worker_second_stage( optimize_with_retry(second_stage_model) treat_termination_status(second_stage_model, policy_training_options, t, scenario, iteration) coefs, rhs, obj = get_cut(second_stage_model, state) - future_cost = get_future_cost(second_stage_model, policy_training_options) - UB = (JuMP.objective_value(second_stage_model) - future_cost) / policy_training_options.num_scenarios return SecondStageAnswer( - UB, coefs, rhs, obj, diff --git a/src/training_strategies/benders_serial.jl b/src/training_strategies/benders_serial.jl index 1a7b934..2192ac2 100644 --- a/src/training_strategies/benders_serial.jl +++ b/src/training_strategies/benders_serial.jl @@ -40,10 +40,10 @@ function serial_benders_train(; coefs, rhs, obj = get_cut(second_stage_model, state) # Store the opening cut in a temporary cut pool store_cut!(local_pools, coefs, state, rhs, obj) - future_cost = get_future_cost(second_stage_model, policy_training_options) - progress.UB[progress.current_iteration] += - (JuMP.objective_value(second_stage_model) - future_cost) / policy_training_options.num_scenarios end + progress.UB[progress.current_iteration] += second_stage_upper_bound_contribution( + policy_training_options, local_pools.obj + ) # Store the (stage, scenario) cut(s) in a persitent pool. # Cuts here can be following the single cut strategy or # the multi cut strategy diff --git a/test/script_newsvendor_benders_job_queue.jl b/test/script_newsvendor_benders_job_queue.jl index 0ed5574..03d664b 100644 --- a/test/script_newsvendor_benders_job_queue.jl +++ b/test/script_newsvendor_benders_job_queue.jl @@ -49,7 +49,10 @@ function second_stage_modifier(sp, inputs, s) return nothing end -function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) +function newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.MultiCut, + risk_measure = LightBenders.RiskNeutral(), + ) inputs = Inputs(5, 10, 1, 100, [10, 20, 30]) num_scenarios = length(inputs.demand) @@ -62,6 +65,7 @@ function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) min_iterations = 2, ), cut_strategy = cut_strategy, + risk_measure = risk_measure, ) policy = LightBenders.train(; @@ -77,9 +81,6 @@ function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) return nothing end - @test LightBenders.lower_bound(policy) ≈ -70 - @test LightBenders.upper_bound(policy) ≈ -70 - results = LightBenders.simulate(; state_variables_builder, first_stage_builder, @@ -93,15 +94,55 @@ function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) ), ) - @test results["objective", 0] ≈ -70 atol = 1e-2 + return policy, results end function test_newsvendor_benders() @testset "[Job Queue] Benders Newsvendor single cut" begin - newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.SingleCut) + v = newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.SingleCut + ) + if v !== nothing + policy, results = v + @test LightBenders.lower_bound(policy) ≈ -70 + @test LightBenders.upper_bound(policy) ≈ -70 + @test results["objective", 0] ≈ -70 atol = 1e-2 + end end @testset "[Job Queue] Benders Newsvendor multi cut" begin - newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) + v = newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.MultiCut + ) + if v !== nothing + policy, results = v + @test LightBenders.lower_bound(policy) ≈ -70 + @test LightBenders.upper_bound(policy) ≈ -70 + @test results["objective", 0] ≈ -70 atol = 1e-2 + end + end + @testset "[Job Queue] Benders Newsvendor single cut CVaR" begin + v = newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.SingleCut, + risk_measure = LightBenders.CVaR(alpha = 0.9, lambda = 0.5) + ) + if v !== nothing + policy, results = v + @test LightBenders.lower_bound(policy) ≈ -50 + @test LightBenders.upper_bound(policy) ≈ -50 + @test results["objective", 0] ≈ -50 atol = 1e-2 + end + end + @testset "[Job Queue] Benders Newsvendor multi cut CVaR" begin + v = newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.MultiCut, + risk_measure = LightBenders.CVaR(alpha = 0.9, lambda = 0.5) + ) + if v !== nothing + policy, results = v + @test LightBenders.lower_bound(policy) ≈ -50 + @test LightBenders.upper_bound(policy) ≈ -50 + @test results["objective", 0] ≈ -50 atol = 1e-2 + end end end diff --git a/test/test_newsvendor_benders.jl b/test/test_newsvendor_benders.jl index e9374ec..fefd004 100644 --- a/test/test_newsvendor_benders.jl +++ b/test/test_newsvendor_benders.jl @@ -49,7 +49,10 @@ function second_stage_modifier(sp, inputs, s) return nothing end -function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) +function newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.MultiCut, + risk_measure = LightBenders.RiskNeutral(), +) inputs = Inputs(5, 10, 1, 100, [10, 20, 30]) num_scenarios = length(inputs.demand) @@ -62,6 +65,7 @@ function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) min_iterations = 2, ), cut_strategy = cut_strategy, + risk_measure = risk_measure, ) policy = LightBenders.train(; @@ -73,9 +77,6 @@ function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) policy_training_options, ) - @test LightBenders.lower_bound(policy) ≈ -70 - @test LightBenders.upper_bound(policy) ≈ -70 - results = LightBenders.simulate(; state_variables_builder, first_stage_builder, @@ -89,7 +90,7 @@ function newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) ), ) - @test results["objective", 0] ≈ -70 atol = 1e-2 + return policy, results end function newsvendor_deterministic() @@ -111,11 +112,35 @@ function newsvendor_deterministic() end function test_newsvendor_benders() - @testset "Benders Newsvendor single cut" begin - newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.SingleCut) + @testset "Benders Newsvendor single cut risk neutral" begin + policy, results = newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.SingleCut) + @test LightBenders.lower_bound(policy) ≈ -70 + @test LightBenders.upper_bound(policy) ≈ -70 + @test results["objective", 0] ≈ -70 atol = 1e-2 + end + @testset "Benders Newsvendor multi cut risk neutral" begin + policy, results = newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) + @test LightBenders.lower_bound(policy) ≈ -70 + @test LightBenders.upper_bound(policy) ≈ -70 + @test results["objective", 0] ≈ -70 atol = 1e-2 + end + @testset "Benders Newsvendor single cut CVaR" begin + policy, results = newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.SingleCut, + risk_measure = LightBenders.CVaR(alpha = 0.9, lambda = 0.5), + ) + @test LightBenders.lower_bound(policy) ≈ -50 + @test LightBenders.upper_bound(policy) ≈ -50 + @test results["objective", 0] ≈ -50 atol = 1e-2 end - @testset "Benders Newsvendor multi cut" begin - newsvendor_benders(; cut_strategy = LightBenders.CutStrategy.MultiCut) + @testset "Benders Newsvendor multi cut CVaR" begin + policy, results = newsvendor_benders(; + cut_strategy = LightBenders.CutStrategy.MultiCut, + risk_measure = LightBenders.CVaR(alpha = 0.9, lambda = 0.5), + ) + @test LightBenders.lower_bound(policy) ≈ -50 + @test LightBenders.upper_bound(policy) ≈ -50 + @test results["objective", 0] ≈ -50 atol = 1e-2 end @testset "Deterministic equivalent Newsvendor" begin newsvendor_deterministic() From 800247fba64c21679e737e89438193b91b088778 Mon Sep 17 00:00:00 2001 From: guilhermebodin Date: Wed, 11 Dec 2024 11:59:58 -0300 Subject: [PATCH 2/3] remove HiGHS as a dependency --- Project.toml | 1 - 1 file changed, 1 deletion(-) diff --git a/Project.toml b/Project.toml index 2a8c5a9..718393a 100644 --- a/Project.toml +++ b/Project.toml @@ -5,7 +5,6 @@ version = "0.2.0" [deps] EnumX = "4e289a0a-7415-4d19-859d-a7e5c4648b56" -HiGHS = "87dc4568-4c63-4d18-b0c0-bb2238e4078b" JobQueueMPI = "32d208e1-246e-420c-b6ff-18b71b410923" JuMP = "4076af6c-e467-56ae-b986-b466b2749572" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" From c1a7281cb68c40b706d24bfb345566febcfc1bc9 Mon Sep 17 00:00:00 2001 From: guilhermebodin Date: Wed, 11 Dec 2024 12:23:08 -0300 Subject: [PATCH 3/3] update logs --- .../benders_training_iterations.jl | 23 +++++++++++-------- src/progress_logs/deterministic_equivalent.jl | 8 +++---- src/simulation_strategies/benders_serial.jl | 4 ++-- src/training_strategies/benders_job_queue.jl | 3 +-- src/training_strategies/benders_serial.jl | 3 +-- 5 files changed, 21 insertions(+), 20 deletions(-) diff --git a/src/progress_logs/benders_training_iterations.jl b/src/progress_logs/benders_training_iterations.jl index 4250dab..fbb04e8 100644 --- a/src/progress_logs/benders_training_iterations.jl +++ b/src/progress_logs/benders_training_iterations.jl @@ -7,15 +7,14 @@ Base.@kwdef mutable struct BendersTrainingIterationsLog <: AbstractProgressLog end function BendersTrainingIterationsLog(policy_training_options::PolicyTrainingOptions) - println(" ") - println("Benders Training") - println(" ") - println("Training options:") - println("Number of scenarios: ", policy_training_options.num_scenarios) - println("Cut strategy: ", policy_training_options.cut_strategy) - println("Risk measure: ", policy_training_options.risk_measure) - println("Stopping rule: ", policy_training_options.stopping_rule) - # TODO add more prints + @info(" ") + @info("Benders Training") + @info(" ") + @info("Training options:") + @info("Number of scenarios: ", policy_training_options.num_scenarios) + @info("Cut strategy: ", policy_training_options.cut_strategy) + @info("Risk measure: ", policy_training_options.risk_measure) + @info("Stopping rule: ", policy_training_options.stopping_rule) progress_table = ProgressTable( header = ["Iteration", "Lower bound", "Upper bound", "Gap", "Time [s]"], @@ -70,7 +69,11 @@ function report_current_bounds(progress::BendersTrainingIterationsLog) return nothing end -function finish_training!(progress::BendersTrainingIterationsLog) +function finish_training!( + progress::BendersTrainingIterationsLog, + convergence_result::ConvergenceResult, +) finalize(progress.progress_table) + @info(results_message(convergence_result)) return nothing end diff --git a/src/progress_logs/deterministic_equivalent.jl b/src/progress_logs/deterministic_equivalent.jl index f206a4a..ddc2778 100644 --- a/src/progress_logs/deterministic_equivalent.jl +++ b/src/progress_logs/deterministic_equivalent.jl @@ -1,7 +1,7 @@ function DeterministicEquivalentLog(num_scenarios::Int) - println(" ") - println("Deterministic Equivalent") - println(" ") - println("Number of scenarios: ", num_scenarios) + @info(" ") + @info("Deterministic Equivalent") + @info(" ") + @info("Number of scenarios: ", num_scenarios) return nothing end diff --git a/src/simulation_strategies/benders_serial.jl b/src/simulation_strategies/benders_serial.jl index 1c47647..311da74 100644 --- a/src/simulation_strategies/benders_serial.jl +++ b/src/simulation_strategies/benders_serial.jl @@ -14,7 +14,7 @@ function serial_benders_simulate(; results = Dict{Tuple{String, Int}, Any}() # (variable_name, scenario) => value # first stage - println("Simulating first stage...") + @info("Simulating first stage...") state_variables_model = state_variables_builder(inputs) model = first_stage_builder(state_variables_model, inputs) @@ -31,7 +31,7 @@ function serial_benders_simulate(; end # second stage - println("Simulating second stage...") + @info("Simulating second stage...") state = if simulation_options.state_handling == SimulationStateHandling.StatesRecalculatedInSimulation get_state(model) diff --git a/src/training_strategies/benders_job_queue.jl b/src/training_strategies/benders_job_queue.jl index 00d1a80..22173b8 100644 --- a/src/training_strategies/benders_job_queue.jl +++ b/src/training_strategies/benders_job_queue.jl @@ -116,8 +116,7 @@ function job_queue_benders_train(; convergence_result = convergence_test(progress, policy_training_options.stopping_rule) if has_converged(convergence_result) - println(results_message(convergence_result)) - finish_training!(progress) + finish_training!(progress, convergence_result) JQM.send_termination_message() break end diff --git a/src/training_strategies/benders_serial.jl b/src/training_strategies/benders_serial.jl index 2192ac2..b4b5739 100644 --- a/src/training_strategies/benders_serial.jl +++ b/src/training_strategies/benders_serial.jl @@ -54,8 +54,7 @@ function serial_benders_train(; convergence_result = convergence_test(progress, policy_training_options.stopping_rule) if has_converged(convergence_result) - finish_training!(progress) - println(results_message(convergence_result)) + finish_training!(progress, convergence_result) break end end