diff --git a/docs/src/guides/add_a_risk_measure.md b/docs/src/guides/add_a_risk_measure.md index 3463b52d8..2631d768d 100644 --- a/docs/src/guides/add_a_risk_measure.md +++ b/docs/src/guides/add_a_risk_measure.md @@ -304,3 +304,32 @@ round.(risk_adjusted_probability, digits = 1) 0.8 0.0 ``` + +### Entropic + +```@docs +SDDP.Entropic +``` + +```jldoctest intermediate_risk +risk_measure = SDDP.Entropic(0.1) + +SDDP.adjust_probability( + risk_measure, + risk_adjusted_probability, + nominal_probability, + noise_supports, + cost_realizations, + is_minimization +) + +round.(risk_adjusted_probability, digits = 4) + +# output + +4-element Array{Float64,1}: + 0.11 + 0.1991 + 0.3648 + 0.326 +``` diff --git a/src/plugins/risk_measures.jl b/src/plugins/risk_measures.jl index 6feef8fbf..089e6aa94 100644 --- a/src/plugins/risk_measures.jl +++ b/src/plugins/risk_measures.jl @@ -541,3 +541,46 @@ function adjust_probability( copyto!(risk_adjusted_probability, JuMP.value.(p)) return 0.0 end + +# ================================= Entropic ============================== # + +""" + Entropic(γ::Float64) + +The entropic risk measure as described by Dowson, Morton, and Pagnoncelli +(2020). Multistage stochastic programs with the entropic risk measure. +http://www.optimization-online.org/DB_HTML/2020/08/7984.html +""" +mutable struct Entropic <: AbstractRiskMeasure + γ::Float64 +end + +function Base.show(io::IO, ent::Entropic) + return print(io, "Entropic risk measure with γ = $(ent.γ)") +end + +function adjust_probability( + measure::Entropic, + Q::Vector{Float64}, + p::Vector{Float64}, + ::Vector, + X::Vector{Float64}, + is_min::Bool, +) + if measure.γ == 0.0 # Special case of entropic: if γ = 0, F[X] = E[X]. + Q .= p + return 0.0 + end + # Handle maximization problems by negating γ. Usually we would negate X, but + # with convex RM's, we also have to negate the α(q) term, so it's easier to + # just negate γ. + γ = (is_min ? 1.0 : -1.0) * measure.γ + # Use `BigFloat` to avoid overflow that occurs when calculating `exp(x)`. + y = p .* exp.(big.(γ .* X)) + Q .= y / sum(y) + α = sum( + # To avoid numerical issues with the log, skip elements that are `≈ 0`. + qi * log(qi / pi) for (pi, qi) in zip(p, Q) if pi > 1e-14 && qi > 1e-14 + ) + return -α / γ +end diff --git a/test/plugins/risk_measures.jl b/test/plugins/risk_measures.jl index 602fe8de9..e9480767b 100644 --- a/test/plugins/risk_measures.jl +++ b/test/plugins/risk_measures.jl @@ -211,7 +211,7 @@ end ) @test risk_adjusted_probability ≈ [0.1, 0.2, 0.3, 0.2, 0.2] atol = 1e-6 end - + @testset "Non-uniform Worst Case" begin risk_adjusted_probability = Vector{Float64}(undef, 5) SDDP.adjust_probability( @@ -224,7 +224,7 @@ end ) @test risk_adjusted_probability ≈ [0.0, 1.0, 0.0, 0.0, 0.0] atol = 1e-6 end - + @testset "Non-uniform Worst Case Max" begin risk_adjusted_probability = Vector{Float64}(undef, 5) SDDP.adjust_probability( @@ -237,7 +237,7 @@ end ) @test risk_adjusted_probability ≈ [0.0, 0.0, 0.0, 0.0, 1.0] atol = 1e-6 end - + @testset "Non-uniform" begin risk_adjusted_probability = Vector{Float64}(undef, 5) SDDP.adjust_probability( @@ -410,3 +410,40 @@ end @test risk_adjusted_probability ≈ [0.0, 1 / 6, 5 / 6, 0.0] end end + +@testset "Entropic" begin + @test sprint(show, SDDP.Entropic(0.1)) == "Entropic risk measure with γ = 0.1" + @testset "Minimization" begin + # Test that increasing values of θ lead to larger values for F[X]. + X = [1.0, 2.0, 3.0] + P = [0.5, 0.5, 0.0] + Q = [NaN, NaN, NaN] + last, last_q2 = -Inf, 0.0 + for i = -4:4 + θ = 10.0^i + α = SDDP.adjust_probability(SDDP.Entropic(θ), Q, P, [], X, true) + current = Q' * X + α + @test current >= last + @test Q[2] >= last_q2 + last, last_q2 = current, Q[2] + end + end + @testset "Maximization" begin + # Test that increasing values of θ lead to smaller values for F[X]. + X = [1.0, 2.0, 3.0] + P = [0.5, 0.5, 0.0] + Q = [NaN, NaN, NaN] + last, last_q1 = Inf, 0.0 + for i = -4:4 + θ = 10.0^i + α = SDDP.adjust_probability(SDDP.Entropic(θ), Q, P, [], X, false) + current = Q' * X + α + if Q[1] < 1 - 1e-6 + # If Q[1] ≈ 1, this test can fail due to numerical error. + @test current <= last + end + @test Q[1] >= last_q1 + last, last_q1 = current, Q[1] + end + end +end