Add min_slope and max_slope params, update monopoly criteria

database/migrations/38/down.sql

@@ -0,0 +1,7 @@
+PRAGMA user_version = 37;
+
+ALTER TABLE Configuration DROP COLUMN supply_function_equilibrium_min_slope;
+ALTER TABLE Configuration DROP COLUMN supply_function_equilibrium_max_slope;
+
+ALTER TABLE Configuration ADD COLUMN supply_function_equilibrium_extra_bid_quantity REAL DEFAULT 1.0;
+ALTER TABLE Configuration ADD COLUMN supply_function_equilibrium_tolerance REAL DEFAULT 0.000001;

database/migrations/38/up.sql

@@ -0,0 +1,7 @@
+PRAGMA user_version = 38;
+
+ALTER TABLE Configuration DROP COLUMN supply_function_equilibrium_extra_bid_quantity;
+ALTER TABLE Configuration DROP COLUMN supply_function_equilibrium_tolerance;
+
+ALTER TABLE Configuration ADD COLUMN supply_function_equilibrium_min_slope REAL NOT NULL DEFAULT 0.02;
+ALTER TABLE Configuration ADD COLUMN supply_function_equilibrium_max_slope REAL NOT NULL DEFAULT 50.0;

docs/src/supply_function_equilibrium.md

