Fix of objective_value when using TwoLevelTree.

Project.toml

@@ -1,7 +1,7 @@
 name = "TimeStruct"
 uuid = "f9ed5ce0-9f41-4eaa-96da-f38ab8df101c"
 authors = ["Lars Hellemo <Lars.Hellemo@sintef.no>, Truls.Flatberg <Truls.Flatberg@sintef.no>"]
-version = "0.9.10"
+version = "0.9.11"
 
 [deps]
 Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"

docs/src/manual/discount.md

@@ -27,3 +27,11 @@ df_avg = [discount(sp, ts, 0.05; type = "avg_year") for sp in sps]
 This approach results in a slightly higher discount factor.
 
 To help set up the objective function in a typical optimization problem, there is a utility function [`objective_weight`](@ref) that returns the weight to give a time period in the objective, considering both discount factor, probability, and possible multiplicity.
+
+!!! note "`objective_weight` for `TwoLevelTree`
+    You must be careful when you utilize the function [`objective_weight`](@ref) for a [`TimeStruct.StratNode`](@ref).
+    The function includes the multiplication with the branch probability.
+    This implies that you must divide the function by [`probability_branch`](@ref) if you have a variable that is indexed over strategic periods and already includes the branch probability.
+
+    You do not have to do this if your variable is not calculated from variables indexed over operational periods.
+    This is frequently the case for, *e.g.*,  investment variables with one example being [`EnergyModelsInvestments`](https://energymodelsx.github.io/EnergyModelsInvestments.jl/stable/manual/optimization-variables/).

src/discount.jl

@@ -135,12 +135,17 @@ end
     objective_weight(t::Union{TimePeriod,TimeStructurePeriod}, ts::TimeStructure, discount_rate; type = "start", timeunit_to_year = 1.0)
     objective_weight(t::Union{TimePeriod,TimeStructurePeriod}, disc::Discounter; type = "start")
 
-
 Calculates the overall objective weight for a time period `t` using a fixed `discount_rate`.
 The weight considers both discounting, the probability, and potential multiplicity of `t`.
 The function can be either called using a [`Discounter`](@ref) type or by specifying the
 parameters (time structure `ts`, `discount_rate` and potentially `timeunit_to_year`) directly.
 
+!!! note "Probabilities"
+    The objective weight includes the probability, both given through [`OperationalScenarios`](@ref)
+    or through [`TwoLevelTree`](@ref). This implies that if you use it on, *e.g.*, variables
+    indexed over strategic periods and calculate these variables through using the [`probability`](@ref)
+    function on operational periods, you have to take care not to have the calculation twice.
+
 There are two types of discounting available:
 
 1. Discounting to the start of the strategic period containing the time period:\n
@@ -195,5 +200,5 @@ function _objective_value(
     type,
     timeunit_to_year,
 )
-    return discount(t, ts, discount_rate; type, timeunit_to_year)
+    return probability_branch(t) * discount(t, ts, discount_rate; type, timeunit_to_year)
 end

src/profiles.jl

@@ -426,8 +426,9 @@ end
 -(a::StrategicProfile{T}) where {T} = StrategicProfile(-a.vals)
 -(a::ScenarioProfile{T}) where {T} = ScenarioProfile(-a.vals)
 -(a::RepresentativeProfile{T}) where {T} = RepresentativeProfile(-a.vals)
--(a::StrategicStochasticProfile{T}) where {T} =
-    StrategicStochasticProfile([.-v for v in a.vals])
+function -(a::StrategicStochasticProfile{T}) where {T}
+    return StrategicStochasticProfile([.-v for v in a.vals])
+end
 
 +(a::TimeProfile{T}) where {T} = a
 

src/structures.jl

@@ -85,14 +85,15 @@ multiple(t::TimePeriod) = 1.0
 """
     probability(t::TimePeriod)
 
-Returns the probability associated with the time period.
+Returns the probability associated with the time period. This include both the probability
+through [`OperationalScenarios`](@ref) and through [`TwoLevelTree`](@ref) structures.
 """
 probability(t::TimePeriod) = 1.0
 
 """
     probability_branch(t::Union{TimePeriod, TimeStructurePeriod})
 
-Returns the branch probability associated with the time period or time structure period
+Returns the branch probability associated with the time period or time structure period.
 """
 probability_branch(t::Union{TimePeriod,TimeStructurePeriod}) = 1.0
 

test/runtests.jl

@@ -1784,7 +1784,9 @@ end
     using Unitful
 
     uniform_years = TwoLevel(10, 1, SimpleTimes(1, 1))  # 10 years with duration of 1 year
+    uniform_years_tlt = TwoLevelTree(1, ones(Int, 9)*2, SimpleTimes(1, 1))
     disc = Discounter(0.04, 1, uniform_years)
+    disc_tlt = Discounter(0.04, 1, uniform_years_tlt)
 
     δ = 1 / 1.04
     for (i, t) in enumerate(uniform_years)
@@ -1797,46 +1799,101 @@ end
     ) ≈ 7.8017 atol = 1e-3
 
     @test sum(objective_weight(t, disc) for t in uniform_years) ≈ 8.435 atol = 1e-3
