From c008cf3c2e21584a756fdf7a92d1c507af80f54b Mon Sep 17 00:00:00 2001 From: Cameron Wade Date: Fri, 9 Oct 2026 10:55:24 -0300 Subject: [PATCH] fix: limit storage charging by energy drawn (#376) Co-Authored-By: Claude Opus 5.5 --- temoa/components/storage.py | 19 ++++-- tests/conftest.py | 1 + tests/legacy_test_values.py | 3 +- tests/test_full_runs.py | 6 +- tests/test_storage.py | 66 +++++++++++++++++++ .../config_storage_charge_rate.toml | 11 ++++ tests/testing_data/storage_charge_rate.sql | 61 +++++++++++++++++ 7 files changed, 158 insertions(+), 9 deletions(-) create mode 100644 tests/testing_configs/config_storage_charge_rate.toml create mode 100644 tests/testing_data/storage_charge_rate.sql diff --git a/temoa/components/storage.py b/temoa/components/storage.py index c7415c96e..545ee1cd6 100644 --- a/temoa/components/storage.py +++ b/temoa/components/storage.py @@ -429,10 +429,16 @@ def storage_charge_rate_constraint( This constraint ensures that the charge rate of the storage unit is limited by the power capacity (typically GW) of the storage unit. + The limit applies to the energy drawn from the input commodity, before the + storage efficiency is applied. Capacity therefore limits charging and + discharging on the same side of the storage losses: what is drawn from the + input commodity and what is delivered to the output commodity. Losses are + still taken on charging, in :code:`storage_energy_constraint`. + .. math:: :label: StorageChargeRate - \sum_{I, O} \textbf{FIS}_{r, p, s, d, i, t, v, o} \cdot EFF_{r,i,t,v,o} + \sum_{I, O} \textbf{FIS}_{r, p, s, d, i, t, v, o} \le \textbf{CAP}_{r,p,t,v} \cdot C2A_{r,t} \cdot SEG_{s,d} @@ -440,10 +446,9 @@ def storage_charge_rate_constraint( \forall \{r, p, s, d, t, v\} \in \Theta_{\text{StorageChargeRate}} """ - # Calculate energy charge in each time slice + # Energy drawn from the input commodity in each time slice slice_charge = quicksum( model.v_flow_in[r, p, s, d, S_i, t, v, S_o] - * get_variable_efficiency(model, r, p, s, d, S_i, t, v, S_o) for S_i in model.process_inputs[r, p, t, v] for S_o in model.process_outputs_by_input[r, p, t, v, S_i] ) @@ -488,14 +493,16 @@ def storage_throughput_constraint( It is not enough to only limit the charge and discharge rate separately. We also need to ensure that the maximum throughput (charge + discharge) does not exceed - the capacity (typically GW) of the storage unit. + the capacity (typically GW) of the storage unit. As in + :code:`storage_charge_rate_constraint`, charging is counted as the energy drawn + from the input commodity, before the storage efficiency is applied. .. math:: :label: StorageThroughput \sum_{I, O} \textbf{FO}_{r, p, s, d, i, t, v, o} + - \sum_{I, O} \textbf{FIS}_{r, p, s, d, i, t, v, o} \cdot EFF_{r,i,t,v,o} + \sum_{I, O} \textbf{FIS}_{r, p, s, d, i, t, v, o} \le \textbf{CAP}_{r,p,t,v} \cdot C2A_{r,t} \cdot SEG_{s,d} @@ -508,9 +515,9 @@ def storage_throughput_constraint( for S_i in model.process_inputs_by_output[r, p, t, v, S_o] ) + # Energy drawn from the input commodity, as in storage_charge_rate_constraint charge = quicksum( model.v_flow_in[r, p, s, d, S_i, t, v, S_o] - * get_variable_efficiency(model, r, p, s, d, S_i, t, v, S_o) for S_i in model.process_inputs[r, p, t, v] for S_o in model.process_outputs_by_input[r, p, t, v, S_i] ) diff --git a/tests/conftest.py b/tests/conftest.py index 7cb63eb9a..f35331ad8 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -103,6 +103,7 @@ def refresh_databases() -> None: ('materials.sql', 'materials.sqlite'), ('simple_linked_tech.sql', 'simple_linked_tech.sqlite'), ('storageville.sql', 'storageville.sqlite'), + ('storage_charge_rate.sql', 'storage_charge_rate.sqlite'), ('test_week.sql', 'test_week.sqlite'), (None, 'utopia_stochastic.sqlite'), ] diff --git a/tests/legacy_test_values.py b/tests/legacy_test_values.py index 3fb7d2e13..46015446a 100644 --- a/tests/legacy_test_values.py +++ b/tests/legacy_test_values.py @@ -45,7 +45,8 @@ class ExpectedVals(Enum): # reduced after changing fixed costs from MLP to PL # reduced by <1 after changing season definition (segfrac no longer rounded) # increased by 143 after activating CF constraints for storage techs - ExpectedVals.OBJ_VALUE: 34853.6835, + # increased by 162 after limiting storage charging at the draw (#376) + ExpectedVals.OBJ_VALUE: 35015.9790, ExpectedVals.EFF_DOMAIN_SIZE: 12312, ExpectedVals.EFF_INDEX_SIZE: 64, # reduced 3/27: unlim_cap techs now employed. diff --git a/tests/test_full_runs.py b/tests/test_full_runs.py index 05d4f586a..88fb5f4d4 100644 --- a/tests/test_full_runs.py +++ b/tests/test_full_runs.py @@ -171,7 +171,8 @@ def test_myopic_utopia( # increased after rework of inter-season sequencing # reduced by <1 after changing season definition (segfrac no longer rounded) # decreased by 41 after activating CF constraints for storage techs - assert invest_sum == pytest.approx(10963.1018), ( + # decreased by 27 after limiting storage charging at the draw (#376) + assert invest_sum == pytest.approx(10935.8768), ( 'sum of investment costs did not match expected' ) @@ -275,7 +276,8 @@ def test_stochastic_utopia( # Stochastic Expected Value for current utopia configuration # reduced by <1 after changing season definition (segfrac no longer rounded) # increased by 145 after activating CF constraints for storage techs - expected_obj = 34534.1822 + # increased by 171 after limiting storage charging at the draw (#376) + expected_obj = 34705.3051 assert sequencer.stochastic_sequencer is not None assert sequencer.stochastic_sequencer.objective_value == pytest.approx(expected_obj, rel=1e-5) diff --git a/tests/test_storage.py b/tests/test_storage.py index 30a2827fe..6c14fe55b 100644 --- a/tests/test_storage.py +++ b/tests/test_storage.py @@ -7,8 +7,10 @@ from typing import Any import pytest +from pyomo.environ import value from temoa.core.model import TemoaModel +from temoa.types.core_types import Period, Region, Season, Technology, TimeOfDay, Vintage logger = logging.getLogger(__name__) # suitable scenarios for storage testing....singleton for now. @@ -153,6 +155,70 @@ def test_storage_flow_balance(system_test_run: tuple[str, Any, TemoaModel, Any]) ) +charge_rate_config_files = [ + {'name': 'storage_charge_rate', 'filename': 'config_storage_charge_rate.toml'}, +] + + +@pytest.mark.parametrize( + 'system_test_run', + argvalues=charge_rate_config_files, + indirect=True, + ids=[d['name'] for d in charge_rate_config_files], +) +def test_storage_charge_limited_at_draw(system_test_run: tuple[str, Any, TemoaModel, Any]) -> None: + """ + Capacity limits the energy a storage process draws while charging, not the energy + that enters storage after the efficiency is applied (issue #376). + + The battery has capacity 1, efficiency 0.8 and C2A 1, and each of the four time slices + is a quarter of the year, so each slice's limit is 0.25. Free solar is available only in + d1 and gas is expensive, so the battery charges as fast as it is allowed in d1. + """ + model: TemoaModel + _, _, model, _ = system_test_run + r, p, t, v = Region('R'), Period(2025), Technology('batt'), Vintage(2025) + + assert model.v_capacity[r, p, t, v].value == pytest.approx(1.0) + + total_drawn = 0.0 + total_delivered = 0.0 + for s in model.time_season: + for d in model.time_of_day: + drawn = sum( + model.v_flow_in[r, p, s, d, S_i, t, v, S_o].value + for S_i in model.process_inputs[r, p, t, v] + for S_o in model.process_outputs_by_input[r, p, t, v, S_i] + ) + delivered = sum( + model.v_flow_out[r, p, s, d, S_i, t, v, S_o].value + for S_o in model.process_outputs[r, p, t, v] + for S_i in model.process_inputs_by_output[r, p, t, v, S_o] + ) + limit = ( + model.v_capacity[r, p, t, v].value + * value(model.capacity_to_activity[r, t]) + * value(model.segment_fraction[s, d]) + ) + assert drawn <= limit + 1e-6, f'draw {drawn} exceeds the limit {limit} in {s, d}' + assert drawn + delivered <= limit + 1e-6, ( + f'draw plus discharge {drawn + delivered} exceeds the limit {limit} in {s, d}' + ) + total_drawn += drawn + total_delivered += delivered + + # The battery charges at its limit in d1 (0.3125 = 0.25 / 0.8 before the change) + d1_drawn = sum( + model.v_flow_in[r, p, Season('s1'), TimeOfDay('d1'), S_i, t, v, S_o].value + for S_i in model.process_inputs[r, p, t, v] + for S_o in model.process_outputs_by_input[r, p, t, v, S_i] + ) + assert d1_drawn == pytest.approx(0.25, abs=1e-6) + # Losses are still taken on charging + assert total_drawn == pytest.approx(0.25, abs=1e-6) + assert total_delivered == pytest.approx(0.8 * 0.25, abs=1e-6) + + # devnote: the storage_init constraint was reworked into LimitStorageLevelFraction # @pytest.mark.skip('not ready for primetime') # def test_hard_initialization(): diff --git a/tests/testing_configs/config_storage_charge_rate.toml b/tests/testing_configs/config_storage_charge_rate.toml new file mode 100644 index 000000000..4223bb2c2 --- /dev/null +++ b/tests/testing_configs/config_storage_charge_rate.toml @@ -0,0 +1,11 @@ +scenario = "storage charge rate" +scenario_mode = "perfect_foresight" +input_database = "tests/testing_outputs/storage_charge_rate.sqlite" +output_database = "tests/testing_outputs/storage_charge_rate.sqlite" +neos = false +solver_name = "appsi_highs" +save_excel = false +save_duals = false +save_lp_file = false +time_sequencing = "seasonal_timeslices" +days_per_period = 365 diff --git a/tests/testing_data/storage_charge_rate.sql b/tests/testing_data/storage_charge_rate.sql new file mode 100644 index 000000000..4b66a595b --- /dev/null +++ b/tests/testing_data/storage_charge_rate.sql @@ -0,0 +1,61 @@ +-- One battery charged from free solar available only in d1 (issue #376). +-- Capacity 1, efficiency 0.8, C2A 1 and four slices of a quarter year each: +-- the battery may draw at most 0.25 in d1. +REPLACE INTO "metadata" VALUES('DB_MAJOR',4,''); +REPLACE INTO "metadata" VALUES('DB_MINOR',0,''); +REPLACE INTO "metadata_real" VALUES('default_loan_rate',0.05,''); +REPLACE INTO "metadata_real" VALUES('global_discount_rate',0.05,''); +REPLACE INTO "operator" VALUES('e','equal to'); +REPLACE INTO "operator" VALUES('le','less than or equal to'); +REPLACE INTO "operator" VALUES('ge','greater than or equal to'); +REPLACE INTO "commodity_type" VALUES('p','physical commodity'); +REPLACE INTO "commodity_type" VALUES('d','demand commodity'); +REPLACE INTO "commodity_type" VALUES('s','source commodity'); +REPLACE INTO "commodity_type" VALUES('e','emissions commodity'); +REPLACE INTO "technology_type" VALUES('p','production technology'); +REPLACE INTO "technology_type" VALUES('ps','storage production technology'); +REPLACE INTO "sector_label" VALUES('electric',NULL); +REPLACE INTO "region" VALUES('R',NULL); +REPLACE INTO "time_period_type" VALUES('e','existing vintages'); +REPLACE INTO "time_period_type" VALUES('f','future'); +REPLACE INTO "time_period" VALUES(1,2020,'e'); +REPLACE INTO "time_period" VALUES(2,2025,'f'); +REPLACE INTO "time_period" VALUES(3,2030,'f'); +REPLACE INTO "time_season" VALUES(1,'s1',1.0,NULL); +REPLACE INTO "time_of_day" VALUES(1,'d1',6.0,NULL); +REPLACE INTO "time_of_day" VALUES(2,'d2',6.0,NULL); +REPLACE INTO "time_of_day" VALUES(3,'d3',6.0,NULL); +REPLACE INTO "time_of_day" VALUES(4,'d4',6.0,NULL); +REPLACE INTO "commodity" VALUES('earth','s','source',NULL); +REPLACE INTO "commodity" VALUES('ELC','p','electricity',NULL); +REPLACE INTO "commodity" VALUES('RL','d','demand',NULL); +REPLACE INTO "technology" VALUES('sun','p','electric','solar','',0,0,0,0,0,0,0,0,'free, d1 only'); +REPLACE INTO "technology" VALUES('gas','p','electric','gas','',0,0,0,0,0,0,0,0,'expensive'); +REPLACE INTO "technology" VALUES('batt','ps','electric','battery','',0,0,0,0,0,0,0,0,'battery'); +REPLACE INTO "technology" VALUES('load','p','electric','load','',0,0,0,0,0,0,0,0,'to demand'); +REPLACE INTO "efficiency" VALUES('R','earth','sun',2025,'ELC',1.0,NULL,NULL); +REPLACE INTO "efficiency" VALUES('R','earth','gas',2025,'ELC',1.0,NULL,NULL); +REPLACE INTO "efficiency" VALUES('R','ELC','batt',2025,'ELC',0.8,NULL,NULL); +REPLACE INTO "efficiency" VALUES('R','ELC','load',2025,'RL',1.0,NULL,NULL); +REPLACE INTO "capacity_to_activity" VALUES('R','sun',1.0,NULL,''); +REPLACE INTO "capacity_to_activity" VALUES('R','gas',1.0,NULL,''); +REPLACE INTO "capacity_to_activity" VALUES('R','batt',1.0,NULL,''); +REPLACE INTO "capacity_to_activity" VALUES('R','load',1.0,NULL,''); +REPLACE INTO "capacity_factor_tech" VALUES('R','s1','d1','sun',1.0,NULL); +REPLACE INTO "capacity_factor_tech" VALUES('R','s1','d2','sun',0.0,NULL); +REPLACE INTO "capacity_factor_tech" VALUES('R','s1','d3','sun',0.0,NULL); +REPLACE INTO "capacity_factor_tech" VALUES('R','s1','d4','sun',0.0,NULL); +REPLACE INTO "storage_duration" VALUES('R','batt',24.0,NULL); +REPLACE INTO "cost_invest" VALUES('R','sun',2025,0.001,NULL,NULL); +REPLACE INTO "cost_invest" VALUES('R','gas',2025,0.001,NULL,NULL); +REPLACE INTO "cost_invest" VALUES('R','batt',2025,0.001,NULL,NULL); +REPLACE INTO "cost_invest" VALUES('R','load',2025,0.001,NULL,NULL); +REPLACE INTO "cost_variable" VALUES('R',2025,'gas',2025,100.0,'',''); +REPLACE INTO "limit_capacity" VALUES('R',2025,'sun','le',10.0,'',''); +REPLACE INTO "limit_capacity" VALUES('R',2025,'batt','ge',1.0,'',''); +REPLACE INTO "limit_capacity" VALUES('R',2025,'batt','le',1.0,'',''); +REPLACE INTO "demand" VALUES('R',2025,'RL',4.0,'',''); +REPLACE INTO "demand_specific_distribution" VALUES('R',2025,'s1','d1','RL',0.25,NULL); +REPLACE INTO "demand_specific_distribution" VALUES('R',2025,'s1','d2','RL',0.25,NULL); +REPLACE INTO "demand_specific_distribution" VALUES('R',2025,'s1','d3','RL',0.25,NULL); +REPLACE INTO "demand_specific_distribution" VALUES('R',2025,'s1','d4','RL',0.25,NULL);