@@ -26,15 +26,15 @@ IARA.update_configuration!(db;
 
 ### SFE Parameters
 
-- **`supply_function_equilibrium_extra_bid_quantity`** (Float64, default: 1.0)
-  Determines the quantity value for the artificial bid added during SFE preprocessing.
+- **`supply_function_equilibrium_min_slope`** (Float64, default: 0.02)
+  Minimum slope allowed for the bid curves.
+
+- **`supply_function_equilibrium_max_slope`** (Float64, default: 50.0)
+  Maximum slope allowed for the bid curves. Must be greater than `supply_function_equilibrium_min_slope`.
 
 - **`supply_function_equilibrium_max_iterations`** (Int, default: 20)
   Maximum iterations for equilibrium computation within each period/scenario.
 
-- **`supply_function_equilibrium_tolerance`** (Float64, default: 0.000001)
-  Minimum slope tolerance. Curves with slopes below this value trigger errors.
-
 - **`supply_function_equilibrium_max_cost_multiplier`** (Float64, default: 2.0)
   Maximum price cap as multiplier of deficit cost.
 
@@ -43,9 +43,9 @@ IARA.update_configuration!(db;
 IARA.update_configuration!(db;
     bid_processing = IARA.Configurations_BidProcessing.ITERATED_BIDS_FROM_SUPPLY_FUNCTION_EQUILIBRIUM,
     reference_curve_number_of_segments = 10,
-    supply_function_equilibrium_extra_bid_quantity = 1.0,
+    supply_function_equilibrium_min_slope = 0.02,
+    supply_function_equilibrium_max_slope = 50.0,
     supply_function_equilibrium_max_iterations = 20,
-    supply_function_equilibrium_tolerance = 0.000001,
     supply_function_equilibrium_max_cost_multiplier = 2.0,
     demand_deficit_cost = 3000.0,
 )

src/collections/configurations.jl

@@ -101,8 +101,8 @@ Configurations for the problem.
     bid_price_limit_high_reference::Float64 = 0.0
     reference_curve_number_of_segments::Int = 0
     reference_curve_final_segment_price_markup::Float64 = 0.0
-    supply_function_equilibrium_extra_bid_quantity::Float64 = 0.0
-    supply_function_equilibrium_tolerance::Float64 = 0.0
+    supply_function_equilibrium_min_slope::Float64 = 0.0
+    supply_function_equilibrium_max_slope::Float64 = 0.0
     supply_function_equilibrium_max_iterations::Int = 0
     supply_function_equilibrium_max_cost_multiplier::Float64 = 0.0
     supply_function_equilibrium_price_taker_weight::Float64 = 0.0
@@ -330,10 +330,10 @@ function initialize!(configurations::Configurations, inputs::AbstractInputs)
         PSRI.get_parms(inputs.db, "Configuration", "reference_curve_number_of_segments")[1]
     configurations.reference_curve_final_segment_price_markup =
         PSRI.get_parms(inputs.db, "Configuration", "reference_curve_final_segment_price_markup")[1]
-    configurations.supply_function_equilibrium_extra_bid_quantity =
-        PSRI.get_parms(inputs.db, "Configuration", "supply_function_equilibrium_extra_bid_quantity")[1]
-    configurations.supply_function_equilibrium_tolerance =
-        PSRI.get_parms(inputs.db, "Configuration", "supply_function_equilibrium_tolerance")[1]
+    configurations.supply_function_equilibrium_min_slope =
+        PSRI.get_parms(inputs.db, "Configuration", "supply_function_equilibrium_min_slope")[1]
+    configurations.supply_function_equilibrium_max_slope =
+        PSRI.get_parms(inputs.db, "Configuration", "supply_function_equilibrium_max_slope")[1]
     configurations.supply_function_equilibrium_max_iterations =
         PSRI.get_parms(inputs.db, "Configuration", "supply_function_equilibrium_max_iterations")[1]
     configurations.supply_function_equilibrium_max_cost_multiplier =
@@ -539,6 +539,20 @@ function validate(configurations::Configurations)
             num_errors += 1
         end
     end
+    if configurations.supply_function_equilibrium_min_slope < DEFAULT_TOLERANCE
+        @error(
+            "supply_function_equilibrium_min_slope must be at least $(DEFAULT_TOLERANCE). " *
+            "Current value: $(configurations.supply_function_equilibrium_min_slope)."
+        )
+        num_errors += 1
+    end
+    if configurations.supply_function_equilibrium_max_slope <= configurations.supply_function_equilibrium_min_slope
+        @error(
+            "supply_function_equilibrium_max_slope must be greater than supply_function_equilibrium_min_slope. " *
+            "Current values: max_slope = $(configurations.supply_function_equilibrium_max_slope), min_slope = $(configurations.supply_function_equilibrium_min_slope)."
+        )
+        num_errors += 1
+    end
     if configurations.cvar_alpha <= 0.0 || configurations.cvar_alpha > 1.0
         @error("cvar_alpha must be in (0.0, 1.0].")
         num_errors += 1
@@ -1654,20 +1668,20 @@ reference_curve_final_segment_price_markup(inputs::AbstractInputs) =
     inputs.collections.configurations.reference_curve_final_segment_price_markup
 
 """
-    supply_function_equilibrium_extra_bid_quantity(inputs::AbstractInputs)
+    supply_function_equilibrium_min_slope(inputs::AbstractInputs)
 
-Return the extra bid quantity for the Supply Function Equilibrium.
+Return the minimum slope for the bid curves in the Supply Function Equilibrium.
 """
-supply_function_equilibrium_extra_bid_quantity(inputs::AbstractInputs) =
-    inputs.collections.configurations.supply_function_equilibrium_extra_bid_quantity
+supply_function_equilibrium_min_slope(inputs::AbstractInputs) =
+    inputs.collections.configurations.supply_function_equilibrium_min_slope
 
 """
-    supply_function_equilibrium_tolerance(inputs)
+    supply_function_equilibrium_max_slope(inputs::AbstractInputs)
 
-Return the tolerance for the bid slopes in the Supply Function Equilibrium.
+Return the maximum slope for the bid curves in the Supply Function Equilibrium.
 """
-supply_function_equilibrium_tolerance(inputs::AbstractInputs) =
-    inputs.collections.configurations.supply_function_equilibrium_tolerance
+supply_function_equilibrium_max_slope(inputs::AbstractInputs) =
+    inputs.collections.configurations.supply_function_equilibrium_max_slope
 
 """
     supply_function_equilibrium_max_iterations(inputs)

src/inputs.jl

@@ -144,8 +144,8 @@ Optional arguments:
   - `language::String`
   - `train_mincost_time_limit_sec::Int64`
   - `bid_price_limit_low_reference::Float64`
-  - `supply_function_equilibrium_extra_bid_quantity::Float64`
-  - `supply_function_equilibrium_tolerance::Float64`
+  - `supply_function_equilibrium_min_slope::Float64`
+  - `supply_function_equilibrium_max_slope::Float64`
   - `supply_function_equilibrium_max_iterations::Int64`
   - `supply_function_equilibrium_max_cost_multiplier::Float64`
   - `subperiod_duration_in_hours::Vector{Float64}`: Subperiod duration in hours of the configuration

src/supply_function_equilibrium.jl

@@ -341,12 +341,11 @@ function treat_reference_curve_data(
             price[vr_index, ao, :],
         )
         treated_quantity_bids[i] = quantity_points_from_segments(treated_quantity_bids[i])
-        treated_quantity_bids[i], treated_price_bids[i] = reverse_bid_order_and_add_point(
+        treated_quantity_bids[i], treated_price_bids[i], treated_slopes[i] = reverse_bid_order_and_add_points(
             inputs,
             treated_quantity_bids[i],
             treated_price_bids[i],
         )
-        treated_slopes[i] = diff(treated_price_bids[i]) ./ diff(treated_quantity_bids[i])
     end
 
     validate_nash_input_data(
@@ -395,12 +394,11 @@ function treat_bidding_group_data(
             aggregated_price[i, :],
         )
         treated_quantity_bids[i] = quantity_points_from_segments(treated_quantity_bids[i])
-        treated_quantity_bids[i], treated_price_bids[i] = reverse_bid_order_and_add_point(
+        treated_quantity_bids[i], treated_price_bids[i], treated_slopes[i] = reverse_bid_order_and_add_points(
             inputs,
             treated_quantity_bids[i],
             treated_price_bids[i],
         )
-        treated_slopes[i] = diff(treated_price_bids[i]) ./ diff(treated_quantity_bids[i])
     end
 
     validate_nash_input_data(
@@ -493,15 +491,64 @@ function quantity_points_from_segments(
     return new_quantity
 end
 
-function reverse_bid_order_and_add_point(
+function reverse_bid_order_and_add_points(
     inputs::AbstractInputs,
     quantity::Vector{Float64},
     price::Vector{Float64},
 )
-    new_quantity = vcat((quantity[end] + supply_function_equilibrium_extra_bid_quantity(inputs)), reverse(quantity))
-    new_price = vcat(demand_deficit_cost(inputs), reverse(price))
+    new_quantity = Float64[]
+    new_price = Float64[]
+    new_slope = Float64[]
 
-    return new_quantity, new_price
+    @assert length(quantity) == length(price)
+
+    reference_quantity = vcat(0, quantity)
+    reference_price = vcat(price, demand_deficit_cost(inputs))
+    number_of_points = length(reference_quantity)
+
+    min_slope = supply_function_equilibrium_min_slope(inputs)
+    max_slope = supply_function_equilibrium_max_slope(inputs)
+
+    #(q0, p0)
+    push!(new_quantity, reference_quantity[1])
+    push!(new_price, reference_price[1])
+
+    for i in 1:number_of_points-1
+        q0 = reference_quantity[i]
+        q1 = reference_quantity[i+1]
+        p0 = reference_price[i]
+        p1 = reference_price[i+1]
+
+        slope = (p1 - p0) / (q1 - q0)
+
+        # If the slope is outside the [min_slope, max_slope] range, no new point is added
+        if slope <= min_slope || slope >= max_slope
+            push!(new_quantity, q1)
+            push!(new_price, p1)
+            push!(new_slope, slope)
+            continue
+        end
+
+        # (q_new, p_new) is the intersection of a line passing through (q0, p0) with slope eta_min
+        # and a line passing through (p1, q1) with slope eta_max
+        q_new = q0 + (max_slope * (q1 - q0) - (p1 - p0)) / (max_slope - min_slope)
+        p_new = p0 + min_slope * (q_new - q0)
+
+        # new point
+        push!(new_quantity, q_new)
+        push!(new_price, p_new)
+        push!(new_slope, min_slope)
+        # (q1, p1)
+        push!(new_quantity, q1)
+        push!(new_price, p1)
+        push!(new_slope, max_slope)
+    end
+
+    reverse!(new_quantity)
+    reverse!(new_price)
+    reverse!(new_slope)
+
+    return new_quantity, new_price, new_slope
 end
 
 function validate_nash_input_data(
@@ -514,7 +561,7 @@ function validate_nash_input_data(
     @assert length(agent_labels) == length(slope)
 
     for i in eachindex(agent_labels)
-        @assert all(slope[i] .> supply_function_equilibrium_tolerance(inputs)) "Reference bid curve for $(agent_labels[i]) has a segment with slope below the tolerance: $(slope[i])"
+        @assert all(slope[i] .> DEFAULT_TOLERANCE^2) "Reference bid curve for $(agent_labels[i]) has a segment with slope below the tolerance: $(slope[i])"
         @assert all(price[i] .<= supply_function_equilibrium_max_cost_multiplier(inputs) * demand_deficit_cost(inputs)) "Reference bid curve for $(agent_labels[i]) has a price point above the demand deficit cost: $(price[i])"
     end
 
@@ -527,8 +574,8 @@ end
 Convert the equilibrium curve of a single agent into incremental quantity and price bids, in the IARA convention.
 
 The equilibrium curve is a list of cumulative quantity points in descending price order, whose first point is the
-synthetic point that `reverse_bid_order_and_add_point` added at the demand deficit cost. This function inverts that
-representation: it drops the synthetic point, restores the ascending price order, and differentiates the cumulative
+synthetic point that `reverse_bid_order_and_add_points` added at the demand deficit cost. This function inverts that
+representation: it drops the last synthetic point, restores the ascending price order, and differentiates the cumulative
 quantities back into incremental bids (inverting `quantity_points_from_segments`). The result is trimmed or zero-padded
 to `number_of_bid_segments`.
 """
@@ -796,7 +843,7 @@ function run_supply_function_equilibrium_iteration(
         minimum_quantities = [minimum(original_quantity[i]) for i in 1:number_of_asset_owners]
         if maximum(
             [new_quantity[i][segment] for i in 1:number_of_asset_owners] - minimum_quantities,
-        ) == 0
+        ) <= DEFAULT_TOLERANCE
             break
         end
 
@@ -889,8 +936,11 @@ function update_slope(
 
     agent_indexes = findall(isfinite, original_slope_in_segment)
 
-    if length(agent_indexes) < 3
-        new_slope[agent_indexes] .= supply_function_equilibrium_tolerance(inputs)
+    if sum(agent_price_type_weight[agent_indexes]) < 3
+        for agent_idx in agent_indexes
+            new_slope[agent_idx] =
+                min(supply_function_equilibrium_min_slope(inputs), original_slope_in_segment[agent_idx])
+        end
     end
 
     return new_slope
@@ -992,10 +1042,12 @@ function get_current_segment(
     minimum_quantities = [minimum(original_quantity[i]) for i in 1:number_of_asset_owners]
 
     for i in 1:number_of_asset_owners
-        if current_quantity_in_segment[i] > minimum_quantities[i]
-            idx = findfirst(reference_quantity[i] .< current_quantity_in_segment[i])
+        if current_quantity_in_segment[i] > minimum_quantities[i] + DEFAULT_TOLERANCE
+            idx = findfirst(
+                reference_quantity[i] .< current_quantity_in_segment[i] - DEFAULT_TOLERANCE,
+            )
             if isnothing(idx)
-                segments[i] = number_of_segments_for_vr_in_nash_equilibrium(inputs, number_of_asset_owners)
+                segments[i] = length(reference_quantity[i]) - 1
             else
                 segments[i] = idx - 1
             end
@@ -1036,7 +1088,7 @@ function test_inversion(
         (original_quantity_in_segment .- quantity_in_segment[agent_indexes]) .* original_slope_in_segment
     price_delta = price_in_segment .- reference_price
 
-    if any(price_delta .< 0)
+    if any(price_delta .< -DEFAULT_TOLERANCE)
         @warn("Curve inversion")
     end
 

test/case_07/supply_function_equilibrium/expected_outputs/deficit_ex_ante_commercial.csv

@@ -1,25 +1,25 @@
 period,scenario,subperiod,dem_1
 1,1,1,300
-1,1,2,150
-1,1,3,278.378
+1,1,2,289.266
+1,1,3,300
 1,1,4,300
 1,2,1,400
-1,2,2,250
-1,2,3,378.378
+1,2,2,389.266
+1,2,3,400
 1,2,4,400
-2,1,1,100
-2,1,2,100
-2,1,3,100
-2,1,4,100
-2,2,1,200
-2,2,2,200
-2,2,3,200
-2,2,4,200
-3,1,1,100
-3,1,2,85.8109
-3,1,3,100
-3,1,4,100
-3,2,1,200
-3,2,2,185.811
-3,2,3,200
-3,2,4,200
+2,1,1,300
+2,1,2,226.431
+2,1,3,300
+2,1,4,300
+2,2,1,400
+2,2,2,326.431
+2,2,3,400
+2,2,4,400
+3,1,1,178.808
+3,1,2,300
+3,1,3,300
+3,1,4,300
+3,2,1,278.808
+3,2,2,400
+3,2,3,400
+3,2,4,400

test/case_07/supply_function_equilibrium/expected_outputs/deficit_ex_ante_physical.csv

@@ -1,25 +1,25 @@
 period,scenario,subperiod,dem_1
 1,1,1,300
-1,1,2,150
-1,1,3,278.378
+1,1,2,289.266
+1,1,3,300
 1,1,4,300
 1,2,1,400
-1,2,2,250
-1,2,3,378.378
+1,2,2,389.266
+1,2,3,400
 1,2,4,400
-2,1,1,100
-2,1,2,100
-2,1,3,100
-2,1,4,100
-2,2,1,200
-2,2,2,200
-2,2,3,200
-2,2,4,200
-3,1,1,100
-3,1,2,85.8109
-3,1,3,100
-3,1,4,100
-3,2,1,200
-3,2,2,185.811
-3,2,3,200
-3,2,4,200
+2,1,1,300
+2,1,2,226.431
+2,1,3,300
+2,1,4,300
+2,2,1,400
+2,2,2,326.431
+2,2,3,400
+2,2,4,400
+3,1,1,178.808
+3,1,2,300
+3,1,3,300
+3,1,4,300
+3,2,1,278.808
+3,2,2,400
+3,2,3,400
+3,2,4,400