+    @test sum(objective_weight(t, disc_tlt) for t in uniform_years_tlt) ≈ 8.435 atol = 1e-3
 
     uniform_day = SimpleTimes(24, 1)
     periods = TwoLevel(10, 5 * 8760, uniform_day)
+    periods_tlt = TwoLevelTree(5 * 8760, ones(Int, 9)*2, uniform_day)
 
     @test sum(
         objective_weight(sp, periods, 0.04; timeunit_to_year = 1 / 8760, type = "start") for
         sp in strat_periods(periods)
     ) ≈ 4.825 atol = 1e-3
+    @test sum(
+        objective_weight(
+            sp,
+            periods_tlt,
+            0.04;
+            timeunit_to_year = 1 / 8760,
+            type = "start",
+        ) for sp in strat_periods(periods_tlt)
+    ) ≈ 4.825 atol = 1e-3
 
     @test sum(
         objective_weight(t, periods, 0.04; timeunit_to_year = 1 / 8760, type = "start") /
         (length(uniform_day) * probability(t) * multiple(t)) for t in periods
     ) ≈ 4.825 atol = 1e-3
+    @test sum(
+        objective_weight(
+            t,
+            periods_tlt,
+            0.04;
+            timeunit_to_year = 1 / 8760,
+            type = "start",
+        ) / (length(uniform_day) * probability(t) * multiple(t) / probability_branch(t)) for
+        t in periods_tlt
+    ) ≈ 4.825 atol = 1e-3
 
     @test sum(
         objective_weight(sp, periods, 0.04; timeunit_to_year = 1 / 8760, type = "avg_year")
         for sp in strat_periods(periods)
     ) ≈ 4.468 atol = 1e-3
+    @test sum(
+        objective_weight(
+            sp,
+            periods_tlt,
+            0.04;
+            timeunit_to_year = 1 / 8760,
+            type = "avg_year",
+        ) for sp in strat_periods(periods_tlt)
+    ) ≈ 4.468 atol = 1e-3
 
     @test sum(
         objective_weight(sp, periods, 0.04; timeunit_to_year = 1 / 8760, type = "avg") for
         sp in strat_periods(periods)
     ) ≈ 4.382 atol = 1e-3
+    @test sum(
+        objective_weight(sp, periods_tlt, 0.04; timeunit_to_year = 1 / 8760, type = "avg")
+        for sp in strat_periods(periods_tlt)
+    ) ≈ 4.382 atol = 1e-3
 
     oscs = OperationalScenarios(2, uniform_day)
     periods = TwoLevel(10, 5 * 8760, oscs)
+    periods_tlt = TwoLevelTree(5 * 8760, ones(Int, 9)*2, oscs)
 
     @test sum(
         objective_weight(t, periods, 0.04; timeunit_to_year = 1 / 8760, type = "start") /
         (length(oscs) * probability(t) * multiple(t)) for t in periods
     ) ≈ 4.825 atol = 1e-3
-    oscs
+    @test sum(
+        objective_weight(
+            t,
+            periods_tlt,
+            0.04;
+            timeunit_to_year = 1 / 8760,
+            type = "start",
+        ) / (length(oscs) * probability(t) * multiple(t) / probability_branch(t)) for
+        t in periods_tlt
+    ) ≈ 4.825 atol = 1e-3
 
     rps = RepresentativePeriods(2, 1, uniform_day)
     periods = TwoLevel(10, 5 * 8760, rps)
+    periods_tlt = TwoLevelTree(5 * 8760, ones(Int, 9)*2, rps)
 
     @test sum(
         objective_weight(t, periods, 0.04; timeunit_to_year = 1 / 8760, type = "start") /
         (length(rps) * probability(t) * multiple(t)) for t in periods
     ) ≈ 4.825 atol = 1e-3
+    @test sum(
+        objective_weight(
+            t,
+            periods_tlt,
+            0.04;
+            timeunit_to_year = 1 / 8760,
+            type = "start",
+        ) / (length(rps) * probability(t) * multiple(t) / probability_branch(t)) for
+        t in periods_tlt
+    ) ≈ 4.825 atol = 1e-3
 
     uniform_day = SimpleTimes(24, 1u"hr")
     periods_unit = TwoLevel(10, 365.125u"d", uniform_day)