From cc7f4e818c9dd2790803ec964124268e82752321 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Wed, 4 Feb 2026 14:50:44 +0100 Subject: [PATCH 01/23] Wrote basic code infrastructure. --- src/mesido/esdl/esdl_heat_model.py | 7 +- .../pycml/component_library/milp/__init__.py | 2 + .../milp/heat/geothermal_source_elec.py | 39 ++++ .../model/sourcesink_with_geo_elec.esdl | 70 +++++++ tests/test_geothermal_source_elec.py | 187 ++++++++++++++++++ 5 files changed, 304 insertions(+), 1 deletion(-) create mode 100644 src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py create mode 100644 tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl create mode 100644 tests/test_geothermal_source_elec.py diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 40b1f35a..7c9502ac 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -40,6 +40,7 @@ GasSubstation, GasTankStorage, GeothermalSource, + GeothermalSourceElec, HeatBuffer, HeatDemand, HeatExchanger, @@ -1313,7 +1314,11 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS f"'{asset.name}' will not be actuated in a constant manner" ) - return GeothermalSource, modifiers + if len(asset.in_ports) == 2: + modifiers["cop"] = asset.attributes["COP"] + return GeothermalSourceElec, modifiers + else: + return GeothermalSource, modifiers elif asset.asset_type == "HeatPump": modifiers["cop"] = asset.attributes["COP"] return AirWaterHeatPump, modifiers diff --git a/src/mesido/pycml/component_library/milp/__init__.py b/src/mesido/pycml/component_library/milp/__init__.py index 6557874c..77e4a40c 100644 --- a/src/mesido/pycml/component_library/milp/__init__.py +++ b/src/mesido/pycml/component_library/milp/__init__.py @@ -21,6 +21,7 @@ from .heat.cold_demand import ColdDemand from .heat.control_valve import ControlValve from .heat.geothermal_source import GeothermalSource +from .heat.geothermal_source_elec import GeothermalSourceElec from .heat.heat_buffer import HeatBuffer from .heat.heat_demand import HeatDemand from .heat.heat_exchanger import HeatExchanger @@ -68,6 +69,7 @@ "GasSubstation", "GasTankStorage", "GeothermalSource", + "GeothermalSourceElec", "HeatExchanger", "HeatFourPort", "HeatPipe", diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py new file mode 100644 index 00000000..58d7f736 --- /dev/null +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -0,0 +1,39 @@ +from mesido.pycml import Variable +from mesido.pycml.component_library.milp.electricity.electricity_base import ElectricityPort +from mesido.pycml.component_library.milp.heat.geothermal_source import GeothermalSource +from mesido.pycml.pycml_mixin import add_variables_documentation_automatically + +from numpy import nan + + +@add_variables_documentation_automatically +class GeothermalSourceElec(GeothermalSource): + """ + TODO: document + + Variables created: + {add_variable_names_for_documentation_here} + + Parameters: + name : The name of the asset. \n + modifiers : Dictionary with asset information. + """ + + def __init__(self, name, **modifiers): + super().__init__( + name, + **modifiers, + ) + + self.elec_power_nominal = nan + self.cop = modifiers["cop"] + self.component_subtype = "geothermal_source_elec" + + self.add_variable(ElectricityPort, "ElectricityIn") + self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) + + self.add_equation(((self.ElectricityIn.Power - self.Power_elec) / self.elec_power_nominal)) + + self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) + + diff --git a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl new file mode 100644 index 00000000..5d785a3d --- /dev/null +++ b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl @@ -0,0 +1,70 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/test_geothermal_source_elec.py b/tests/test_geothermal_source_elec.py new file mode 100644 index 00000000..016eaee6 --- /dev/null +++ b/tests/test_geothermal_source_elec.py @@ -0,0 +1,187 @@ +from pathlib import Path +from unittest import TestCase + +from mesido.esdl.esdl_parser import ESDLFileParser +from mesido.esdl.profile_parser import ProfileReaderFromFile +from mesido.util import run_esdl_mesido_optimization + +import numpy as np + +from utils_tests import ( + demand_matching_test, + electric_power_conservation_test, + energy_conservation_test, + heat_to_discharge_test, +) + + +class TestGeothermalSourceElec(TestCase): + + # def asset_cost_calculation_tests(self, solution, results): + # # Check the cost components + # for asset in [ + # *solution.energy_system_components.get("heat_source", []), + # ]: + # esdl_asset = solution.esdl_assets[solution.esdl_asset_name_to_id_map[f"{asset}"]] + # costs_esdl_asset = esdl_asset.attributes["costInformation"] + + # # investment cost + # investment_cost_info = costs_esdl_asset.investmentCosts.value + # investment_cost = investment_cost_info * results[f"{asset}__max_size"] / 1.0e6 + # np.testing.assert_allclose( + # investment_cost, results[f"{asset}__investment_cost"], atol=1.0e-7 + # ) + + # # installation cost + # installation_cost = costs_esdl_asset.installationCosts.value + # np.testing.assert_allclose(installation_cost, results[f"{asset}__installation_cost"]) + + # # variable operational cost + # timesteps_hr = np.diff(solution.times()) / 3600.0 + # variable_operational_cost = 0.0 + # var_op_costs = costs_esdl_asset.variableOperationalCosts.value + # for ii in range(1, len(solution.times())): + # variable_operational_cost += ( + # var_op_costs + # * results[f"{asset}.Power_consumed"][ii] + # * timesteps_hr[ii - 1] + # / 1e6 + # ) + # np.testing.assert_allclose( + # variable_operational_cost, results[f"{asset}__variable_operational_cost"] + # ) + + # def test_elec_heat_source_elec(self): + # """ + # This tests checks the elec heat sourc elec for the standard checks and the energy + # conservation over the commodity change. + + # Checks: + # 1. demand is matched + # 2. energy conservation in the network + # 3. heat to discharge + # 4. energy conservation over the heat and electricity commodity + # 5. ElectricBoiler cost components + # """ + # import models.source_pipe_sink.src.double_pipe_heat as example + # from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + # base_folder = Path(example.__file__).resolve().parent.parent + + # heat_problem = run_esdl_mesido_optimization( + # SourcePipeSink, + # base_folder=base_folder, + # esdl_file_name="sourcesink_with_eboiler.esdl", + # esdl_parser=ESDLFileParser, + # profile_reader=ProfileReaderFromFile, + # input_timeseries_file="timeseries_import.csv", + # ) + # results = heat_problem.extract_results() + # parameters = heat_problem.parameters(0) + + # demand_matching_test(heat_problem, results) + # energy_conservation_test(heat_problem, results) + # heat_to_discharge_test(heat_problem, results) + # electric_power_conservation_test(heat_problem, results) + + # np.testing.assert_array_less(0.0, results["ElectricBoiler_9aab.Heat_source"]) + # np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) + # np.testing.assert_allclose( + # parameters["ElectricBoiler_9aab.efficiency"] + # * results["ElectricBoiler_9aab.Power_consumed"], + # results["ElectricBoiler_9aab.Heat_source"], + # ) + + # # Check the cost components of ElectricBoiler + # self.asset_cost_calculation_tests(heat_problem, results) + + # def test_heat_source_elec(self): + # """ + # This tests checks the heat sourc elec for the standard checks and the energy conservation + # over the commodity change. + + # Checks: + # 1. demand is matched + # 2. energy conservation in the network + # 3. heat to discharge + # 4. energy conservation over the heat and electricity commodity + # 5. ElectricBoiler cost components + # """ + # import models.source_pipe_sink.src.double_pipe_heat as example + # from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + # base_folder = Path(example.__file__).resolve().parent.parent + + # heat_problem = run_esdl_mesido_optimization( + # SourcePipeSink, + # base_folder=base_folder, + # esdl_file_name="sourcesink_with_eboiler_no_elec.esdl", + # esdl_parser=ESDLFileParser, + # profile_reader=ProfileReaderFromFile, + # input_timeseries_file="timeseries_import.csv", + # ) + # results = heat_problem.extract_results() + # parameters = heat_problem.parameters(0) + + # demand_matching_test(heat_problem, results) + # energy_conservation_test(heat_problem, results) + # heat_to_discharge_test(heat_problem, results) + + # np.testing.assert_array_less(0.0, results["ElectricBoiler_9aab.Heat_source"]) + # np.testing.assert_allclose( + # parameters["ElectricBoiler_9aab.efficiency"] + # * results["ElectricBoiler_9aab.Power_consumed"], + # results["ElectricBoiler_9aab.Heat_source"], + # ) + + # # Check the cost components of ElectricBoiler + # self.asset_cost_calculation_tests(heat_problem, results) + + def test_geothermal_source_elec(self): + """ + TODO: document. + + Checks: + 1. TODO + """ + # import models.source_pipe_sink.src.double_pipe_heat as example + # from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + # base_folder = Path(example.__file__).resolve().parent.parent + + # heat_problem = run_esdl_mesido_optimization( + # SourcePipeSink, + # base_folder=base_folder, + # esdl_file_name="sourcesink_withHP.esdl", + # esdl_parser=ESDLFileParser, + # profile_reader=ProfileReaderFromFile, + # input_timeseries_file="timeseries_import.csv", + # ) + # results = heat_problem.extract_results() + # parameters = heat_problem.parameters(0) + + # demand_matching_test(heat_problem, results) + # energy_conservation_test(heat_problem, results) + # heat_to_discharge_test(heat_problem, results) + # electric_power_conservation_test(heat_problem, results) + + # np.testing.assert_array_less(0.0, results["HeatPump_d8fd.Heat_source"]) + # np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) + # np.testing.assert_array_less( + # parameters["HeatPump_d8fd.cop"] * results["HeatPump_d8fd.Power_elec"], + # results["HeatPump_d8fd.Heat_source"] + 1.0e-6, + # ) + + # # Check how variable operation cost is calculated + # np.testing.assert_allclose( + # parameters["HeatPump_d8fd.variable_operational_cost_coefficient"] + # * sum(results["HeatPump_d8fd.Heat_source"][1:]) + # / parameters["HeatPump_d8fd.cop"], + # results["HeatPump_d8fd__variable_operational_cost"], + # ) + pass + +if __name__ == "__main__": + + TestElecBoiler = TestGeothermalSourceElec() + TestElecBoiler.test_geothermal_source_elec() From cadf89fe580ec33c2a7de249e9a525997c8cfc36 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Wed, 4 Feb 2026 15:11:33 +0100 Subject: [PATCH 02/23] Created test and wrote esdl parsin routine. Equations need to be fixed. --- src/mesido/esdl/esdl_model_base.py | 22 ++++++++++ .../model/sourcesink_with_geo_elec.esdl | 4 +- tests/test_geothermal_source_elec.py | 41 +++++++++---------- 3 files changed, 44 insertions(+), 23 deletions(-) diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index 66b06d94..49e8492b 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -258,6 +258,28 @@ def __set_primary_secondary_heat_ports(): raise Exception( f"{asset.name} must have one inport for electricity and one outport for gas" ) + + elif ( + asset.asset_type == "GeothermalSource" + and len(asset.out_ports) == 1 + and len(asset.in_ports) == 2 + ): + for p in [*asset.in_ports, *asset.out_ports]: + + if isinstance(p, InPort) and isinstance( + p.carrier, esdl.ElectricityCommodity + ): + port_map[p.id] = getattr(component, elec_in_suf) + elif isinstance(p, InPort) and isinstance(p.carrier, esdl.HeatCommodity): + port_map[p.id] = getattr(component, in_suf) + elif isinstance(p, OutPort): # OutPort + port_map[p.id] = getattr(component, out_suf) + else: + raise Exception( + f"{asset.name} has does not have (1 electricity in_port) 1 heat " + f"in port and 1 Heat out_ports " + ) + elif ( asset.in_ports is None and isinstance(asset.out_ports[0].carrier, esdl.ElectricityCommodity) diff --git a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl index 5d785a3d..790bd4e6 100644 --- a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl +++ b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl @@ -1,5 +1,5 @@ - + @@ -52,7 +52,7 @@ - + diff --git a/tests/test_geothermal_source_elec.py b/tests/test_geothermal_source_elec.py index 016eaee6..69cb62bd 100644 --- a/tests/test_geothermal_source_elec.py +++ b/tests/test_geothermal_source_elec.py @@ -144,26 +144,26 @@ def test_geothermal_source_elec(self): Checks: 1. TODO """ - # import models.source_pipe_sink.src.double_pipe_heat as example - # from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink - - # base_folder = Path(example.__file__).resolve().parent.parent - - # heat_problem = run_esdl_mesido_optimization( - # SourcePipeSink, - # base_folder=base_folder, - # esdl_file_name="sourcesink_withHP.esdl", - # esdl_parser=ESDLFileParser, - # profile_reader=ProfileReaderFromFile, - # input_timeseries_file="timeseries_import.csv", - # ) - # results = heat_problem.extract_results() - # parameters = heat_problem.parameters(0) - - # demand_matching_test(heat_problem, results) - # energy_conservation_test(heat_problem, results) - # heat_to_discharge_test(heat_problem, results) - # electric_power_conservation_test(heat_problem, results) + import models.source_pipe_sink.src.double_pipe_heat as example + from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + base_folder = Path(example.__file__).resolve().parent.parent + + heat_problem = run_esdl_mesido_optimization( + SourcePipeSink, + base_folder=base_folder, + esdl_file_name="sourcesink_with_geo_elec.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_import.csv", + ) + results = heat_problem.extract_results() + parameters = heat_problem.parameters(0) + + demand_matching_test(heat_problem, results) + energy_conservation_test(heat_problem, results) + heat_to_discharge_test(heat_problem, results) + electric_power_conservation_test(heat_problem, results) # np.testing.assert_array_less(0.0, results["HeatPump_d8fd.Heat_source"]) # np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) @@ -179,7 +179,6 @@ def test_geothermal_source_elec(self): # / parameters["HeatPump_d8fd.cop"], # results["HeatPump_d8fd__variable_operational_cost"], # ) - pass if __name__ == "__main__": From b7a2c1dd8730150a7c9f753adf4c10c043ea0ec8 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Wed, 4 Feb 2026 16:55:08 +0100 Subject: [PATCH 03/23] Fixed some esdl implementation errors. Remaining issue with conservation of electricity. --- src/mesido/esdl/asset_to_component_base.py | 1 + src/mesido/esdl/esdl_heat_model.py | 58 ++++++++++++++++++- .../milp/heat/geothermal_source_elec.py | 3 + 3 files changed, 61 insertions(+), 1 deletion(-) diff --git a/src/mesido/esdl/asset_to_component_base.py b/src/mesido/esdl/asset_to_component_base.py index b739209b..c402bc51 100644 --- a/src/mesido/esdl/asset_to_component_base.py +++ b/src/mesido/esdl/asset_to_component_base.py @@ -898,6 +898,7 @@ def _get_connected_i_nominal_and_max(self, asset: Asset) -> Tuple[float, float]: or asset.asset_type == "Electrolyzer" or asset.asset_type == "ElectricBoiler" or asset.asset_type == "HeatPump" + or asset.asset_type == "GeothermalSource" ): for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 7c9502ac..c2a341cf 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1315,8 +1315,9 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS ) if len(asset.in_ports) == 2: - modifiers["cop"] = asset.attributes["COP"] + _, modifiers = self.convert_geothermal_source_elec(asset) return GeothermalSourceElec, modifiers + else: return GeothermalSource, modifiers elif asset.asset_type == "HeatPump": @@ -2692,6 +2693,61 @@ def convert_air_water_heat_pump_elec( return AirWaterHeatPumpElec, modifiers + def convert_geothermal_source_elec( + self, asset: Asset + ) -> Tuple[AirWaterHeatPumpElec, MODIFIERS]: + """ + TODO: document + + """ + assert asset.asset_type in {"GeothermalSource"} + + get_potential_errors().add_potential_issue( + MesidoAssetIssueType.ELECT_ASSET_TYPE_NOT_SUPPORTED, + asset.id, + f"Asset named {asset.name}: This is an asset that includes electricity and it should be" + " replaced with a supported asset", + ) + + max_supply = asset.attributes["power"] + + if not max_supply: + logger.error(f"{asset.asset_type} '{asset.name}' has no max power specified. ") + assert max_supply > 0.0 + + id_mapping = asset.global_properties["carriers"][asset.in_ports[0].carrier.id][ + "id_number_mapping" + ] + + # TODO: CO2 coefficient + + q_nominal = self._get_connected_q_nominal(asset) + for port in asset.in_ports: + if isinstance(port.carrier, esdl.ElectricityCommodity): + min_voltage = port.carrier.voltage + i_max, i_nom = self._get_connected_i_nominal_and_max(asset) + cop = asset.attributes["COP"] if asset.attributes["COP"] else 1.0 + + modifiers = dict( + Heat_source=dict(min=0.0, max=max_supply, nominal=max_supply / 2.0), + id_mapping_carrier=id_mapping, + ElectricityIn=dict( + Power=dict(min=0.0, max=max_supply, nominal=max_supply / 2.0), + I=dict(min=0.0, max=i_max, nominal=i_nom), + V=dict(min=min_voltage, nominal=min_voltage), + ), + min_voltage=min_voltage, + elec_power_nominal=max_supply, + cop=cop, + **self._generic_modifiers(asset), + **self._generic_heat_modifiers(0.0, max_supply, q_nominal["Q_nominal"]), + **self._supply_return_temperature_modifiers(asset), + **self._rho_cp_modifiers, + **self._get_cost_figure_modifiers(asset), + ) + + return GeothermalSourceElec, modifiers + class ESDLHeatModel(_ESDLModelBase): """ diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index 58d7f736..d650f790 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -28,9 +28,12 @@ def __init__(self, name, **modifiers): self.elec_power_nominal = nan self.cop = modifiers["cop"] self.component_subtype = "geothermal_source_elec" + self.id_mapping_carrier = -1 + self.min_voltage = nan self.add_variable(ElectricityPort, "ElectricityIn") self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) + #self.add_variable(Variable, "Power_elec", min=0.0, nominal=10.0) self.add_equation(((self.ElectricityIn.Power - self.Power_elec) / self.elec_power_nominal)) From d49c6d0d1cdd464e91218918fe90bf96504a907d Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Wed, 4 Feb 2026 17:34:18 +0100 Subject: [PATCH 04/23] Implementation finished. Still need to document, format, lint etc. --- .../milp/heat/geothermal_source_elec.py | 2 -- tests/test_geothermal_source_elec.py | 15 --------------- tests/utils_tests.py | 1 + 3 files changed, 1 insertion(+), 17 deletions(-) diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index d650f790..8217378c 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -33,10 +33,8 @@ def __init__(self, name, **modifiers): self.add_variable(ElectricityPort, "ElectricityIn") self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) - #self.add_variable(Variable, "Power_elec", min=0.0, nominal=10.0) self.add_equation(((self.ElectricityIn.Power - self.Power_elec) / self.elec_power_nominal)) - self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) diff --git a/tests/test_geothermal_source_elec.py b/tests/test_geothermal_source_elec.py index 69cb62bd..71bba989 100644 --- a/tests/test_geothermal_source_elec.py +++ b/tests/test_geothermal_source_elec.py @@ -165,21 +165,6 @@ def test_geothermal_source_elec(self): heat_to_discharge_test(heat_problem, results) electric_power_conservation_test(heat_problem, results) - # np.testing.assert_array_less(0.0, results["HeatPump_d8fd.Heat_source"]) - # np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) - # np.testing.assert_array_less( - # parameters["HeatPump_d8fd.cop"] * results["HeatPump_d8fd.Power_elec"], - # results["HeatPump_d8fd.Heat_source"] + 1.0e-6, - # ) - - # # Check how variable operation cost is calculated - # np.testing.assert_allclose( - # parameters["HeatPump_d8fd.variable_operational_cost_coefficient"] - # * sum(results["HeatPump_d8fd.Heat_source"][1:]) - # / parameters["HeatPump_d8fd.cop"], - # results["HeatPump_d8fd__variable_operational_cost"], - # ) - if __name__ == "__main__": TestElecBoiler = TestGeothermalSourceElec() diff --git a/tests/utils_tests.py b/tests/utils_tests.py index b4211892..bd9d6d7b 100644 --- a/tests/utils_tests.py +++ b/tests/utils_tests.py @@ -370,6 +370,7 @@ def electric_power_conservation_test(solution, results, atol=1e-2): "elec_heat_source_elec", "heat_pump_elec", "air_water_heat_pump_elec", + "geothermal_source_elec", ] ) producers = solution.energy_system_components_get(["electricity_source"]) From 4d52db291c3c5630f6a234ed3a1cf11e6a1efcd7 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 5 Feb 2026 15:54:55 +0100 Subject: [PATCH 05/23] Wrote missing tests, some documentation. --- .../model/sourcesink_with_geo_elec.esdl | 2 +- tests/test_geothermal_source_elec.py | 152 ++++-------------- 2 files changed, 29 insertions(+), 125 deletions(-) diff --git a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl index 790bd4e6..c437cc6d 100644 --- a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl +++ b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl @@ -1,5 +1,5 @@ - + diff --git a/tests/test_geothermal_source_elec.py b/tests/test_geothermal_source_elec.py index 71bba989..47bf64e0 100644 --- a/tests/test_geothermal_source_elec.py +++ b/tests/test_geothermal_source_elec.py @@ -17,132 +17,18 @@ class TestGeothermalSourceElec(TestCase): - # def asset_cost_calculation_tests(self, solution, results): - # # Check the cost components - # for asset in [ - # *solution.energy_system_components.get("heat_source", []), - # ]: - # esdl_asset = solution.esdl_assets[solution.esdl_asset_name_to_id_map[f"{asset}"]] - # costs_esdl_asset = esdl_asset.attributes["costInformation"] - - # # investment cost - # investment_cost_info = costs_esdl_asset.investmentCosts.value - # investment_cost = investment_cost_info * results[f"{asset}__max_size"] / 1.0e6 - # np.testing.assert_allclose( - # investment_cost, results[f"{asset}__investment_cost"], atol=1.0e-7 - # ) - - # # installation cost - # installation_cost = costs_esdl_asset.installationCosts.value - # np.testing.assert_allclose(installation_cost, results[f"{asset}__installation_cost"]) - - # # variable operational cost - # timesteps_hr = np.diff(solution.times()) / 3600.0 - # variable_operational_cost = 0.0 - # var_op_costs = costs_esdl_asset.variableOperationalCosts.value - # for ii in range(1, len(solution.times())): - # variable_operational_cost += ( - # var_op_costs - # * results[f"{asset}.Power_consumed"][ii] - # * timesteps_hr[ii - 1] - # / 1e6 - # ) - # np.testing.assert_allclose( - # variable_operational_cost, results[f"{asset}__variable_operational_cost"] - # ) - - # def test_elec_heat_source_elec(self): - # """ - # This tests checks the elec heat sourc elec for the standard checks and the energy - # conservation over the commodity change. - - # Checks: - # 1. demand is matched - # 2. energy conservation in the network - # 3. heat to discharge - # 4. energy conservation over the heat and electricity commodity - # 5. ElectricBoiler cost components - # """ - # import models.source_pipe_sink.src.double_pipe_heat as example - # from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink - - # base_folder = Path(example.__file__).resolve().parent.parent - - # heat_problem = run_esdl_mesido_optimization( - # SourcePipeSink, - # base_folder=base_folder, - # esdl_file_name="sourcesink_with_eboiler.esdl", - # esdl_parser=ESDLFileParser, - # profile_reader=ProfileReaderFromFile, - # input_timeseries_file="timeseries_import.csv", - # ) - # results = heat_problem.extract_results() - # parameters = heat_problem.parameters(0) - - # demand_matching_test(heat_problem, results) - # energy_conservation_test(heat_problem, results) - # heat_to_discharge_test(heat_problem, results) - # electric_power_conservation_test(heat_problem, results) - - # np.testing.assert_array_less(0.0, results["ElectricBoiler_9aab.Heat_source"]) - # np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) - # np.testing.assert_allclose( - # parameters["ElectricBoiler_9aab.efficiency"] - # * results["ElectricBoiler_9aab.Power_consumed"], - # results["ElectricBoiler_9aab.Heat_source"], - # ) - - # # Check the cost components of ElectricBoiler - # self.asset_cost_calculation_tests(heat_problem, results) - - # def test_heat_source_elec(self): - # """ - # This tests checks the heat sourc elec for the standard checks and the energy conservation - # over the commodity change. - - # Checks: - # 1. demand is matched - # 2. energy conservation in the network - # 3. heat to discharge - # 4. energy conservation over the heat and electricity commodity - # 5. ElectricBoiler cost components - # """ - # import models.source_pipe_sink.src.double_pipe_heat as example - # from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink - - # base_folder = Path(example.__file__).resolve().parent.parent - - # heat_problem = run_esdl_mesido_optimization( - # SourcePipeSink, - # base_folder=base_folder, - # esdl_file_name="sourcesink_with_eboiler_no_elec.esdl", - # esdl_parser=ESDLFileParser, - # profile_reader=ProfileReaderFromFile, - # input_timeseries_file="timeseries_import.csv", - # ) - # results = heat_problem.extract_results() - # parameters = heat_problem.parameters(0) - - # demand_matching_test(heat_problem, results) - # energy_conservation_test(heat_problem, results) - # heat_to_discharge_test(heat_problem, results) - - # np.testing.assert_array_less(0.0, results["ElectricBoiler_9aab.Heat_source"]) - # np.testing.assert_allclose( - # parameters["ElectricBoiler_9aab.efficiency"] - # * results["ElectricBoiler_9aab.Power_consumed"], - # results["ElectricBoiler_9aab.Heat_source"], - # ) - - # # Check the cost components of ElectricBoiler - # self.asset_cost_calculation_tests(heat_problem, results) - def test_geothermal_source_elec(self): """ - TODO: document. + This tests checks the electric geothermal source for the standard checks and the energy conservation + over the commodity change. Checks: - 1. TODO + 1. Demand is matched. + 2. Energy is conserved in the network. + 3. Heat to discharge. + 4. Electric power is conserved. + 5. Checks that the asset equations behave as expected. + 6. Variable operational cost check. """ import models.source_pipe_sink.src.double_pipe_heat as example from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink @@ -160,12 +46,30 @@ def test_geothermal_source_elec(self): results = heat_problem.extract_results() parameters = heat_problem.parameters(0) + # Standard checks. demand_matching_test(heat_problem, results) energy_conservation_test(heat_problem, results) heat_to_discharge_test(heat_problem, results) electric_power_conservation_test(heat_problem, results) -if __name__ == "__main__": - + # Equations check + np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) + np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) + np.testing.assert_array_less( + parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], + results["GeothermalSource_a77b.Heat_source"] + 1.0e-6, + ) + + # Variable operational cost check. + np.testing.assert_allclose( + parameters["GeothermalSource_a77b.variable_operational_cost_coefficient"] + * sum(results["GeothermalSource_a77b.Heat_source"][1:]) + / parameters["GeothermalSource_a77b.cop"], + results["GeothermalSource_a77b__variable_operational_cost"], + ) + + +if __name__ == "__main__": TestElecBoiler = TestGeothermalSourceElec() TestElecBoiler.test_geothermal_source_elec() + From c7354b3fc0271b06855ad0e96c2085970a82cabe Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 5 Feb 2026 16:02:19 +0100 Subject: [PATCH 06/23] Minor fixes. Need to document, format and lint to complete PR submission. --- src/mesido/esdl/esdl_heat_model.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index c2a341cf..d432b2b6 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1313,8 +1313,8 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS f"{asset.asset_type} '{asset.name}' has no desired flow rate specified. " f"'{asset.name}' will not be actuated in a constant manner" ) - - if len(asset.in_ports) == 2: + # TODO: A better check might be needed here. Perhaps check for an elec asset. + if len(asset.in_ports) == 2: # Geothermal source with an electricity in port. _, modifiers = self.convert_geothermal_source_elec(asset) return GeothermalSourceElec, modifiers From 6d1c63b5b2f590766d73958f7a821f2af5db6b3f Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 26 Feb 2026 15:30:10 +0100 Subject: [PATCH 07/23] Implemented comments from the PR. --- src/mesido/esdl/esdl_heat_model.py | 73 ++++--------------- src/mesido/esdl/esdl_model_base.py | 16 ++-- .../milp/heat/geothermal_source.py | 7 ++ .../milp/heat/geothermal_source_elec.py | 5 -- tests/test_geothermal_source_elec.py | 1 + 5 files changed, 33 insertions(+), 69 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index d432b2b6..4c1e2caf 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1313,11 +1313,25 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS f"{asset.asset_type} '{asset.name}' has no desired flow rate specified. " f"'{asset.name}' will not be actuated in a constant manner" ) + max_supply=asset.attributes["power"] + modifiers["elec_power_nominal"] = max_supply + modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 1.0 + # TODO: A better check might be needed here. Perhaps check for an elec asset. if len(asset.in_ports) == 2: # Geothermal source with an electricity in port. - _, modifiers = self.convert_geothermal_source_elec(asset) + for port in asset.in_ports: + if isinstance(port.carrier, esdl.ElectricityCommodity): + min_voltage = port.carrier.voltage + i_max, i_nom = self._get_connected_i_nominal_and_max(asset) + modifiers.update( + min_voltage=min_voltage, + ElectricityIn=dict( + Power=dict(min=0.0, max=max_supply, nominal=max_supply / 2.0), + I=dict(min=0.0, max=i_max, nominal=i_nom), + V=dict(min=min_voltage, nominal=min_voltage), + ) + ) return GeothermalSourceElec, modifiers - else: return GeothermalSource, modifiers elif asset.asset_type == "HeatPump": @@ -2693,61 +2707,6 @@ def convert_air_water_heat_pump_elec( return AirWaterHeatPumpElec, modifiers - def convert_geothermal_source_elec( - self, asset: Asset - ) -> Tuple[AirWaterHeatPumpElec, MODIFIERS]: - """ - TODO: document - - """ - assert asset.asset_type in {"GeothermalSource"} - - get_potential_errors().add_potential_issue( - MesidoAssetIssueType.ELECT_ASSET_TYPE_NOT_SUPPORTED, - asset.id, - f"Asset named {asset.name}: This is an asset that includes electricity and it should be" - " replaced with a supported asset", - ) - - max_supply = asset.attributes["power"] - - if not max_supply: - logger.error(f"{asset.asset_type} '{asset.name}' has no max power specified. ") - assert max_supply > 0.0 - - id_mapping = asset.global_properties["carriers"][asset.in_ports[0].carrier.id][ - "id_number_mapping" - ] - - # TODO: CO2 coefficient - - q_nominal = self._get_connected_q_nominal(asset) - for port in asset.in_ports: - if isinstance(port.carrier, esdl.ElectricityCommodity): - min_voltage = port.carrier.voltage - i_max, i_nom = self._get_connected_i_nominal_and_max(asset) - cop = asset.attributes["COP"] if asset.attributes["COP"] else 1.0 - - modifiers = dict( - Heat_source=dict(min=0.0, max=max_supply, nominal=max_supply / 2.0), - id_mapping_carrier=id_mapping, - ElectricityIn=dict( - Power=dict(min=0.0, max=max_supply, nominal=max_supply / 2.0), - I=dict(min=0.0, max=i_max, nominal=i_nom), - V=dict(min=min_voltage, nominal=min_voltage), - ), - min_voltage=min_voltage, - elec_power_nominal=max_supply, - cop=cop, - **self._generic_modifiers(asset), - **self._generic_heat_modifiers(0.0, max_supply, q_nominal["Q_nominal"]), - **self._supply_return_temperature_modifiers(asset), - **self._rho_cp_modifiers, - **self._get_cost_figure_modifiers(asset), - ) - - return GeothermalSourceElec, modifiers - class ESDLHeatModel(_ESDLModelBase): """ diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index 49e8492b..a09c2771 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -180,7 +180,7 @@ def __set_primary_secondary_heat_ports(): f"milp(4) and electricity (1) ports" ) elif ( - asset.asset_type == "HeatPump" + (asset.asset_type == "HeatPump" or asset.asset_type == "GeothermalSource") and len(asset.out_ports) == 1 and len(asset.in_ports) in [1, 2] ): @@ -200,12 +200,14 @@ def __set_primary_secondary_heat_ports(): f"in port and 1 Heat out_ports " ) else: - raise Exception( - f"{asset.name} has incorrect number of in/out ports. HeatPumps are allows " - f"to have 1 in and 1 out port for air-water HP, 2 in ports and 2 out ports " - f"when modelling a water-water HP, or 3 in ports and 2 out ports when the " - f"electricity connection of the water-water HP is modelled." - ) + if asset.asset_type == "HeatPump": + raise Exception( + f"{asset.name} has incorrect number of in/out ports. HeatPumps are allows " + f"to have 1 in and 1 out port for air-water HP, 2 in ports and 2 out ports " + f"when modelling a water-water HP, or 3 in ports and 2 out ports when the " + f"electricity connection of the water-water HP is modelled." + ) + elif ( asset.asset_type == "GasHeater" and len(asset.out_ports) == 1 diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py index a91a4099..9f955211 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py @@ -1,3 +1,4 @@ +from mesido.pycml import Variable from mesido.pycml.pycml_mixin import add_variables_documentation_automatically from numpy import nan @@ -31,4 +32,10 @@ def __init__(self, name, **modifiers): self.target_flow_rate = nan self.single_doublet_power = nan + self.cop = nan + self.elec_power_nominal = nan self.nr_of_doublets = 1.0 + + self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) + + self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index 8217378c..72d2e8d5 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -1,4 +1,3 @@ -from mesido.pycml import Variable from mesido.pycml.component_library.milp.electricity.electricity_base import ElectricityPort from mesido.pycml.component_library.milp.heat.geothermal_source import GeothermalSource from mesido.pycml.pycml_mixin import add_variables_documentation_automatically @@ -25,16 +24,12 @@ def __init__(self, name, **modifiers): **modifiers, ) - self.elec_power_nominal = nan - self.cop = modifiers["cop"] self.component_subtype = "geothermal_source_elec" self.id_mapping_carrier = -1 self.min_voltage = nan self.add_variable(ElectricityPort, "ElectricityIn") - self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) self.add_equation(((self.ElectricityIn.Power - self.Power_elec) / self.elec_power_nominal)) - self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) diff --git a/tests/test_geothermal_source_elec.py b/tests/test_geothermal_source_elec.py index 47bf64e0..ea7cd06d 100644 --- a/tests/test_geothermal_source_elec.py +++ b/tests/test_geothermal_source_elec.py @@ -68,6 +68,7 @@ def test_geothermal_source_elec(self): results["GeothermalSource_a77b__variable_operational_cost"], ) + # TODO: move test to a better location, improve the equations checks, check both the geothermal asset and the one with the elec port. if __name__ == "__main__": TestElecBoiler = TestGeothermalSourceElec() From 291f69b2cae2b7fdf294c1e49867521731fab1a8 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 26 Feb 2026 16:35:04 +0100 Subject: [PATCH 08/23] Refactored and modified tests for geothermal source assets. --- .../model/sourcesink_with_geo.esdl | 57 ++++++++++ tests/test_geothermal_source_elec.py | 76 ------------- tests/test_multicommodity.py | 100 ++++++++++++++++++ 3 files changed, 157 insertions(+), 76 deletions(-) create mode 100644 tests/models/source_pipe_sink/model/sourcesink_with_geo.esdl delete mode 100644 tests/test_geothermal_source_elec.py diff --git a/tests/models/source_pipe_sink/model/sourcesink_with_geo.esdl b/tests/models/source_pipe_sink/model/sourcesink_with_geo.esdl new file mode 100644 index 00000000..f69ec905 --- /dev/null +++ b/tests/models/source_pipe_sink/model/sourcesink_with_geo.esdl @@ -0,0 +1,57 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/test_geothermal_source_elec.py b/tests/test_geothermal_source_elec.py deleted file mode 100644 index ea7cd06d..00000000 --- a/tests/test_geothermal_source_elec.py +++ /dev/null @@ -1,76 +0,0 @@ -from pathlib import Path -from unittest import TestCase - -from mesido.esdl.esdl_parser import ESDLFileParser -from mesido.esdl.profile_parser import ProfileReaderFromFile -from mesido.util import run_esdl_mesido_optimization - -import numpy as np - -from utils_tests import ( - demand_matching_test, - electric_power_conservation_test, - energy_conservation_test, - heat_to_discharge_test, -) - - -class TestGeothermalSourceElec(TestCase): - - def test_geothermal_source_elec(self): - """ - This tests checks the electric geothermal source for the standard checks and the energy conservation - over the commodity change. - - Checks: - 1. Demand is matched. - 2. Energy is conserved in the network. - 3. Heat to discharge. - 4. Electric power is conserved. - 5. Checks that the asset equations behave as expected. - 6. Variable operational cost check. - """ - import models.source_pipe_sink.src.double_pipe_heat as example - from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink - - base_folder = Path(example.__file__).resolve().parent.parent - - heat_problem = run_esdl_mesido_optimization( - SourcePipeSink, - base_folder=base_folder, - esdl_file_name="sourcesink_with_geo_elec.esdl", - esdl_parser=ESDLFileParser, - profile_reader=ProfileReaderFromFile, - input_timeseries_file="timeseries_import.csv", - ) - results = heat_problem.extract_results() - parameters = heat_problem.parameters(0) - - # Standard checks. - demand_matching_test(heat_problem, results) - energy_conservation_test(heat_problem, results) - heat_to_discharge_test(heat_problem, results) - electric_power_conservation_test(heat_problem, results) - - # Equations check - np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) - np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) - np.testing.assert_array_less( - parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], - results["GeothermalSource_a77b.Heat_source"] + 1.0e-6, - ) - - # Variable operational cost check. - np.testing.assert_allclose( - parameters["GeothermalSource_a77b.variable_operational_cost_coefficient"] - * sum(results["GeothermalSource_a77b.Heat_source"][1:]) - / parameters["GeothermalSource_a77b.cop"], - results["GeothermalSource_a77b__variable_operational_cost"], - ) - - # TODO: move test to a better location, improve the equations checks, check both the geothermal asset and the one with the elec port. - -if __name__ == "__main__": - TestElecBoiler = TestGeothermalSourceElec() - TestElecBoiler.test_geothermal_source_elec() - diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index d38742a9..f0ce4b64 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -323,9 +323,109 @@ def solver_options(self): ) np.testing.assert_allclose(var_opex_hp_calc, var_opex_hp) +class TestGeothermalSourceElec(TestCase): + + def test_geothermal_source(self): + """ + This tests the electric behavior of the regular geothermal source asset. + + It first does the standard checks, and then the equations that are added to the geothermal + asset, including the electricity power consumption. + + """ + import models.source_pipe_sink.src.double_pipe_heat as example + from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + base_folder = Path(example.__file__).resolve().parent.parent + + heat_problem = run_esdl_mesido_optimization( + SourcePipeSink, + base_folder=base_folder, + esdl_file_name="sourcesink_with_geo.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_import.csv", + ) + results = heat_problem.extract_results() + parameters = heat_problem.parameters(0) + + # Standard checks. + demand_matching_test(heat_problem, results) + energy_conservation_test(heat_problem, results) + heat_to_discharge_test(heat_problem, results) + electric_power_conservation_test(heat_problem, results) + + # Equations check + np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) + np.testing.assert_allclose( + parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], + results["GeothermalSource_a77b.Heat_source"], + ) + + # Variable operational cost check. + np.testing.assert_allclose( + parameters["GeothermalSource_a77b.variable_operational_cost_coefficient"] + * sum(results["GeothermalSource_a77b.Heat_source"][1:]) + / parameters["GeothermalSource_a77b.cop"], + results["GeothermalSource_a77b__variable_operational_cost"], + ) + + def test_geothermal_source_elec(self): + """ + This tests checks the electric geothermal producer asset with an electricity port. + + It does all the same checks as with the regular geothermal asset, plus tests on the + electricity in port. + """ + import models.source_pipe_sink.src.double_pipe_heat as example + from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + base_folder = Path(example.__file__).resolve().parent.parent + + heat_problem = run_esdl_mesido_optimization( + SourcePipeSink, + base_folder=base_folder, + esdl_file_name="sourcesink_with_geo_elec.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_import.csv", + ) + results = heat_problem.extract_results() + parameters = heat_problem.parameters(0) + + # Standard checks. + demand_matching_test(heat_problem, results) + energy_conservation_test(heat_problem, results) + heat_to_discharge_test(heat_problem, results) + electric_power_conservation_test(heat_problem, results) + + # Equations check + np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) + np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) + np.testing.assert_allclose( + parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], + results["GeothermalSource_a77b.Heat_source"], + ) + + # Test electricity port + np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.ElectricityIn.Power"]) + np.testing.assert_allclose(results["GeothermalSource_a77b.ElectricityIn.Power"], + results["GeothermalSource_a77b.Power_elec"]) + + # Variable operational cost check. + np.testing.assert_allclose( + parameters["GeothermalSource_a77b.variable_operational_cost_coefficient"] + * sum(results["GeothermalSource_a77b.Heat_source"][1:]) + / parameters["GeothermalSource_a77b.cop"], + results["GeothermalSource_a77b__variable_operational_cost"], + ) + if __name__ == "__main__": test_cold_demand = TestMultiCommodityHeatPump() test_cold_demand.test_air_to_water_heat_pump_elec_min_elec() # test_cold_demand.test_heat_pump_elec_min_elec() + test_geothermal = TestGeothermalSourceElec() + test_geothermal.test_geothermal_source_elec() + test_geothermal.test_geothermal_source() \ No newline at end of file From 75df94d5bab02786b6cf756db176e5413d29a2b9 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 15:23:01 +0100 Subject: [PATCH 09/23] Fixed comments from the PR. --- src/mesido/esdl/esdl_heat_model.py | 6 +++--- .../milp/heat/geothermal_source.py | 13 ++++++++----- .../milp/heat/geothermal_source_elec.py | 5 ++++- tests/test_multicommodity.py | 6 +++--- 4 files changed, 18 insertions(+), 12 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 4c1e2caf..01726223 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1316,9 +1316,9 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS max_supply=asset.attributes["power"] modifiers["elec_power_nominal"] = max_supply modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 1.0 - - # TODO: A better check might be needed here. Perhaps check for an elec asset. - if len(asset.in_ports) == 2: # Geothermal source with an electricity in port. + # Check to see if there is an electricity carrier at the in ports. + in_port_carriers = [port.carrier for port in asset.in_ports] + if any([isinstance(carrier, esdl.esdl.ElectricityCommodity) for carrier in in_port_carriers]): for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): min_voltage = port.carrier.voltage diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py index 9f955211..a398c94d 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py @@ -10,11 +10,14 @@ class GeothermalSource(HeatSource): """ The geothermal source component is used to model geothermal doublets. It is equivilent to a - normal source with the only difference being in the modelling of doublets. The main reason for - this component instead of using just a regular source is that to have the integer behaviour of - increasing the amount of doublets. In the HeatMixin an integer is created _aggregation_count to - model the amount of doublets and the maximum power will scale with this integer instead of - continuous. This will also ensure that the cost will scale with this integer. + normal source with the only difference being the modelling of doublets and power consumption. + The main reason for this component instead of using just a regular source is that to have the + integer behaviour of increasing the amount of doublets. In the HeatMixin an integer is created + _aggregation_count to model the amount of doublets and the maximum power will scale with this + integer instead of continuous. This will also ensure that the cost will scale with this integer. + The power consumption is computed through a COP calculation, directly linked to the heat source + production. COP is set to a default value of 0 in order to ensure power and its associated costs + are only inlcuded in the computations intentionally. Variables created: {add_variable_names_for_documentation_here} diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index 72d2e8d5..918f4552 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -8,7 +8,10 @@ @add_variables_documentation_automatically class GeothermalSourceElec(GeothermalSource): """ - TODO: document + The geothermal source electric asset is almost identical to the geothermal source one. The only + difference is that this one includes an electricity in port. The electricity power calculation that + it inherits from the geothermal source asset needs to be satisfied by an electricity carrier supplied + through this new in port. Variables created: {add_variable_names_for_documentation_here} diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index f0ce4b64..02db88be 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -423,9 +423,9 @@ def test_geothermal_source_elec(self): if __name__ == "__main__": - test_cold_demand = TestMultiCommodityHeatPump() - test_cold_demand.test_air_to_water_heat_pump_elec_min_elec() + #test_cold_demand = TestMultiCommodityHeatPump() + #test_cold_demand.test_air_to_water_heat_pump_elec_min_elec() # test_cold_demand.test_heat_pump_elec_min_elec() test_geothermal = TestGeothermalSourceElec() test_geothermal.test_geothermal_source_elec() - test_geothermal.test_geothermal_source() \ No newline at end of file + #test_geothermal.test_geothermal_source() \ No newline at end of file From 6e97d915423d392d24f8287f85824d623effa181 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 15:23:31 +0100 Subject: [PATCH 10/23] Merged main. --- CHANGELOG.md | 30 +- examples/PoCTutorial/src/run_grow_tutorial.py | 4 +- examples/municipality/src/run_municipality.py | 4 +- manual_run/testing_bugs/src/run_grow.py | 6 +- setup.py | 5 +- src/mesido/asset_sizing_mixin.py | 28 +- src/mesido/esdl/esdl_mixin.py | 41 +- src/mesido/esdl/profile_parser.py | 86 +-- src/mesido/financial_mixin.py | 127 ++--- src/mesido/workflows/gas_elect_workflow.py | 4 +- src/mesido/workflows/grow_workflow.py | 17 +- src/mesido/workflows/io/write_output.py | 9 +- .../multicommodity_simulator_workflow.py | 4 +- src/mesido/workflows/simulator_workflow.py | 4 +- .../model/heat_exchanger_with_costs.esdl | 313 +++++++++++ .../case_2a_ensemble/input/ensemble.csv | 3 + .../forecast_high_prob/timeseries_import.csv | 4 + .../forecast_low_prob/timeseries_import.csv | 4 + .../unit_cases/case_2a_ensemble/model/2a.esdl | 521 ++++++++++++++++++ .../case_2a_ensemble/output/.gitignore | 2 + .../unit_cases/case_2a_ensemble/src/run_2a.py | 183 ++++++ tests/test_end_scenario_sizing.py | 37 ++ tests/test_head_loss.py | 4 +- ...est_multiple_in_and_out_port_components.py | 7 +- tests/test_profile_parsing.py | 47 ++ tests/utils_tests.py | 14 +- tox.ini | 9 +- 27 files changed, 1347 insertions(+), 170 deletions(-) create mode 100644 tests/models/heat_exchange/model/heat_exchanger_with_costs.esdl create mode 100644 tests/models/unit_cases/case_2a_ensemble/input/ensemble.csv create mode 100644 tests/models/unit_cases/case_2a_ensemble/input/forecast_high_prob/timeseries_import.csv create mode 100644 tests/models/unit_cases/case_2a_ensemble/input/forecast_low_prob/timeseries_import.csv create mode 100644 tests/models/unit_cases/case_2a_ensemble/model/2a.esdl create mode 100644 tests/models/unit_cases/case_2a_ensemble/output/.gitignore create mode 100644 tests/models/unit_cases/case_2a_ensemble/src/run_2a.py diff --git a/CHANGELOG.md b/CHANGELOG.md index b62ac1c3..bae68ff7 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,4 +1,29 @@ -# [Unreleased-main] - 2025-12-08 +# [Unreleased-main] - 2026-02-25 + +## Added +- Parsing of ensemble profiles when using input profiles from csv. + +## Changed +- Speed-up timeseries check in from InfluxDB + +## Fixed +- xxx + + +# [0.1.17] - 2026-02-25 + +## Added +- xxx + +## Changed +- xxx + +## Fixed +- DBAccesType naming convention +- eboiler KPI + + +# [0.1.16] - 2026-02-05 ## Added - ATES variable cost calculation utilizing a split ates charging and discharging variable. @@ -22,11 +47,14 @@ - Costs of available pipe classes are updated based on the asset measures and templates if they are provided. - The charging and discharging variable for electricity storage is created without a binary variable using the convex hull description. - Gas Boiler asset is renamed as HeatSourceGas +- Improvement: Loops over timesteps are vectorized. ## Fixed - Bug: Write updated esdl for 2 port heat pump - Bug: 2 port heatpump write result profiles to database - Bug: setting of self._pipe_heat_loss_nominals was not accounting for negative values when T_ground > carrier temperature +- Bug: heat exchanger state and capacity are updated in optimized esdl file +- Additional slicing of priceprofile is required # [0.1.15] - 2025-11-19 diff --git a/examples/PoCTutorial/src/run_grow_tutorial.py b/examples/PoCTutorial/src/run_grow_tutorial.py index ce9184ef..6b02e8f4 100644 --- a/examples/PoCTutorial/src/run_grow_tutorial.py +++ b/examples/PoCTutorial/src/run_grow_tutorial.py @@ -1,6 +1,6 @@ from pathlib import Path -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_parser import ESDLFileParser from mesido.workflows import EndScenarioSizingStaged, run_end_scenario_sizing @@ -16,7 +16,7 @@ class EndScenarioSizingStagedHighs(EndScenarioSizingStaged): "write_result_db_profiles": True, "database_connections": [ { - "access_type": DBAccesType.WRITE, # DBAccesType.READ or DBAccesType.READ_WRITE + "access_type": DBAccessType.WRITE, # DBAccessType.READ or DBAccessType.READ_WRITE "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/examples/municipality/src/run_municipality.py b/examples/municipality/src/run_municipality.py index fab2efcf..dc529d66 100644 --- a/examples/municipality/src/run_municipality.py +++ b/examples/municipality/src/run_municipality.py @@ -1,6 +1,6 @@ from pathlib import Path -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_parser import ESDLFileParser from mesido.workflows import EndScenarioSizingStaged, run_end_scenario_sizing @@ -15,7 +15,7 @@ kwargs = { "database_connections": [ { - "access_type": DBAccesType.READ, # DBAccesType.WRITE or DBAccesType.READ_WRITE + "access_type": DBAccessType.READ, # DBAccessType.WRITE or DBAccessType.READ_WRITE "influxdb_host": "required_user_input", "influxdb_port": 1234, "influxdb_username": "required_user_input", diff --git a/manual_run/testing_bugs/src/run_grow.py b/manual_run/testing_bugs/src/run_grow.py index 7e92dc40..1d2f7662 100644 --- a/manual_run/testing_bugs/src/run_grow.py +++ b/manual_run/testing_bugs/src/run_grow.py @@ -2,7 +2,7 @@ from mesido.esdl.esdl_parser import ESDLFileParser from mesido.workflows import EndScenarioSizingStaged, run_end_scenario_sizing -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType if __name__ == "__main__": import time @@ -13,7 +13,7 @@ kwargs = { "database_connections": [ { - "access_type": DBAccesType.READ, # or DBAccesType.WRITE or DBAccesType.READ_WRITE + "access_type": DBAccessType.READ, # or DBAccessType.WRITE or DBAccessType.READ_WRITE "influxdb_host": "required_user_input", "influxdb_port": 1234, "influxdb_username": "required_user_input", @@ -22,7 +22,7 @@ "influxdb_verify_ssl": False, }, { - "access_type": DBAccesType.WRITE, # or DBAccesType.WRITE or DBAccesType.READ_WRITE + "access_type": DBAccessType.WRITE, # or DBAccessType.WRITE or DBAccessType.READ_WRITE "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/setup.py b/setup.py index 82a0baec..c913a76e 100644 --- a/setup.py +++ b/setup.py @@ -59,13 +59,16 @@ "pyecore >= 0.13.2", "pymoca >= 0.9.0", "rtc-tools-gil-comp == 2.6.1", + # setuptools version limitations currently: + # < 81.0.0 needed for pandapipes (still to be removed) + # < 82.0.0 needed for pkg_resources (used in rtctools) + "setuptools <= 80.9.0", "pyesdl == 25.7", "pandas >= 1.3.1, < 2.0", "casadi-gil-comp == 3.6.7", "StrEnum == 0.4.15", "CoolProp==6.6.0", ], - tests_require=["pytest", "pytest-runner", "numpy"], include_package_data=True, python_requires=">=3.10,<3.12", cmdclass=versioneer.get_cmdclass(), diff --git a/src/mesido/asset_sizing_mixin.py b/src/mesido/asset_sizing_mixin.py index 95d8c5fb..709c1504 100644 --- a/src/mesido/asset_sizing_mixin.py +++ b/src/mesido/asset_sizing_mixin.py @@ -1956,15 +1956,13 @@ def __max_size_constraints(self, ensemble_member): logger.error(f"Unexpected state: {state_val}") sys.exit(1) - for i in range(0, len(self.times())): - constraints.append( - ( - (profile_scaled[i] * max_heat_var - heat_source[i]) - / constraint_nominal, - 0.0, - np.inf, - ) + constraints.append( + ( + (profile_scaled * max_heat_var - heat_source) / constraint_nominal, + 0.0, + np.inf, ) + ) # Option 2: Normalised profile (0.0-1.0) shape that scales with maximum size of the # producer # Note: If the asset is not optional then the profile will be scaled to the @@ -1993,15 +1991,13 @@ def __max_size_constraints(self, ensemble_member): # for ProfileContraint. Future addition can be to use a different unit/quantity # etc. so that the profile is used in a normalised way and scale to max_size - for i in range(0, len(self.times())): - constraints.append( - ( - (profile_scaled[i] * max_heat - heat_source[i]) - / constraint_nominal, - 0.0, - np.inf, - ) + constraints.append( + ( + (profile_scaled * max_heat - heat_source) / constraint_nominal, + 0.0, + np.inf, ) + ) else: RuntimeError(f"{s}: Unforeseen error in adding a profile contraint") else: diff --git a/src/mesido/esdl/esdl_mixin.py b/src/mesido/esdl/esdl_mixin.py index 208a9f37..52cc77ad 100644 --- a/src/mesido/esdl/esdl_mixin.py +++ b/src/mesido/esdl/esdl_mixin.py @@ -2,6 +2,7 @@ import copy import dataclasses import logging +import os import sys from datetime import timedelta from pathlib import Path @@ -45,7 +46,7 @@ DEFAULT_END_TIMESTAMP = "2018-01-01T00:00:00+00:00" -class DBAccesType(StrEnum): +class DBAccessType(StrEnum): """ Enumeration for database access types """ @@ -71,6 +72,8 @@ class ESDLMixin( for example demand profiles. """ + csv_ensemble_mode: bool = False + esdl_run_info_path: Path = None esdl_pi_validate_timeseries: bool = False @@ -118,8 +121,8 @@ def __init__(self, *args, **kwargs) -> None: self.__energy_system_handler: esdl.esdl_handler.EnergySystemHandler = esdl_parser.get_esh() self._esdl_measures: Dict[str, Asset] = esdl_parser.get_measures() self._database_credentials: Optional[Dict[str, Tuple[str, str]]] = { - DBAccesType.READ: [], - DBAccesType.WRITE: [], + DBAccessType.READ: [], + DBAccessType.WRITE: [], } profile_reader_class = kwargs.get("profile_reader", InfluxDBProfileReader) @@ -131,7 +134,7 @@ def __init__(self, *args, **kwargs) -> None: database_connection_info = kwargs.get("database_connections", {}) read_only_dbase_credentials: Dict[str, Tuple[str, str]] = {} # for profile reader for dbconnection in database_connection_info: - if dbconnection["access_type"] != DBAccesType.WRITE: + if dbconnection["access_type"] != DBAccessType.WRITE: database_host_port = "{}:{}".format( dbconnection["influxdb_host"], dbconnection["influxdb_port"], @@ -140,7 +143,7 @@ def __init__(self, *args, **kwargs) -> None: dbconnection["influxdb_username"], dbconnection["influxdb_password"], ) - if dbconnection["access_type"] != DBAccesType.READ_WRITE: + if dbconnection["access_type"] != DBAccessType.READ_WRITE: self._database_credentials[dbconnection["access_type"]].append( { "influxdb_host": dbconnection["influxdb_host"], @@ -151,8 +154,8 @@ def __init__(self, *args, **kwargs) -> None: "influxdb_verify_ssl": dbconnection["influxdb_verify_ssl"], } ) - elif dbconnection["access_type"] == DBAccesType.READ_WRITE: - both_read_and_write = [DBAccesType.READ, DBAccesType.WRITE] + elif dbconnection["access_type"] == DBAccessType.READ_WRITE: + both_read_and_write = [DBAccessType.READ, DBAccessType.WRITE] for rw in both_read_and_write: self._database_credentials[rw].append( { @@ -167,7 +170,7 @@ def __init__(self, *args, **kwargs) -> None: else: logger.error( f"Database access type {dbconnection['access_type']} is not recognized. " - f"Please use DBAccesType.READ, DBAccesType.WRITE or DBAccesType.READ_WRITE." + f"Please use DBAccessType.READ, DBAccessType.WRITE or DBAccessType.READ_WRITE." ) sys.exit(1) @@ -691,6 +694,19 @@ def read(self) -> None: None """ super().read() + ensemble_size = 1 + self.__ensemble = None + if self.csv_ensemble_mode: + self.__ensemble = np.genfromtxt( + os.path.join(self._input_folder, "ensemble.csv"), + delimiter=",", + deletechars="", + dtype=None, + names=True, + encoding=None, + ) + ensemble_size = len(self.__ensemble) + energy_system_components = self.energy_system_components esdl_carriers = self.esdl_carriers self.hot_cold_pipe_relations() @@ -701,7 +717,8 @@ def read(self) -> None: esdl_asset_id_to_name_map=self.esdl_asset_id_to_name_map, esdl_assets=self.esdl_assets, carrier_properties=esdl_carriers, - ensemble_size=self.ensemble_size, + ensemble_size=ensemble_size, + ensemble=self.__ensemble, ) def write(self) -> None: @@ -826,3 +843,9 @@ def filter_asset_measures( filtered_assets[asset_id] = asset_type return filtered_assets + + def ensemble_member_probability(self, ensemble_member): + if self.csv_ensemble_mode: + return self.__ensemble["probability"][ensemble_member] + else: + return 1.0 diff --git a/src/mesido/esdl/profile_parser.py b/src/mesido/esdl/profile_parser.py index d8da9736..0c738508 100644 --- a/src/mesido/esdl/profile_parser.py +++ b/src/mesido/esdl/profile_parser.py @@ -58,6 +58,7 @@ def read_profiles( esdl_assets: Dict[str, Asset], carrier_properties: Dict[str, Dict], ensemble_size: int, + ensemble, ) -> None: """ This function takes a datastore and a dictionary of milp network components and loads a @@ -91,6 +92,7 @@ def read_profiles( esdl_asset_id_to_name_map=esdl_asset_id_to_name_map, carrier_properties=carrier_properties, ensemble_size=ensemble_size, + ensemble=ensemble, ) try: @@ -181,6 +183,7 @@ def _load_profiles_from_source( esdl_asset_id_to_name_map: Dict[str, str], carrier_properties: Dict[str, Dict], ensemble_size: int, + ensemble: None | np.ndarray, ) -> None: """ This function must be implemented by the child. It must load the available @@ -240,6 +243,7 @@ def _load_profiles_from_source( esdl_asset_id_to_name_map: Dict[str, str], carrier_properties: Dict[str, Dict], ensemble_size: int, + ensemble: None | np.ndarray, ) -> None: profiles: Dict[str, np.ndarray] = dict() logger.info("Reading profiles from InfluxDB") @@ -315,7 +319,6 @@ def _load_profiles_from_source( ] ) series = unique_series[index_of_unique_profile] - self._check_profile_time_series(profile_time_series=series, profile=profile) container = profile.eContainer() asset = container.eContainer() converted_dataframe = self._convert_profile_to_correct_unit( @@ -506,12 +509,15 @@ def _check_profile_time_series( # Error check: ensure that the profile data has a time resolution of 3600s (1hour) as # expected - for d1, d2 in zip(profile_time_series.index, profile_time_series.index[1:]): - if d2 - d1 != pd.Timedelta(hours=1): - raise RuntimeError( - f"The timestep for variable {profile.field} between {d1} and {d2} isn't " - f"exactly 1 hour" - ) + idx = profile_time_series.index + dt = idx.to_series().diff()[1:] + problem_timesteps = idx[:-1][dt != pd.Timedelta(hours=1)] + if not problem_timesteps.empty: + raise RuntimeError( + f"The timestep for variable {profile.field} at timestamp(s) {problem_timesteps} " + f"isn't " + f"exactly 1 hour" + ) # Check if any NaN values exist if profile_time_series.isnull().any().any(): raise Exception( @@ -614,6 +620,7 @@ def _load_profiles_from_source( esdl_asset_id_to_name_map: Dict[str, str], carrier_properties: Dict[str, Dict], ensemble_size: int, + ensemble: None | np.ndarray, ) -> None: if self._file_path.suffix == ".xml": logger.warning( @@ -628,6 +635,7 @@ def _load_profiles_from_source( energy_system_components=energy_system_components, carrier_properties=carrier_properties, ensemble_size=ensemble_size, + ensemble=ensemble, ) else: raise _ProfileParserException( @@ -639,25 +647,27 @@ def _load_csv( energy_system_components: Dict[str, Set[str]], carrier_properties: Dict[str, Dict], ensemble_size: int, + ensemble, ) -> None: - data = pd.read_csv(self._file_path) - if len(data.filter(like="Unnamed").columns) > 0: - raise Exception( - f"An unnamed column has been found in profile source file: {self._file_path}" - ) + data_dict = {} + if ensemble_size > 1: + input_folder = self._file_path.parent + file_name = self._file_path.name + for i, ensemble_name, _ in ensemble: + data_dict[i] = pd.read_csv(Path(input_folder / ensemble_name / file_name)) + else: + data_dict[0] = pd.read_csv(self._file_path) - try: - timeseries_import_times = [ - datetime.datetime.strptime(entry.replace("Z", ""), "%Y-%m-%d %H:%M:%S").replace( - tzinfo=datetime.timezone.utc + for data in data_dict.values(): + if len(data.filter(like="Unnamed").columns) > 0: + raise Exception( + f"An unnamed column has been found in profile source file: {self._file_path}" ) - for entry in data["DateTime"].to_numpy() - ] - except ValueError: + try: timeseries_import_times = [ - datetime.datetime.strptime(entry.replace("Z", ""), "%Y-%m-%dT%H:%M:%S").replace( + datetime.datetime.strptime(entry.replace("Z", ""), "%Y-%m-%d %H:%M:%S").replace( tzinfo=datetime.timezone.utc ) for entry in data["DateTime"].to_numpy() @@ -666,48 +676,58 @@ def _load_csv( try: timeseries_import_times = [ datetime.datetime.strptime( - entry.replace("Z", ""), "%d-%m-%Y %H:%M" + entry.replace("Z", ""), "%Y-%m-%dT%H:%M:%S" ).replace(tzinfo=datetime.timezone.utc) for entry in data["DateTime"].to_numpy() ] except ValueError: - raise _ProfileParserException("Date time string is not in a supported format") + try: + timeseries_import_times = [ + datetime.datetime.strptime( + entry.replace("Z", ""), "%d-%m-%Y %H:%M" + ).replace(tzinfo=datetime.timezone.utc) + for entry in data["DateTime"].to_numpy() + ] + except ValueError: + raise _ProfileParserException( + "Date time string is not in a supported format" + ) - logger.warning("Timezone specification not supported yet: default UTC has been used") + logger.warning("Timezone specification not supported yet: default UTC has been used") - self._reference_datetimes = timeseries_import_times + self._reference_datetimes = timeseries_import_times - for ensemble_member in range(ensemble_size): + for e_m in range(ensemble_size): + data_em = data_dict[e_m] for component_type, var_name in self.component_type_to_var_name_map.items(): for component_name in energy_system_components.get(component_type, []): try: column_name = f"{component_name.replace(' ', '')}" - values = data[column_name].to_numpy() + values = data_em[column_name].to_numpy() if np.isnan(values).any(): raise Exception( f"Column name: {column_name}, NaN exists in the profile source" f" file {self._file_path}." - f" Detials: {data[data[column_name].isnull()]}" + f" Detials: {data_em[data_em[column_name].isnull()]}" ) except KeyError: pass else: - self._profiles[ensemble_member][component_name + var_name] = values + self._profiles[e_m][component_name + var_name] = values for properties in carrier_properties.values(): carrier_name = properties.get("name") try: - values = data[carrier_name].to_numpy() + values = data_em[carrier_name].to_numpy() if np.isnan(values).any(): raise Exception( f"Carrier name: {carrier_name}, NaN exists in the profile source file" - f" {self._file_path}. Details: {data[data[carrier_name].isnull()]}" + f" {self._file_path}. Details: " + f"{data_em[data_em[carrier_name].isnull()]}" ) except KeyError: pass else: - self._profiles[ensemble_member][ - carrier_name + self.carrier_profile_var_name - ] = values + self._profiles[e_m][carrier_name + self.carrier_profile_var_name] = values def _load_xml(self, energy_system_components, esdl_asset_id_to_name_map): timeseries_import_basename = self._file_path.stem diff --git a/src/mesido/financial_mixin.py b/src/mesido/financial_mixin.py index 9a3368d7..09783d79 100644 --- a/src/mesido/financial_mixin.py +++ b/src/mesido/financial_mixin.py @@ -975,16 +975,14 @@ def __variable_operational_cost_constraints(self, ensemble_member): else: price_profile = Timeseries(self.times(), np.zeros(len(self.times()))) - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient - * (heat_charge[i] + heat_discharge[i]) - * timesteps[i - 1] - ) - sum += price_profile.values[i] * pump_power[i] * timesteps[i - 1] / eff + sum_ = ca.sum1( + variable_operational_cost_coefficient + * (heat_charge[1:] + heat_discharge[1:]) + * timesteps + ) + sum_ += ca.sum1(price_profile.values[1:] * pump_power[1:] * timesteps / eff) - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for asset in [ *self.energy_system_components.get("pump", []), @@ -1011,11 +1009,9 @@ def __variable_operational_cost_constraints(self, ensemble_member): else: price_profile = Timeseries(self.times(), np.zeros(len(self.times()))) - sum = 0.0 - for i in range(1, len(self.times())): - sum += price_profile.values[i] * pump_power[i] * timesteps_hr[i - 1] / eff + sum_ = ca.sum1(price_profile.values[1:] * pump_power[1:] * timesteps_hr / eff) - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for s in self.energy_system_components.get("heat_source", []): heat_source = self.__state_vector_scaled(f"{s}.Heat_source", ensemble_member) @@ -1072,17 +1068,14 @@ def __variable_operational_cost_constraints(self, ensemble_member): ) # [W] else: nominator_vector = heat_source - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient - * nominator_vector[i] - * timesteps_hr[i - 1] - / denominator - ) - sum += price_profile.values[i] * pump_power[i] * timesteps_hr[i - 1] / eff - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + sum_ = ( + ca.sum1(variable_operational_cost_coefficient * nominator_vector[1:] * timesteps_hr) + / denominator + ) + sum_ += ca.sum1(price_profile.values[1:] * pump_power[1:] * timesteps_hr / eff) + + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for hp in [ *self.energy_system_components.get("heat_pump", []), @@ -1111,20 +1104,15 @@ def __variable_operational_cost_constraints(self, ensemble_member): else: price_profile = Timeseries(self.times(), np.zeros(len(self.times()))) - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient - * elec_consumption[i] - * timesteps_hr[i - 1] - ) - sum += price_profile.values[i] * pump_power[i] * timesteps_hr[i - 1] / eff - if hp not in self.energy_system_components.get("heat_pump_elec", []): - # assuming that if heatpump has electricity port, the cost for the electricity - # are already made by the electricity producer and transport - sum += price_profile.values[i] * elec_consumption[i] * timesteps_hr[i - 1] - - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + sum_ = ca.sum1( + variable_operational_cost_coefficient * elec_consumption[1:] * timesteps_hr + ) + sum_ += ca.sum1(price_profile.values[1:] * pump_power[1:] * timesteps_hr / eff) + if hp not in self.energy_system_components.get("heat_pump_elec", []): + # assuming that if heatpump has electricity port, the cost for the electricity + # are already made by the electricity producer and transport + sum_ += ca.sum1(price_profile.values[1:] * elec_consumption[1:] * timesteps_hr) + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for ac in self.energy_system_components.get("airco", []): heat_airco = self.__state_vector_scaled(f"{ac}.Heat_airco", ensemble_member) @@ -1136,11 +1124,9 @@ def __variable_operational_cost_constraints(self, ensemble_member): variable_operational_cost_coefficient = parameters[ f"{ac}.variable_operational_cost_coefficient" ] - sum = 0.0 - for i in range(1, len(self.times())): - sum += variable_operational_cost_coefficient * heat_airco[i] * timesteps_hr[i - 1] + sum_ = ca.sum1(variable_operational_cost_coefficient * heat_airco[1:] * timesteps_hr) - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for demand in self.energy_system_components.get("gas_demand", []): gas_mass_flow = self.__state_vector_scaled( @@ -1156,13 +1142,8 @@ def __variable_operational_cost_constraints(self, ensemble_member): f"{demand}.variable_operational_cost_coefficient" ] - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient * gas_mass_flow[i] * timesteps_hr[i - 1] - ) - - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + sum_ = ca.sum1(variable_operational_cost_coefficient * gas_mass_flow[1:] * timesteps_hr) + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for gs in self.energy_system_components.get("gas_source", []): gas_produced_g_s = self.__state_vector_scaled( @@ -1177,14 +1158,11 @@ def __variable_operational_cost_constraints(self, ensemble_member): f"{gs}.variable_operational_cost_coefficient" ] timesteps_sec = np.diff(self.times()) - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient - * gas_produced_g_s[i] - * timesteps_sec[i - 1] - ) # [euro/g] * [g/s] * [s] - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + sum_ = ca.sum1( + variable_operational_cost_coefficient * gas_produced_g_s[1:] * timesteps_sec + ) + # [euro/g] * [g/s] * [s] + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) for es in self.energy_system_components.get("electricity_source", []): elec_produced_w = self.__state_vector_scaled( @@ -1198,12 +1176,10 @@ def __variable_operational_cost_constraints(self, ensemble_member): variable_operational_cost_coefficient = parameters[ # euro / Wh f"{es}.variable_operational_cost_coefficient" ] - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient * elec_produced_w[i] * timesteps_hr[i - 1] - ) # [euro/Wh] * [W] * [hr] - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + sum_ = ca.sum1( + variable_operational_cost_coefficient * elec_produced_w[1:] * timesteps_hr + ) # [euro/Wh] * [W] * [hr] + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) # for a in self.heat_network_components.get("ates", []): # TODO: needs to be replaced with the positive or abs value of this, see varOPEX, @@ -1242,15 +1218,13 @@ def __variable_operational_cost_constraints(self, ensemble_member): f"{electrolyzer}.variable_operational_cost_coefficient" ] - sum = 0.0 - for i in range(1, len(self.times())): - sum += ( - variable_operational_cost_coefficient - * power_consumer[i] - * timesteps_hr[i - 1] # gas_mass_flow unit is g/s - ) + sum_ = ca.sum1( + variable_operational_cost_coefficient + * power_consumer[1:] + * timesteps_hr # gas_mass_flow unit is g/s + ) - constraints.append(((variable_operational_cost - sum) / nominal, 0.0, 0.0)) + constraints.append(((variable_operational_cost - sum_) / nominal, 0.0, 0.0)) # for a in self.heat_network_components.get("ates", []): # TODO: needs to be replaced with the positive or abs value of this, see varOPEX, @@ -1636,7 +1610,12 @@ def __revenue_constraints(self, ensemble_member): carrier_name = attr["name"] cost_multiplier = 1.0 # priceprofile gas is in EUR/g if carrier_name is not None: - price_profile = self.get_timeseries(f"{carrier_name}.price_profile").values + price_profile_timeseries = self.get_timeseries(f"{carrier_name}.price_profile") + # The slicing is required if the timeseries wasn't adapted in the read + mask = (price_profile_timeseries.times >= self.times()[0]) & ( + price_profile_timeseries.times <= self.times()[-1] + ) + price_profile = price_profile_timeseries.values[mask] if demand in self.energy_system_components.get("gas_demand", []): energy_flow = self.__state_vector_scaled( @@ -1652,12 +1631,10 @@ def __revenue_constraints(self, ensemble_member): variable_revenue = self.extra_variable(variable_revenue_var, ensemble_member) nominal = self.variable_nominal(variable_revenue_var) - sum = 0.0 timesteps = np.diff(self.times()) * cost_multiplier - for i in range(1, len(self.times())): - sum += price_profile[i] * energy_flow[i] * timesteps[i - 1] + sum_ = ca.sum1(price_profile[1:] * energy_flow[1:] * timesteps) - constraints.append(((variable_revenue - sum) / (nominal), 0.0, 0.0)) + constraints.append(((variable_revenue - sum_) / (nominal), 0.0, 0.0)) return constraints diff --git a/src/mesido/workflows/gas_elect_workflow.py b/src/mesido/workflows/gas_elect_workflow.py index eeda5d9e..f05ad27a 100644 --- a/src/mesido/workflows/gas_elect_workflow.py +++ b/src/mesido/workflows/gas_elect_workflow.py @@ -2,7 +2,7 @@ import os from mesido.esdl.esdl_additional_vars_mixin import ESDLAdditionalVarsMixin -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_mixin import ESDLMixin from mesido.head_loss_class import HeadLossOption from mesido.network_common import NetworkSettings @@ -200,7 +200,7 @@ def main(runinfo_path, log_level): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/src/mesido/workflows/grow_workflow.py b/src/mesido/workflows/grow_workflow.py index f8bc19ee..23017677 100644 --- a/src/mesido/workflows/grow_workflow.py +++ b/src/mesido/workflows/grow_workflow.py @@ -6,7 +6,7 @@ from typing import Dict from mesido.esdl.esdl_additional_vars_mixin import ESDLAdditionalVarsMixin -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_mixin import ESDLMixin from mesido.head_loss_class import HeadLossOption from mesido.potential_errors import reset_potential_errors @@ -375,9 +375,16 @@ def constraints(self, ensemble_member): vars = self.state_vector(f"{b}.Heat_buffer") symbol_stored_heat = self.state_vector(f"{b}.Stored_heat") constraints.append((symbol_stored_heat[self.__indx_max_peak], 0.0, 0.0)) - for i in range(len(self.times())): - if i < self.__indx_max_peak or i > (self.__indx_max_peak + 23): - constraints.append((vars[i], 0.0, 0.0)) + + ind_peak = int(self.__indx_max_peak) + constraints.append((vars[:ind_peak], 0.0, 0.0)) + constraints.append( + ( + vars[ind_peak + 24 :], + 0.0, + 0.0, + ) + ) # TODO: confirm if volume or heat balance is required over year. This will # determine if cyclic contraint below is for stored_heat or stored_volume @@ -864,7 +871,7 @@ def main(runinfo_path, log_level): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/src/mesido/workflows/io/write_output.py b/src/mesido/workflows/io/write_output.py index 0261b82c..488d3976 100644 --- a/src/mesido/workflows/io/write_output.py +++ b/src/mesido/workflows/io/write_output.py @@ -18,7 +18,7 @@ import mesido.esdl.esdl_parser from mesido.constants import GRAVITATIONAL_CONSTANT from mesido.esdl.edr_pipe_class import EDRPipeClass -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.financial_mixin import calculate_annuity_factor from mesido.network_common import NetworkSettings from mesido.post_processing.post_processing_utils import pipe_pressure, pipe_velocity @@ -60,7 +60,7 @@ def __init__(self, **kwargs): sys.exit(1) if self.write_result_db_profiles: - database_connection_write = self._database_credentials.get(DBAccesType.WRITE, []) + database_connection_write = self._database_credentials.get(DBAccessType.WRITE, []) if len(database_connection_write) == 0: logger.error("The connections settings for writing to a database is empty") @@ -466,8 +466,8 @@ def _add_kpis_to_energy_system(self, energy_system, optimizer_sim: bool = False) or asset.asset_type == "GenericProducer" or asset.asset_type == "ResidualHeatSource" or asset.asset_type == "GeothermalSource" - or asset.asset_type == "ResidualHeatSource" or asset.asset_type == "GasHeater" + or asset.asset_type == "ElectricBoiler" ): heat_source_energy_wh[asset.name] = np.sum( results[f"{asset.name}.Heat_source"][1:] * diff_times / 3600 @@ -1031,6 +1031,7 @@ def _write_updated_esdl( *self.energy_system_components.get("heat_buffer", []), *self.energy_system_components.get("heat_pump", []), *self.energy_system_components.get("airco", []), + *self.energy_system_components.get("heat_exchanger", []), ]: asset = self._name_to_asset(energy_system, name) asset_placement_var = self._asset_aggregation_count_var_map[name] @@ -1050,6 +1051,8 @@ def _write_updated_esdl( * parameters[f"{name}.rho"] * parameters[f"{name}.dT"] ) + elif asset.name in self.energy_system_components.get("heat_exchanger", []): + asset.capacity = max_size elif asset.name in [ *self.energy_system_components.get("heat_pump", []), *self.energy_system_components.get("airco", []), diff --git a/src/mesido/workflows/multicommodity_simulator_workflow.py b/src/mesido/workflows/multicommodity_simulator_workflow.py index 21e27464..a078625d 100644 --- a/src/mesido/workflows/multicommodity_simulator_workflow.py +++ b/src/mesido/workflows/multicommodity_simulator_workflow.py @@ -5,7 +5,7 @@ import esdl -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_mixin import ESDLMixin from mesido.esdl.esdl_parser import ESDLFileParser from mesido.esdl.profile_parser import ProfileReaderFromFile @@ -833,7 +833,7 @@ def main(runinfo_path, log_level): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/src/mesido/workflows/simulator_workflow.py b/src/mesido/workflows/simulator_workflow.py index 934aee2f..4b9f3274 100644 --- a/src/mesido/workflows/simulator_workflow.py +++ b/src/mesido/workflows/simulator_workflow.py @@ -3,7 +3,7 @@ import esdl -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_mixin import ESDLMixin from mesido.head_loss_class import HeadLossOption from mesido.techno_economic_mixin import TechnoEconomicMixin @@ -359,7 +359,7 @@ def main(runinfo_path, log_level): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/tests/models/heat_exchange/model/heat_exchanger_with_costs.esdl b/tests/models/heat_exchange/model/heat_exchanger_with_costs.esdl new file mode 100644 index 00000000..8148e8ea --- /dev/null +++ b/tests/models/heat_exchange/model/heat_exchanger_with_costs.esdl @@ -0,0 +1,313 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/unit_cases/case_2a_ensemble/input/ensemble.csv b/tests/models/unit_cases/case_2a_ensemble/input/ensemble.csv new file mode 100644 index 00000000..c26cb83e --- /dev/null +++ b/tests/models/unit_cases/case_2a_ensemble/input/ensemble.csv @@ -0,0 +1,3 @@ +Unnamed: 0,name,probability +0,forecast_high_prob,0.7 +1,forecast_low_prob,0.3 diff --git a/tests/models/unit_cases/case_2a_ensemble/input/forecast_high_prob/timeseries_import.csv b/tests/models/unit_cases/case_2a_ensemble/input/forecast_high_prob/timeseries_import.csv new file mode 100644 index 00000000..1e9ee86c --- /dev/null +++ b/tests/models/unit_cases/case_2a_ensemble/input/forecast_high_prob/timeseries_import.csv @@ -0,0 +1,4 @@ +DateTime,HeatingDemand_7484,HeatingDemand_c6c8,HeatingDemand_6f99 +2013-05-19 22:00:00,350000,350000,350000 +2013-05-19 23:00:00,350000,350000,350000 +2013-05-20 00:00:00,350000,350000,350000 diff --git a/tests/models/unit_cases/case_2a_ensemble/input/forecast_low_prob/timeseries_import.csv b/tests/models/unit_cases/case_2a_ensemble/input/forecast_low_prob/timeseries_import.csv new file mode 100644 index 00000000..74e7a9e1 --- /dev/null +++ b/tests/models/unit_cases/case_2a_ensemble/input/forecast_low_prob/timeseries_import.csv @@ -0,0 +1,4 @@ +DateTime,HeatingDemand_7484,HeatingDemand_c6c8,HeatingDemand_6f99 +2013-05-19 22:00:00,350000,350000,300000 +2013-05-19 23:00:00,350000,350000,300000 +2013-05-20 00:00:00,350000,350000,300000 diff --git a/tests/models/unit_cases/case_2a_ensemble/model/2a.esdl b/tests/models/unit_cases/case_2a_ensemble/model/2a.esdl new file mode 100644 index 00000000..1410c252 --- /dev/null +++ b/tests/models/unit_cases/case_2a_ensemble/model/2a.esdl @@ -0,0 +1,521 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/unit_cases/case_2a_ensemble/output/.gitignore b/tests/models/unit_cases/case_2a_ensemble/output/.gitignore new file mode 100644 index 00000000..c4697e99 --- /dev/null +++ b/tests/models/unit_cases/case_2a_ensemble/output/.gitignore @@ -0,0 +1,2 @@ +diag.xml +timeseries_export.xml diff --git a/tests/models/unit_cases/case_2a_ensemble/src/run_2a.py b/tests/models/unit_cases/case_2a_ensemble/src/run_2a.py new file mode 100644 index 00000000..ed69bd82 --- /dev/null +++ b/tests/models/unit_cases/case_2a_ensemble/src/run_2a.py @@ -0,0 +1,183 @@ +from mesido.asset_sizing_mixin import AssetSizingMixin +from mesido.esdl.esdl_mixin import ESDLMixin +from mesido.esdl.esdl_parser import ESDLFileParser +from mesido.esdl.profile_parser import ProfileReaderFromFile +from mesido.physics_mixin import PhysicsMixin + +import numpy as np + +from rtctools.optimization.collocated_integrated_optimization_problem import ( + CollocatedIntegratedOptimizationProblem, +) +from rtctools.optimization.control_tree_mixin import ControlTreeMixin +from rtctools.optimization.goal_programming_mixin import Goal, GoalProgrammingMixin +from rtctools.optimization.linearized_order_goal_programming_mixin import ( + LinearizedOrderGoalProgrammingMixin, +) +from rtctools.util import run_optimization_problem + + +class TargetDemandGoal(Goal): + priority = 1 + + order = 1 + + def __init__(self, state, target, index): + self.state = state + self.index = index + + self.targets = target + self.function_nominal = np.median(target[0].values[index]) + + def function(self, optimization_problem, ensemble_member): + """Note that for path goals, the ensemble member index is not passed to the call + to :func:`OptimizationProblem.state`. This call returns a time-independent symbol + that is also independent of the active ensemble member. Path goals are + applied to all times and all ensemble members simultaneously.""" + nom = optimization_problem.variable_nominal(self.state) + vector_state = ( + nom * optimization_problem.state_vector(self.state, ensemble_member)[self.index] + ) + + return self.targets[ensemble_member].values[self.index] - vector_state + + +class TargetDemandPathGoal(Goal): + priority = 1 + + order = 2 + + def __init__(self, state, target): + self.state = state + + self.target_min = target + self.target_max = target + self.function_range = (0.0, 2.0 * max(target.values)) + self.function_nominal = np.median(target.values) + + def function(self, optimization_problem, ensemble_member): + """Path goals are applied to all timeseries and ensemble members simultaneously, + therefore the goals ensemble independent and only one target for all the ensembles can be + applied.""" + + return optimization_problem.state(self.state) + + +class MinimizeSourcesSizeGoal(Goal): + priority = 2 + + order = 2 + + def __init__(self, source): + self.target_max = 0.0 + self.function_range = (0.0, 10e6) + self.source = source + self.function_nominal = 1e5 + + def function(self, optimization_problem, ensemble_member): + return optimization_problem.extra_variable(f"{self.source}__max_size", ensemble_member) + + +class _GoalsAndOptions: + + def goals(self): + goals = super().goals().copy() + + for i in range(len(self.times())): + for demand in self.energy_system_components["heat_demand"]: + target = [ + self.get_timeseries( + f"{demand}.target_heat_demand", ensemble_member=ensemble_member + ) + for ensemble_member in range(self.ensemble_size) + ] + state = f"{demand}.Heat_demand" + + goals.append(TargetDemandGoal(state, target, i)) + + return goals + + +class HeatProblem( + _GoalsAndOptions, + PhysicsMixin, + LinearizedOrderGoalProgrammingMixin, + GoalProgrammingMixin, + ESDLMixin, + CollocatedIntegratedOptimizationProblem, +): + + def energy_system_options(self): + options = super().energy_system_options() + options["heat_loss_disconnected_pipe"] = True + self.heat_network_settings["minimum_velocity"] = 0.0001 + + return options + + def solver_options(self): + options = super().solver_options() + options["solver"] = "highs" + return options + + +class HeatProblemEnsemble( + AssetSizingMixin, + HeatProblem, + ControlTreeMixin, +): + csv_ensemble_mode = True + + def goals(self): + goals = super().goals().copy() + for s in self.energy_system_components["heat_source"]: + goals.append(MinimizeSourcesSizeGoal(s)) + + return goals + + def __fixed_max_size(self): + constraints = [] + for heat_source in self.energy_system_components["heat_source"]: + max_size_prev = self.extra_variable(f"{heat_source}__max_size", ensemble_member=0) + for e_m in range(self.ensemble_size - 1): + max_size = self.extra_variable(f"{heat_source}__max_size", ensemble_member=e_m + 1) + constraints.append((max_size - max_size_prev, 0.0, 0.0)) + max_size_prev = max_size + return constraints + + def __update_target_demand_constraint(self, ensemble_member): + """ + This constraint method adds upper limits to the demand heat production based on the + target of a specific ensemble. This requires the first ensemble to always have the + timeseries with the largest values. + """ + constraints = [] + for demand in self.energy_system_components["heat_demand"]: + var = self.state_vector(f"{demand}.Heat_demand", ensemble_member=ensemble_member) + nom = self.variable_nominal(f"{demand}.Heat_demand") + target = np.asarray( + self.get_timeseries(f"{demand}.target_heat_demand", ensemble_member).values + ) + + constraints.append((var * nom / nom, 0.0, target / nom)) + return constraints + + def constraints(self, ensemble_member): + constraints = super().constraints(ensemble_member).copy() + constraints.extend(self.__fixed_max_size()) + constraints.extend(self.__update_target_demand_constraint(ensemble_member)) + return constraints + + def path_constraints(self, ensemble_member): + constraints = super().path_constraints(ensemble_member).copy() + + return constraints + + +if __name__ == "__main__": + run_optimization_problem( + HeatProblem, + esdl_file_name="2a.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_import.xml", + ) diff --git a/tests/test_end_scenario_sizing.py b/tests/test_end_scenario_sizing.py index 1799e360..719e8d6a 100644 --- a/tests/test_end_scenario_sizing.py +++ b/tests/test_end_scenario_sizing.py @@ -24,6 +24,43 @@ class TestEndScenarioSizing(TestCase): + def test_heat_exchanger_sizing(self): + """ + Check heat exchanger can be sized in EndScenarioSizingStaged problem. + After optimization asset state and capacity attributes are changed. + + Checks: + - max_size variable of the asset is calculated + - heat exchanger state attribute is changed + - heat exchanger capacity attribute is updated + """ + import models.heat_exchange.src.run_heat_exchanger as run_heat_exchanger + + base_folder = Path(run_heat_exchanger.__file__).resolve().parent.parent + + solution = run_end_scenario_sizing( + EndScenarioSizingStaged, + base_folder=base_folder, + esdl_file_name="heat_exchanger_with_costs.esdl", + esdl_parser=ESDLFileParser, + ) + results = solution.extract_results() + + # Check heat exchanger is sized + np.testing.assert_allclose( + max(results["HeatExchange_39ed.Secondary_heat"]), results["HeatExchange_39ed__max_size"] + ) + + # Check heat exchanger state attribute is changed from OPTIONAL + # to ENABLED after the optimization + energy_system = solution._ESDLMixin__energy_system_handler.energy_system + asset = solution._name_to_asset(energy_system, "HeatExchange_39ed") + np.testing.assert_equal(esdl.AssetStateEnum.ENABLED, asset.state) + + # Check heat exchanger capacity attribute is updated + # with max_size variable after the optimization + np.testing.assert_allclose(results["HeatExchange_39ed__max_size"], asset.capacity) + @classmethod def setUpClass(cls) -> None: import models.test_case_small_network_ates_buffer_optional_assets.src.run_ates as run_ates diff --git a/tests/test_head_loss.py b/tests/test_head_loss.py index a1af43c4..e54fd7f0 100644 --- a/tests/test_head_loss.py +++ b/tests/test_head_loss.py @@ -3,7 +3,7 @@ import mesido._darcy_weisbach as darcy_weisbach from mesido.constants import GRAVITATIONAL_CONSTANT -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_parser import ESDLFileParser from mesido.esdl.profile_parser import ProfileReaderFromFile from mesido.head_loss_class import HeadLossOption @@ -84,7 +84,7 @@ def energy_system_options(self): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/tests/test_multiple_in_and_out_port_components.py b/tests/test_multiple_in_and_out_port_components.py index 78d46e31..f3dc1a13 100644 --- a/tests/test_multiple_in_and_out_port_components.py +++ b/tests/test_multiple_in_and_out_port_components.py @@ -1,7 +1,7 @@ from pathlib import Path from unittest import TestCase -from mesido.esdl.esdl_mixin import DBAccesType +from mesido.esdl.esdl_mixin import DBAccessType from mesido.esdl.esdl_parser import ESDLFileParser from mesido.esdl.profile_parser import ProfileReaderFromFile from mesido.util import run_esdl_mesido_optimization @@ -12,6 +12,7 @@ class TestHEX(TestCase): + def test_heat_exchanger(self): """ Check the modelling of the heat exchanger component which allows two hydraulically @@ -54,7 +55,7 @@ def energy_system_options(self): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, @@ -352,7 +353,7 @@ def energy_system_options(self): "write_result_db_profiles": False, "database_connections": [ { - "access_type": DBAccesType.WRITE, + "access_type": DBAccessType.WRITE, "influxdb_host": "localhost", "influxdb_port": 8086, "influxdb_username": None, diff --git a/tests/test_profile_parsing.py b/tests/test_profile_parsing.py index fd743918..0184657c 100644 --- a/tests/test_profile_parsing.py +++ b/tests/test_profile_parsing.py @@ -315,6 +315,53 @@ def test_loading_profile_from_esdl(self): atol=1e-2, ) + def test_loading_profiles_ensemble_members(self): + """ + This test constructs multiple ensemble members based on an "ensemble_member" CSV file + that describes the probability of the ensemble member and the name and number. + The profiles related to each ensemble member are read from the respective CSV files and + saved in for each member. + The test checks if the profiles read match the profiles from the CVS files and if the + ensemble_member_size is set accordingly. + """ + import models.unit_cases.case_2a_ensemble.src.run_2a as run_2a + from models.unit_cases.case_2a_ensemble.src.run_2a import HeatProblemEnsemble + + base_folder = Path(run_2a.__file__).resolve().parent.parent + model_folder = base_folder / "model" + input_folder = base_folder / "input" + + problem = HeatProblemEnsemble( + base_folder=base_folder, + model_folder=model_folder, + input_folder=input_folder, + esdl_file_name="2a.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_import.csv", + ) + + problem.pre() + + # check that the ensemble size is set at 2, which is based on the ensemble.csv + np.testing.assert_equal(problem.ensemble_size, 2) + prob_0 = problem.ensemble_member_probability(0) + prob_1 = problem.ensemble_member_probability(1) + np.testing.assert_allclose(prob_0, 0.7) + np.testing.assert_allclose(prob_1, 0.3) + + # check that the timeseries are loaded for all ensemble sizes and that the timeseries are + # equal to a heating demand of 350000 except for HeatingDemand_6f99 at the second + # ensemble, where it is equal to 300000. + timeseries_names = problem.io.get_timeseries_names() + for t_name in timeseries_names: + for e_m in range(problem.ensemble_size): + t_series = problem.get_timeseries(t_name, e_m) + if e_m == 1 and t_name == "HeatingDemand_6f99.target_heat_demand": + np.testing.assert_allclose(t_series.values, [300000] * 3) + else: + np.testing.assert_allclose(t_series.values, [350000] * 3) + if __name__ == "__main__": # unittest.main() diff --git a/tests/utils_tests.py b/tests/utils_tests.py index bd9d6d7b..d50ba5b6 100644 --- a/tests/utils_tests.py +++ b/tests/utils_tests.py @@ -21,21 +21,25 @@ def feasibility_test(solution): ) -def demand_matching_test(solution, results, atol=1.0e-3, rtol=1.0e-6): +def demand_matching_test(solution, results, ensemble_member=0, atol=1.0e-3, rtol=1.0e-6): """ "Test function to check whether the milp demand of each consumer is matched""" for d in solution.energy_system_components.get("heat_demand", []): if len(solution.times()) > 0: len_times = len(solution.times()) else: len_times = len(solution.get_timeseries(f"{d}.target_heat_demand").values) - target = solution.get_timeseries(f"{d}.target_heat_demand").values[0:len_times] + target = solution.get_timeseries(f"{d}.target_heat_demand", ensemble_member).values[ + 0:len_times + ] np.testing.assert_allclose(target, results[f"{d}.Heat_demand"], atol=atol, rtol=rtol) for d in solution.energy_system_components.get("cold_demand", []): if len(solution.times()) > 0: len_times = len(solution.times()) else: len_times = len(solution.get_timeseries(f"{d}.target_cold_demand").values) - target = solution.get_timeseries(f"{d}.target_cold_demand").values[0:len_times] + target = solution.get_timeseries(f"{d}.target_cold_demand", ensemble_member).values[ + 0:len_times + ] np.testing.assert_allclose(target, results[f"{d}.Cold_demand"], atol=atol, rtol=rtol) for d in solution.energy_system_components.get("gas_demand", []): timeseries_name = f"{d}.target_gas_demand" @@ -44,7 +48,7 @@ def demand_matching_test(solution, results, atol=1.0e-3, rtol=1.0e-6): len_times = len(solution.times()) else: len_times = len(solution.get_timeseries(timeseries_name).values) - target = solution.get_timeseries(timeseries_name).values[0:len_times] + target = solution.get_timeseries(timeseries_name, ensemble_member).values[0:len_times] np.testing.assert_allclose( target, results[f"{d}.Gas_demand_mass_flow"], atol=atol, rtol=rtol ) @@ -55,7 +59,7 @@ def demand_matching_test(solution, results, atol=1.0e-3, rtol=1.0e-6): len_times = len(solution.times()) else: len_times = len(solution.get_timeseries(timeseries_name).values) - target = solution.get_timeseries(timeseries_name).values[0:len_times] + target = solution.get_timeseries(timeseries_name, ensemble_member).values[0:len_times] np.testing.assert_allclose( target, results[f"{d}.Electricity_demand"], atol=atol, rtol=rtol ) diff --git a/tox.ini b/tox.ini index 2df141d2..98cf1930 100644 --- a/tox.ini +++ b/tox.ini @@ -15,10 +15,10 @@ deps = pytest-ordering pytest-timeout numpy - pandapipes==0.10.0 - pandapower==2.14.6 - deepdiff==7.0.1 -extras = all + pandapipes == 0.10.0 + pandapower == 2.14.6 + deepdiff == 7.0.1 + pip == 25.3 # Main tests part 1 (typical normal test) [testenv:test_env_main_1] @@ -43,6 +43,7 @@ deps = flake8-comprehensions == 3.15.0 flake8-import-order == 0.18.2 pep8-naming == 0.13.3 + setuptools == 80.9.0 commands = flake8 examples src tests setup.py --exclude=pandapipeesdlparser.py From cf80ec6a1c6079a8f9a48bd0a969921f393a18cb Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 15:51:03 +0100 Subject: [PATCH 11/23] Changed default COP value for geothermal source. --- src/mesido/esdl/esdl_heat_model.py | 2 +- .../models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 01726223..8ef5c4e1 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1315,7 +1315,7 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS ) max_supply=asset.attributes["power"] modifiers["elec_power_nominal"] = max_supply - modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 1.0 + modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 10.0 # Check to see if there is an electricity carrier at the in ports. in_port_carriers = [port.carrier for port in asset.in_ports] if any([isinstance(carrier, esdl.esdl.ElectricityCommodity) for carrier in in_port_carriers]): diff --git a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl index c437cc6d..c5afad36 100644 --- a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl +++ b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec.esdl @@ -52,7 +52,7 @@ - + From 2a163879c7155fcb9ed8c41570a3cd111d468e21 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 16:07:12 +0100 Subject: [PATCH 12/23] Formatting. --- src/mesido/esdl/esdl_heat_model.py | 11 ++++-- src/mesido/esdl/esdl_model_base.py | 38 +++++++++---------- src/mesido/financial_mixin.py | 6 ++- .../milp/heat/geothermal_source.py | 10 ++--- .../milp/heat/geothermal_source_elec.py | 2 - tests/test_multicommodity.py | 13 ++++--- 6 files changed, 43 insertions(+), 37 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 8ef5c4e1..6345f2dd 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1313,12 +1313,17 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS f"{asset.asset_type} '{asset.name}' has no desired flow rate specified. " f"'{asset.name}' will not be actuated in a constant manner" ) - max_supply=asset.attributes["power"] + max_supply = asset.attributes["power"] modifiers["elec_power_nominal"] = max_supply modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 10.0 # Check to see if there is an electricity carrier at the in ports. in_port_carriers = [port.carrier for port in asset.in_ports] - if any([isinstance(carrier, esdl.esdl.ElectricityCommodity) for carrier in in_port_carriers]): + if any( + [ + isinstance(carrier, esdl.esdl.ElectricityCommodity) + for carrier in in_port_carriers + ] + ): for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): min_voltage = port.carrier.voltage @@ -1329,7 +1334,7 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS Power=dict(min=0.0, max=max_supply, nominal=max_supply / 2.0), I=dict(min=0.0, max=i_max, nominal=i_nom), V=dict(min=min_voltage, nominal=min_voltage), - ) + ), ) return GeothermalSourceElec, modifiers else: diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index a09c2771..699a0a75 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -260,28 +260,26 @@ def __set_primary_secondary_heat_ports(): raise Exception( f"{asset.name} must have one inport for electricity and one outport for gas" ) - + elif ( - asset.asset_type == "GeothermalSource" - and len(asset.out_ports) == 1 - and len(asset.in_ports) == 2 - ): - for p in [*asset.in_ports, *asset.out_ports]: + asset.asset_type == "GeothermalSource" + and len(asset.out_ports) == 1 + and len(asset.in_ports) == 2 + ): + for p in [*asset.in_ports, *asset.out_ports]: + + if isinstance(p, InPort) and isinstance(p.carrier, esdl.ElectricityCommodity): + port_map[p.id] = getattr(component, elec_in_suf) + elif isinstance(p, InPort) and isinstance(p.carrier, esdl.HeatCommodity): + port_map[p.id] = getattr(component, in_suf) + elif isinstance(p, OutPort): # OutPort + port_map[p.id] = getattr(component, out_suf) + else: + raise Exception( + f"{asset.name} has does not have (1 electricity in_port) 1 heat " + f"in port and 1 Heat out_ports " + ) - if isinstance(p, InPort) and isinstance( - p.carrier, esdl.ElectricityCommodity - ): - port_map[p.id] = getattr(component, elec_in_suf) - elif isinstance(p, InPort) and isinstance(p.carrier, esdl.HeatCommodity): - port_map[p.id] = getattr(component, in_suf) - elif isinstance(p, OutPort): # OutPort - port_map[p.id] = getattr(component, out_suf) - else: - raise Exception( - f"{asset.name} has does not have (1 electricity in_port) 1 heat " - f"in port and 1 Heat out_ports " - ) - elif ( asset.in_ports is None and isinstance(asset.out_ports[0].carrier, esdl.ElectricityCommodity) diff --git a/src/mesido/financial_mixin.py b/src/mesido/financial_mixin.py index 09783d79..ac91ef31 100644 --- a/src/mesido/financial_mixin.py +++ b/src/mesido/financial_mixin.py @@ -1549,9 +1549,11 @@ def __annualized_capex_constraints(self, ensemble_member): # Input is assumed as as annual percentage discount_percentage = parameters[f"{asset_name}.discount_rate"] if np.isnan(asset_life_years) or np.isnan(discount_percentage): - logger.warning(f"Annualized cost cannot be computed for \ + logger.warning( + f"Annualized cost cannot be computed for \ {asset_name} since technical_life \ - or discount_rate are not set.") + or discount_rate are not set." + ) continue symbol_name = self._annualized_capex_var_map[asset_name] diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py index a398c94d..5ff5ffd6 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py @@ -10,10 +10,10 @@ class GeothermalSource(HeatSource): """ The geothermal source component is used to model geothermal doublets. It is equivilent to a - normal source with the only difference being the modelling of doublets and power consumption. - The main reason for this component instead of using just a regular source is that to have the - integer behaviour of increasing the amount of doublets. In the HeatMixin an integer is created - _aggregation_count to model the amount of doublets and the maximum power will scale with this + normal source with the only difference being the modelling of doublets and power consumption. + The main reason for this component instead of using just a regular source is that to have the + integer behaviour of increasing the amount of doublets. In the HeatMixin an integer is created + _aggregation_count to model the amount of doublets and the maximum power will scale with this integer instead of continuous. This will also ensure that the cost will scale with this integer. The power consumption is computed through a COP calculation, directly linked to the heat source production. COP is set to a default value of 0 in order to ensure power and its associated costs @@ -41,4 +41,4 @@ def __init__(self, name, **modifiers): self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) - self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) + self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index 918f4552..13b87e66 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -34,5 +34,3 @@ def __init__(self, name, **modifiers): self.add_variable(ElectricityPort, "ElectricityIn") self.add_equation(((self.ElectricityIn.Power - self.Power_elec) / self.elec_power_nominal)) - - diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index 02db88be..ca9e9633 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -323,6 +323,7 @@ def solver_options(self): ) np.testing.assert_allclose(var_opex_hp_calc, var_opex_hp) + class TestGeothermalSourceElec(TestCase): def test_geothermal_source(self): @@ -409,8 +410,10 @@ def test_geothermal_source_elec(self): # Test electricity port np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.ElectricityIn.Power"]) - np.testing.assert_allclose(results["GeothermalSource_a77b.ElectricityIn.Power"], - results["GeothermalSource_a77b.Power_elec"]) + np.testing.assert_allclose( + results["GeothermalSource_a77b.ElectricityIn.Power"], + results["GeothermalSource_a77b.Power_elec"], + ) # Variable operational cost check. np.testing.assert_allclose( @@ -423,9 +426,9 @@ def test_geothermal_source_elec(self): if __name__ == "__main__": - #test_cold_demand = TestMultiCommodityHeatPump() - #test_cold_demand.test_air_to_water_heat_pump_elec_min_elec() + # test_cold_demand = TestMultiCommodityHeatPump() + # test_cold_demand.test_air_to_water_heat_pump_elec_min_elec() # test_cold_demand.test_heat_pump_elec_min_elec() test_geothermal = TestGeothermalSourceElec() test_geothermal.test_geothermal_source_elec() - #test_geothermal.test_geothermal_source() \ No newline at end of file + # test_geothermal.test_geothermal_source() From 0d45e48cef302c5b1ad52126dce8ba20ac0d3fa3 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 16:32:03 +0100 Subject: [PATCH 13/23] Linting. --- src/mesido/esdl/esdl_heat_model.py | 11 +++++------ src/mesido/esdl/esdl_model_base.py | 8 ++++---- .../milp/heat/geothermal_source_elec.py | 8 ++++---- 3 files changed, 13 insertions(+), 14 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 6345f2dd..ac77d057 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1318,12 +1318,11 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 10.0 # Check to see if there is an electricity carrier at the in ports. in_port_carriers = [port.carrier for port in asset.in_ports] - if any( - [ - isinstance(carrier, esdl.esdl.ElectricityCommodity) - for carrier in in_port_carriers - ] - ): + elec_port = False + for carrier in in_port_carriers: + if isinstance(carrier, esdl.esdl.ElectricityCommodity) + elec_port = True + if elec_port: for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): min_voltage = port.carrier.voltage diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index 699a0a75..4ee1eace 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -202,10 +202,10 @@ def __set_primary_secondary_heat_ports(): else: if asset.asset_type == "HeatPump": raise Exception( - f"{asset.name} has incorrect number of in/out ports. HeatPumps are allows " - f"to have 1 in and 1 out port for air-water HP, 2 in ports and 2 out ports " - f"when modelling a water-water HP, or 3 in ports and 2 out ports when the " - f"electricity connection of the water-water HP is modelled." + f"{asset.name} has incorrect number of in/out ports. HeatPumps allow " + f"to have 1 in and 1 out port for air-water HP, 2 in ports and 2 out " + f"ports when modelling a water-water HP, or 3 in ports and 2 out ports " + f"when the electricity connection of the water-water HP is modelled." ) elif ( diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index 13b87e66..8b7fa7c1 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -8,10 +8,10 @@ @add_variables_documentation_automatically class GeothermalSourceElec(GeothermalSource): """ - The geothermal source electric asset is almost identical to the geothermal source one. The only - difference is that this one includes an electricity in port. The electricity power calculation that - it inherits from the geothermal source asset needs to be satisfied by an electricity carrier supplied - through this new in port. + The geothermal source electric asset is almost identical to the geothermal source one. The + only difference is that this one includes an electricity in port. The electricity power + calculation that it inherits from the geothermal source asset needs to be satisfied by an + electricity carrier supplied through this new in port. Variables created: {add_variable_names_for_documentation_here} From 0169bed7771fb97eecd2cd50c8549e9fca5efa90 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 16:32:57 +0100 Subject: [PATCH 14/23] Fixed mistake. --- src/mesido/esdl/esdl_heat_model.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index ac77d057..75b60f11 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1320,8 +1320,8 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS in_port_carriers = [port.carrier for port in asset.in_ports] elec_port = False for carrier in in_port_carriers: - if isinstance(carrier, esdl.esdl.ElectricityCommodity) - elec_port = True + if isinstance(carrier, esdl.esdl.ElectricityCommodity): + elec_port = True if elec_port: for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): From dd80124bc0e25a54aa7352fd69c64f839103c21c Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 16:38:06 +0100 Subject: [PATCH 15/23] Linting. --- .../component_library/milp/heat/geothermal_source_elec.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index 8b7fa7c1..e5718333 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -8,9 +8,9 @@ @add_variables_documentation_automatically class GeothermalSourceElec(GeothermalSource): """ - The geothermal source electric asset is almost identical to the geothermal source one. The - only difference is that this one includes an electricity in port. The electricity power - calculation that it inherits from the geothermal source asset needs to be satisfied by an + The geothermal source electric asset is almost identical to the geothermal source one. The + only difference is that this one includes an electricity in port. The electricity power + calculation that it inherits from the geothermal source asset needs to be satisfied by an electricity carrier supplied through this new in port. Variables created: From f3a482baf03827781eed613dbe57e9c47d142afe Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Mon, 2 Mar 2026 16:54:34 +0100 Subject: [PATCH 16/23] Fixed formatting. --- src/mesido/financial_mixin.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/mesido/financial_mixin.py b/src/mesido/financial_mixin.py index ac91ef31..09783d79 100644 --- a/src/mesido/financial_mixin.py +++ b/src/mesido/financial_mixin.py @@ -1549,11 +1549,9 @@ def __annualized_capex_constraints(self, ensemble_member): # Input is assumed as as annual percentage discount_percentage = parameters[f"{asset_name}.discount_rate"] if np.isnan(asset_life_years) or np.isnan(discount_percentage): - logger.warning( - f"Annualized cost cannot be computed for \ + logger.warning(f"Annualized cost cannot be computed for \ {asset_name} since technical_life \ - or discount_rate are not set." - ) + or discount_rate are not set.") continue symbol_name = self._annualized_capex_var_map[asset_name] From e2112bed298ad2f172de8e7800f6ff546bff859a Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 10:37:45 +0100 Subject: [PATCH 17/23] Implemented comments from PR. --- src/mesido/esdl/esdl_heat_model.py | 9 +-------- src/mesido/esdl/esdl_model_base.py | 21 +-------------------- tests/test_multicommodity.py | 3 ++- 3 files changed, 4 insertions(+), 29 deletions(-) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index 75b60f11..f12a05cc 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1313,16 +1313,9 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS f"{asset.asset_type} '{asset.name}' has no desired flow rate specified. " f"'{asset.name}' will not be actuated in a constant manner" ) - max_supply = asset.attributes["power"] modifiers["elec_power_nominal"] = max_supply modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 10.0 - # Check to see if there is an electricity carrier at the in ports. - in_port_carriers = [port.carrier for port in asset.in_ports] - elec_port = False - for carrier in in_port_carriers: - if isinstance(carrier, esdl.esdl.ElectricityCommodity): - elec_port = True - if elec_port: + if len(asset.in_ports) == 2: for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): min_voltage = port.carrier.voltage diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index 4ee1eace..84018874 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -227,7 +227,7 @@ def __set_primary_secondary_heat_ports(): f"Heat out_ports " ) elif ( - asset.asset_type == "ElectricBoiler" + asset.asset_type == "ElectricBoiler" or asset.asset_type == "GeothermalSource" and len(asset.out_ports) == 1 and len(asset.in_ports) == 2 ): @@ -261,25 +261,6 @@ def __set_primary_secondary_heat_ports(): f"{asset.name} must have one inport for electricity and one outport for gas" ) - elif ( - asset.asset_type == "GeothermalSource" - and len(asset.out_ports) == 1 - and len(asset.in_ports) == 2 - ): - for p in [*asset.in_ports, *asset.out_ports]: - - if isinstance(p, InPort) and isinstance(p.carrier, esdl.ElectricityCommodity): - port_map[p.id] = getattr(component, elec_in_suf) - elif isinstance(p, InPort) and isinstance(p.carrier, esdl.HeatCommodity): - port_map[p.id] = getattr(component, in_suf) - elif isinstance(p, OutPort): # OutPort - port_map[p.id] = getattr(component, out_suf) - else: - raise Exception( - f"{asset.name} has does not have (1 electricity in_port) 1 heat " - f"in port and 1 Heat out_ports " - ) - elif ( asset.in_ports is None and isinstance(asset.out_ports[0].carrier, esdl.ElectricityCommodity) diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index ca9e9633..f8f453eb 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -357,7 +357,6 @@ def test_geothermal_source(self): electric_power_conservation_test(heat_problem, results) # Equations check - np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) np.testing.assert_allclose( parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], results["GeothermalSource_a77b.Heat_source"], @@ -407,6 +406,8 @@ def test_geothermal_source_elec(self): parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], results["GeothermalSource_a77b.Heat_source"], ) + np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], + results["GeothermalSource_a77b.ElectricityIn.Power"]) # Test electricity port np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.ElectricityIn.Power"]) From e815dbca01265aa2fc6f099623fcc709fa35dc2e Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 10:39:16 +0100 Subject: [PATCH 18/23] Fixed test. --- tests/test_multicommodity.py | 5 ++--- 1 file changed, 2 insertions(+), 3 deletions(-) diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index f8f453eb..ef19c05a 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -401,13 +401,12 @@ def test_geothermal_source_elec(self): # Equations check np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) - np.testing.assert_array_less(0.0, results["ElectricityProducer_4dde.ElectricityOut.Power"]) + np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], + results["GeothermalSource_a77b.ElectricityIn.Power"]) np.testing.assert_allclose( parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], results["GeothermalSource_a77b.Heat_source"], ) - np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], - results["GeothermalSource_a77b.ElectricityIn.Power"]) # Test electricity port np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.ElectricityIn.Power"]) From 2ff0dfb428968fee440a0fbeeefecd405eb7f660 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 14:46:48 +0100 Subject: [PATCH 19/23] Modified behavior for cases where no cop is provided. Added a test for it. --- src/mesido/esdl/esdl_heat_model.py | 2 +- .../milp/heat/geothermal_source.py | 5 +- .../sourcesink_with_geo_elec_no_cop.esdl | 70 +++++++++++++++++++ tests/test_multicommodity.py | 52 ++++++++++++-- 4 files changed, 120 insertions(+), 9 deletions(-) create mode 100644 tests/models/source_pipe_sink/model/sourcesink_with_geo_elec_no_cop.esdl diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index f12a05cc..f1a997f7 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -1314,7 +1314,7 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS f"'{asset.name}' will not be actuated in a constant manner" ) modifiers["elec_power_nominal"] = max_supply - modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 10.0 + modifiers["cop"] = asset.attributes["COP"] if asset.attributes["COP"] else 0.0 if len(asset.in_ports) == 2: for port in asset.in_ports: if isinstance(port.carrier, esdl.ElectricityCommodity): diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py index 5ff5ffd6..b01049fb 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source.py @@ -41,4 +41,7 @@ def __init__(self, name, **modifiers): self.add_variable(Variable, "Power_elec", min=0.0, nominal=self.elec_power_nominal) - self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) + if self.cop == 0.0: + self.add_equation((self.Power_elec)) + else: + self.add_equation(((self.Power_elec * self.cop - self.Heat_source) / self.Heat_nominal)) diff --git a/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec_no_cop.esdl b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec_no_cop.esdl new file mode 100644 index 00000000..a164e3e6 --- /dev/null +++ b/tests/models/source_pipe_sink/model/sourcesink_with_geo_elec_no_cop.esdl @@ -0,0 +1,70 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index ef19c05a..79cb114c 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -422,13 +422,51 @@ def test_geothermal_source_elec(self): / parameters["GeothermalSource_a77b.cop"], results["GeothermalSource_a77b__variable_operational_cost"], ) + + def test_geothermal_source_elec_no_cop(self): + """ + This tests checks the electric geothermal producer asset when no cop is provided. + + It is the same test case as with the regular one, but in this case, the cop is not + provided and it assigned a value of 0.0, resulting in no power related costs. + """ + import models.source_pipe_sink.src.double_pipe_heat as example + from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink + + base_folder = Path(example.__file__).resolve().parent.parent + + heat_problem = run_esdl_mesido_optimization( + SourcePipeSink, + base_folder=base_folder, + esdl_file_name="sourcesink_with_geo_elec_no_cop.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_import.csv", + ) + results = heat_problem.extract_results() + + # Standard checks. + demand_matching_test(heat_problem, results) + energy_conservation_test(heat_problem, results) + heat_to_discharge_test(heat_problem, results) + electric_power_conservation_test(heat_problem, results) + + # Equations check + np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) + np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], + results["GeothermalSource_a77b.ElectricityIn.Power"]) + np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], + 0.0) + # Test electricity port + np.testing.assert_allclose( + results["GeothermalSource_a77b.ElectricityIn.Power"], + results["GeothermalSource_a77b.Power_elec"], + ) -if __name__ == "__main__": + # Variable operational cost check. + np.testing.assert_allclose( + 0.0, + results["GeothermalSource_a77b__variable_operational_cost"] + ) - # test_cold_demand = TestMultiCommodityHeatPump() - # test_cold_demand.test_air_to_water_heat_pump_elec_min_elec() - # test_cold_demand.test_heat_pump_elec_min_elec() - test_geothermal = TestGeothermalSourceElec() - test_geothermal.test_geothermal_source_elec() - # test_geothermal.test_geothermal_source() From a7d46453304267f0e37ca2a9fd3cfbf98b8894d6 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 15:06:20 +0100 Subject: [PATCH 20/23] Merged main. --- CHANGELOG.md | 18 +- heat_producer) | 0 manual_run/testing_bugs/src/run_grow.py | 1 + setup.py | 2 +- src/mesido/asset_sizing_mixin.py | 253 +++++--- src/mesido/esdl/asset_to_component_base.py | 6 +- src/mesido/esdl/esdl_additional_vars_mixin.py | 25 +- src/mesido/esdl/esdl_heat_model.py | 240 ++++++-- src/mesido/esdl/esdl_mixin.py | 124 +++- src/mesido/esdl/profile_parser.py | 49 +- src/mesido/potential_errors.py | 1 + src/mesido/workflows/gas_elect_workflow.py | 3 +- src/mesido/workflows/grow_workflow.py | 7 +- src/mesido/workflows/utils/error_types.py | 2 + .../model/testing_network_1_all_enabled.esdl | 329 ++++++++++ .../model/testing_network_1_all_optional.esdl | 333 ++++++++++ ...etwork_1_all_optional_excluding_pipes.esdl | 333 ++++++++++ .../model/testing_network_1_supply_only.esdl | 210 +++++++ tests/models/esdl_constraints/src/run_1.py | 31 + ...ork_all_optional_pipe_catalog_pipe_DN.esdl | 574 ++++++++++++++++++ .../input/timeseries_with_pv.csv | 4 + .../model/pv_with_csv_profile.esdl | 37 ++ .../source_sink_cable/src/example.py | 88 ++- tests/test_electric_source_sink.py | 59 +- tests/test_end_scenario_sizing.py | 57 ++ tests/test_esdl_constraints.py | 295 +++++++++ tests/test_potential_errors.py | 2 + tests/test_producer_profiles.py | 4 +- tests/test_profile_parsing.py | 8 +- 29 files changed, 2915 insertions(+), 180 deletions(-) delete mode 100644 heat_producer) create mode 100644 tests/models/esdl_constraints/model/testing_network_1_all_enabled.esdl create mode 100644 tests/models/esdl_constraints/model/testing_network_1_all_optional.esdl create mode 100644 tests/models/esdl_constraints/model/testing_network_1_all_optional_excluding_pipes.esdl create mode 100644 tests/models/esdl_constraints/model/testing_network_1_supply_only.esdl create mode 100644 tests/models/esdl_constraints/src/run_1.py create mode 100644 tests/models/test_case_small_network_ates_buffer_optional_assets/model/test_case_small_network_all_optional_pipe_catalog_pipe_DN.esdl create mode 100644 tests/models/unit_cases_electricity/source_sink_cable/input/timeseries_with_pv.csv create mode 100644 tests/models/unit_cases_electricity/source_sink_cable/model/pv_with_csv_profile.esdl create mode 100644 tests/test_esdl_constraints.py diff --git a/CHANGELOG.md b/CHANGELOG.md index bae68ff7..327abe81 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,10 +1,26 @@ -# [Unreleased-main] - 2026-02-25 +# [Unreleased-main] - 2026-03-12 + +## Added +- xxx + +## Changed +- xxx + +## Fixed +- xxx + + +# [0.1.18] - 2026-03-12 ## Added - Parsing of ensemble profiles when using input profiles from csv. +- Maximum profile constraint for PV asset is considered in PV sizing +- Cater for input via esdl constraints to specify the upper limit for OPTIONAL assets in DTK +- Initial implementation of adaptable pipe DN lower limit per pipe ## Changed - Speed-up timeseries check in from InfluxDB +- Delete solution after stage 1 in a staged workflow ## Fixed - xxx diff --git a/heat_producer) b/heat_producer) deleted file mode 100644 index e69de29b..00000000 diff --git a/manual_run/testing_bugs/src/run_grow.py b/manual_run/testing_bugs/src/run_grow.py index 1d2f7662..06668149 100644 --- a/manual_run/testing_bugs/src/run_grow.py +++ b/manual_run/testing_bugs/src/run_grow.py @@ -11,6 +11,7 @@ base_folder = Path(__file__).resolve().parent.parent kwargs = { + "use_esdl_ranged_constraint": True, # default value in the code is set to False "database_connections": [ { "access_type": DBAccessType.READ, # or DBAccessType.WRITE or DBAccessType.READ_WRITE diff --git a/setup.py b/setup.py index c913a76e..e8af0ab9 100644 --- a/setup.py +++ b/setup.py @@ -63,7 +63,7 @@ # < 81.0.0 needed for pandapipes (still to be removed) # < 82.0.0 needed for pkg_resources (used in rtctools) "setuptools <= 80.9.0", - "pyesdl == 25.7", + "pyesdl == 26.2", "pandas >= 1.3.1, < 2.0", "casadi-gil-comp == 3.6.7", "StrEnum == 0.4.15", diff --git a/src/mesido/asset_sizing_mixin.py b/src/mesido/asset_sizing_mixin.py index 709c1504..e3c6948f 100644 --- a/src/mesido/asset_sizing_mixin.py +++ b/src/mesido/asset_sizing_mixin.py @@ -10,6 +10,8 @@ from mesido.base_component_type_mixin import BaseComponentTypeMixin from mesido.demand_insulation_class import DemandInsulationClass from mesido.esdl.asset_to_component_base import AssetStateEnum +from mesido.esdl.common import Asset +from mesido.esdl.esdl_additional_vars_mixin import get_asset_contraints from mesido.head_loss_class import HeadLossOption from mesido.network_common import NetworkSettings from mesido.pipe_class import CableClass, GasPipeClass, PipeClass @@ -801,13 +803,14 @@ def _make_max_size_var(name, lb, ub, nominal): ub = bounds[f"{asset_name}.Heat_source"][1] # Update bound to account for profile constraint being used instead of 1 value - esdl_asset_attributes = self.esdl_assets[ - self.esdl_asset_name_to_id_map[asset_name] - ].attributes["constraint"] + asset = self.esdl_assets[self.esdl_asset_name_to_id_map[asset_name]] + asset_profile_constraints, qty_asset_profile_constraints = get_asset_contraints( + self, asset, esdl.ProfileConstraint + ) if ( - len(esdl_asset_attributes) > 0 - and hasattr(esdl_asset_attributes.items[0], "maximum") - and esdl_asset_attributes.items[0].maximum.profileQuantityAndUnit.reference.unit + qty_asset_profile_constraints > 0 + and hasattr(asset_profile_constraints[0], "maximum") + and asset_profile_constraints[0].maximum.profileQuantityAndUnit.reference.unit == esdl.UnitEnum.WATT and parameters[f"{asset_name}.state"] == AssetStateEnum.OPTIONAL # Optional asset ): @@ -1857,6 +1860,133 @@ def __electricity_cable_topology_path_constraints(self, ensemble_member): return constraints + def __producer_constraints( + self, + asset: Asset, + variable_suffix: str, + source: ca.MX, + max_power: ca.MX, + constraint_nominal: float, + ensemble_member: int, + ) -> list: + """ + Apply production‑capacity constraints based on an asset’s time‑series profile. + + This method caps heat or electricity output using either: + - an absolute profile in Watts (fixed maximum or optional scaling), or + - a normalized (0–1) profile that scales with installed capacity. + + Parameters + ---------- + asset : Asset object. + variable_suffix : Timeseries suffix indicating the maximum production profile + (e.g., "maximum_heat_source", "maximum_electricity_source"). + source : Vector of path_variables of produced energy to be constrained. + max_power : Installed or decision‑variable capacity of the asset. + constraint_nominal : Normalization constant for constraints. + ensemble_member : Ensemble index used to retrieve asset parameters. + + Returns + ------- + constraints : The function returns a list of profile‑based capacity constraints. + """ + constraints = [] + if f"{asset.name}.{variable_suffix}" in self.io.get_timeseries_names(): + timeseries = self.get_timeseries(f"{asset.name}.{variable_suffix}") + if len(self.times()) < len(timeseries.times): + idx_start = np.where(timeseries.times == self.times()[0])[0][0] + idx_end = np.where(timeseries.times == self.times()[-1])[0][0] + profile_non_scaled = timeseries.values[idx_start : idx_end + 1] + else: + profile_non_scaled = timeseries.values + max_profile_non_scaled = max(profile_non_scaled) + profile_scaled = profile_non_scaled / max_profile_non_scaled + + # Cap the energy produced via a profile. Two profile options below. + # Option 1: Profile specified in absolute values [W] via a ProfileConstraint + asset_profile_constraints, qty_asset_profile_constraints = get_asset_contraints( + self, asset, esdl.ProfileConstraint + ) + if ( + qty_asset_profile_constraints > 0 + and hasattr(asset_profile_constraints[0], "maximum") + and asset_profile_constraints[0].maximum.profileQuantityAndUnit.reference.unit + == esdl.UnitEnum.WATT + ): + parameters = self.parameters(ensemble_member) + asset_state = parameters[f"{asset.name}.state"] + + if asset_state == AssetStateEnum.ENABLED: # Enabled asset + constraints.append( + ( + (max_power - max_profile_non_scaled) / constraint_nominal, + 0.0, + np.inf, + ) + ) + max_power_var = ( + max_profile_non_scaled # maximum heat power or maximum electric power + ) + + elif asset_state == AssetStateEnum.OPTIONAL: # Optional asset + max_power_var = max_power + + else: + state_val = asset_state + logger.error(f"Unexpected state: {state_val}") + sys.exit(1) + + constraints.append( + ( + (profile_scaled * max_power_var - source) / constraint_nominal, + 0.0, + np.inf, + ) + ) + # Option 2: Normalised profile (0.0-1.0) shape that scales with maximum size of the + # producer + # Note: If the asset is not optional then the profile will be scaled to the + # installed capacity + elif ( + # profile is specified without units (xlm/csv) + qty_asset_profile_constraints == 0 + or ( + asset_profile_constraints[ + 0 + ].maximum.profileQuantityAndUnit.reference.physicalQuantity + == esdl.PhysicalQuantityEnum.COEFFICIENT + and ( + asset_profile_constraints[0].maximum.profileQuantityAndUnit.reference.unit + == esdl.UnitEnum.PERCENT + or asset_profile_constraints[ + 0 + ].maximum.profileQuantityAndUnit.reference.unit + == esdl.UnitEnum.NONE + ) + ) # profile from esdl + ): + # TODO: currently this can only be used with a csv file since units must be set + # for ProfileContraint. Future addition can be to use a different unit/quantity + # etc. so that the profile is used in a normalised way and scale to max_size + constraints.append( + ( + (profile_scaled * max_power - source) / constraint_nominal, + 0.0, + np.inf, + ) + ) + else: + RuntimeError(f"{asset.name}: Unforeseen error in adding a profile constraint") + else: + constraints.append( + ( + (np.ones(len(self.times())) * max_power - source) / constraint_nominal, + 0.0, + np.inf, + ) + ) + return constraints + def __max_size_constraints(self, ensemble_member): """ This function makes sure that the __max_size variable is at least as large as needed. For @@ -1919,95 +2049,16 @@ def __max_size_constraints(self, ensemble_member): max_heat = self.extra_variable(max_var, ensemble_member) heat_source = self.__state_vector_scaled(f"{s}.Heat_source", ensemble_member) constraint_nominal = self.variable_nominal(f"{s}.Heat_source") - - if f"{s}.maximum_heat_source" in self.io.get_timeseries_names(): - profile_non_scaled = self.get_timeseries(f"{s}.maximum_heat_source").values - max_profile_non_scaled = max(profile_non_scaled) - profile_scaled = profile_non_scaled / max_profile_non_scaled - - # Cap the heat produced via a profile. Two profile options below. - # Option 1: Profile specified in absolute values [W] via a ProfileConstraint - esdl_asset_attributes = self.esdl_assets[ - self.esdl_asset_name_to_id_map[s] - ].attributes["constraint"] - if ( - len(esdl_asset_attributes) > 0 - and hasattr(esdl_asset_attributes.items[0], "maximum") - and esdl_asset_attributes.items[0].maximum.profileQuantityAndUnit.reference.unit - == esdl.UnitEnum.WATT - ): - parameters = self.parameters(ensemble_member) - - if parameters[f"{s}.state"] == AssetStateEnum.ENABLED: # Enabled asset - constraints.append( - ( - (max_heat - max_profile_non_scaled) / constraint_nominal, - 0.0, - np.inf, - ) - ) - max_heat_var = max_profile_non_scaled - - elif parameters[f"{s}.state"] == AssetStateEnum.OPTIONAL: # Optional asset - max_heat_var = max_heat - - else: - state_val = parameters[f"{s}.state"] - logger.error(f"Unexpected state: {state_val}") - sys.exit(1) - - constraints.append( - ( - (profile_scaled * max_heat_var - heat_source) / constraint_nominal, - 0.0, - np.inf, - ) - ) - # Option 2: Normalised profile (0.0-1.0) shape that scales with maximum size of the - # producer - # Note: If the asset is not optional then the profile will be scaled to the - # installed capacity - elif ( - # profile is specified without units (xlm/csv) - len(esdl_asset_attributes) == 0 - or ( - esdl_asset_attributes.items[ - 0 - ].maximum.profileQuantityAndUnit.reference.physicalQuantity - == esdl.PhysicalQuantityEnum.COEFFICIENT - and ( - esdl_asset_attributes.items[ - 0 - ].maximum.profileQuantityAndUnit.reference.unit - == esdl.UnitEnum.PERCENT - or esdl_asset_attributes.items[ - 0 - ].maximum.profileQuantityAndUnit.reference.unit - == esdl.UnitEnum.NONE - ) - ) # profile from esdl - ): - # TODO: currently this can only be used with a csv file since units must be set - # for ProfileContraint. Future addition can be to use a different unit/quantity - # etc. so that the profile is used in a normalised way and scale to max_size - - constraints.append( - ( - (profile_scaled * max_heat - heat_source) / constraint_nominal, - 0.0, - np.inf, - ) - ) - else: - RuntimeError(f"{s}: Unforeseen error in adding a profile contraint") - else: - constraints.append( - ( - (np_ones * max_heat - heat_source) / constraint_nominal, - 0.0, - np.inf, - ) + constraints.extend( + self.__producer_constraints( + self.esdl_assets[self.esdl_asset_name_to_id_map[s]], + "maximum_heat_source", + heat_source, + max_heat, + constraint_nominal, + ensemble_member, ) + ) # Constraint the aggregation_count to 0.0 when asset is not placed. And also ensure the # max_size_var is a factor of single_doublet_power @@ -2142,12 +2193,14 @@ def __max_size_constraints(self, ensemble_member): f"{d}.Electricity_source", ensemble_member ) constraint_nominal = self.variable_nominal(f"{d}.Electricity_source") - - constraints.append( - ( - (np_ones * max_power - electricity_source) / constraint_nominal, - 0.0, - np.inf, + constraints.extend( + self.__producer_constraints( + self.esdl_assets[self.esdl_asset_name_to_id_map[d]], + "maximum_electricity_source", + electricity_source, + max_power, + constraint_nominal, + ensemble_member, ) ) diff --git a/src/mesido/esdl/asset_to_component_base.py b/src/mesido/esdl/asset_to_component_base.py index c402bc51..532ab7f6 100644 --- a/src/mesido/esdl/asset_to_component_base.py +++ b/src/mesido/esdl/asset_to_component_base.py @@ -141,7 +141,8 @@ def get_density( # pressure and temperature. if carrier is None: logger.warning( - f"Neither gas/hydrogen/heat was used in the carrier at asset named {asset_name}." + f"Neither gas/hydrogen/heat was specified at asset named {asset_name}. " + "Gas properties will be used." ) density = cP.CoolProp.PropsSI( "D", @@ -201,7 +202,8 @@ def get_density( ) else: logger.warning( - f"Neither gas/hydrogen/heat was used in the carrier at asset named {asset_name}." + f"Neither gas/hydrogen/heat was specified in the carrier at asset named {asset_name}. " + "Gas properties will be used." ) density = 6.2 # natural gas at about 8 bar diff --git a/src/mesido/esdl/esdl_additional_vars_mixin.py b/src/mesido/esdl/esdl_additional_vars_mixin.py index 926d43d8..290bd174 100644 --- a/src/mesido/esdl/esdl_additional_vars_mixin.py +++ b/src/mesido/esdl/esdl_additional_vars_mixin.py @@ -1,8 +1,9 @@ import logging -from typing import List +from typing import List, Union import esdl +from mesido.esdl.common import Asset from mesido.network_common import NetworkSettings from mesido.pipe_class import PipeClass @@ -285,3 +286,25 @@ def temperature_regimes(self, carrier): self.__temperature_options[carrier] = temperature_options return temperature_options + + +def get_asset_contraints( + self, asset: Asset, constraint_type: esdl.Constraint +) -> Union[List[esdl.Constraint], int]: + """ + Get the contraints at an asset and the qty thereof. + + Arg: + asset: mesido common asset with all attributes + constraint type: the type of contraint specified (e.g. esdl.RangedConstraint, + esdl.ProfileConstraint) + + Returns: + - Contraints of a specific type + - Number of constraints of specific type that exists + """ + + asset_constraints = [ + cnst for cnst in asset.attributes["constraint"] if isinstance(cnst, constraint_type) + ] + return asset_constraints, len(asset_constraints) diff --git a/src/mesido/esdl/esdl_heat_model.py b/src/mesido/esdl/esdl_heat_model.py index f1a997f7..1af489b5 100644 --- a/src/mesido/esdl/esdl_heat_model.py +++ b/src/mesido/esdl/esdl_heat_model.py @@ -2,6 +2,7 @@ import inspect import logging import math +import sys from typing import Any, Dict, Tuple, Type, Union import esdl @@ -14,6 +15,7 @@ get_energy_content, ) from mesido.esdl.common import Asset +from mesido.esdl.esdl_additional_vars_mixin import get_asset_contraints from mesido.esdl.esdl_model_base import _ESDLModelBase from mesido.potential_errors import MesidoAssetIssueType, get_potential_errors from mesido.pycml.component_library.milp import ( @@ -58,6 +60,9 @@ WindPark, ) +# Importing workflow utilities at module import time can create circular +# imports when workflows import ESDL mixins. Import locally where needed. + from scipy.optimize import fsolve logger = logging.getLogger("mesido") @@ -158,6 +163,8 @@ def __init__( self.primary_port_name_convention = kwargs["primary_port_name_convention"] if "secondary_port_name_convention" in kwargs.keys(): self.secondary_port_name_convention = kwargs["secondary_port_name_convention"] + self.energy_system_esdl_version = kwargs.get("energy_system_esdl_version", None) + self.use_esdl_ranged_constraint = kwargs.get("use_esdl_ranged_constraint", False) @property def _rho_cp_modifiers(self) -> Dict: @@ -324,6 +331,123 @@ def _generic_heat_modifiers(min_heat=None, max_heat=None, q_nominal=None) -> Dic return modifiers + def _validate_attribute_value_not_zero( + self, + asset: Asset, + max_size_attribute: str, + max_value_attribute: float, + constraint_attribute: bool = False, + ) -> None: + """ + This function checks if the asset attribute type value > 0, else raise error if aplicable + + Args: + asset: mesido common asset with all attributes + max_size_attribute: type of attribute e.g. powrr, volume etc. + max_value_attribute: value that of the attribute, + constraint_attribute: Does the atrribute value originate from a constraint + """ + + msg = None + if not constraint_attribute: # value comes from asset attribute directly + msg = f"Asset named {asset.name}: The attribute {max_size_attribute} must be > 0." + else: # value comes from asset constraint attribute + msg = ( + f"Asset named {asset.name}: The maximum value in the range " + f"constraint for attribute {max_size_attribute} must be > 0." + ) + if not max_value_attribute > 0.0: + get_potential_errors().add_potential_issue( + MesidoAssetIssueType.ASSET_UPPER_LIMIT, + asset.id, + msg, + ) + # Raise the potential error here if applicable, with feedback to user + # Else a normal error exit might occer which will not give feedback to the + # user + # Import occurs here to prevent circular reference + from mesido.workflows.utils.error_types import potential_error_to_error + + potential_error_to_error(self._error_type_check) + + def _get_asset_max_size_input(self, asset: Asset, max_size_attribute: str) -> float: + """ + This function gets the max size value that is used as an upper limit for an asset's size + + Args: + asset: mesido common asset with all attributes + max_size_attribute: type of attribute (power or volume etc.) used as max size + + Returns: value that should be used as the max size value for an asset + """ + + asset_range_constraints, qty_asset_range_constraints = get_asset_contraints( + self, asset, esdl.RangedConstraint + ) + + if ( + self.energy_system_esdl_version is not None + and self.energy_system_esdl_version >= "v2602" + and self.use_esdl_ranged_constraint + # "v2602" contains the items needed for the ranged constraint implementation in MESIDO + ): + if asset.attributes["state"] == esdl.AssetStateEnum.OPTIONAL: + if qty_asset_range_constraints > 1: + logger.error( + f"Asset named {asset.name}: The code currently does not cater for more than" + " 1 RangedConstraint" + ) + sys.exit(1) + elif qty_asset_range_constraints == 0: + logger.warning( + "Expected a range contraint (upper size limit) for asset named " + f"{asset.name}, but none has been specified." + ) + get_potential_errors().add_potential_issue( + MesidoAssetIssueType.ASSET_UPPER_LIMIT, + asset.id, + f"Asset named {asset.name}: The upper limit of the asset size has to be " + "specified via a maximum value in a range constraint.", + ) + # Raise the potential error here if applicable, with feedback to user + # Else a normal error exit might occer which will not give feedback to the user + # Import occurs here to prevent circular reference + from mesido.workflows.utils.error_types import potential_error_to_error + + potential_error_to_error(self._error_type_check) + else: + logger.warning( + f"For asset named {asset.name}, the range constraint value is used for the " + f"asset's upper limit for the attribute {max_size_attribute}." + ) + max_value_range = asset_range_constraints[0].range.maxValue + + self._validate_attribute_value_not_zero( + asset, max_size_attribute, max_value_range, True + ) + + return max_value_range + + elif asset.attributes["state"] == esdl.AssetStateEnum.ENABLED: + if qty_asset_range_constraints > 0: + logger.warning( + f"The constraint that has been assigned to asset name {asset.name} is not " + "being used because the asset state has been specified as ENABLED." + ) + + max_value_attribute = asset.attributes[max_size_attribute] + + self._validate_attribute_value_not_zero( + asset, max_size_attribute, max_value_attribute + ) + + return max_value_attribute + + else: + exit(f"{asset.name}: asset state DISABLED is not supported yet") + else: # Catering for backwards compatibility + return asset.attributes[max_size_attribute] + def convert_heat_buffer(self, asset: Asset) -> Tuple[Type[HeatBuffer], MODIFIERS]: """ This function converts the buffer object in esdl to a set of modifiers that can be used in @@ -386,35 +510,34 @@ def convert_heat_buffer(self, asset: Asset) -> Tuple[Type[HeatBuffer], MODIFIERS f"Volume with value of {asset.attributes['volume']} m3 will be used." ) - capacity = 0.0 - if asset.attributes["volume"]: - capacity = ( - asset.attributes["volume"] - * self.rho - * self.cp - * (supply_temperature - return_temperature) + asset_capacity_joule = 0.0 + asset_volume_m3 = self._get_asset_max_size_input(asset, "volume") + + if asset_volume_m3: + asset_capacity_joule = ( + asset_volume_m3 * self.rho * self.cp * (supply_temperature - return_temperature) ) elif asset.attributes["capacity"]: - capacity = asset.attributes["capacity"] + asset_capacity_joule = self._get_asset_max_size_input(asset, "capacity") else: logger.error( f"{asset.asset_type} '{asset.name}' has both not capacity and volume specified. " f"Please specify one of the two" ) - assert capacity > 0.0 + assert asset_capacity_joule > 0.0 min_fraction_tank_volume = self.min_fraction_tank_volume if self.get_state(asset) == 0 or self.get_state(asset) == 2: min_fraction_tank_volume = 0.0 # We assume that the height equals the radius of the buffer. r = ( - capacity + asset_capacity_joule * (1 + min_fraction_tank_volume) / (self.rho * self.cp * (supply_temperature - return_temperature) * math.pi) ) ** (1.0 / 3.0) - min_heat = capacity * min_fraction_tank_volume - max_heat = capacity * (1 + min_fraction_tank_volume) + min_heat = asset_capacity_joule * min_fraction_tank_volume + max_heat = asset_capacity_joule * (1 + min_fraction_tank_volume) assert max_heat > 0.0 # default is set to 10MW @@ -765,9 +888,7 @@ def convert_gas_pipe(self, asset: Asset) -> Tuple[Type[GasPipe], MODIFIERS]: return GasPipe, modifiers - def convert_heat_pipe( - self, asset: Asset - ) -> Tuple[Union[Type[HeatPipe], Type[GasPipe]], MODIFIERS]: + def convert_heat_pipe(self, asset: Asset) -> Tuple[Type[HeatPipe], MODIFIERS]: """ This function converts the pipe object in esdl to a set of modifiers that can be used in a pycml object. Most important: @@ -1031,13 +1152,18 @@ def convert_heat_exchanger(self, asset: Asset) -> Tuple[Type[HeatExchanger], MOD f"transfer heat from primary to secondary.", ) + asset_power = None + asset_capacity = None if asset.asset_type == "GenericConversion": - max_power = asset.attributes["power"] if asset.attributes["power"] else math.inf + asset_power = self._get_asset_max_size_input(asset, "power") + max_power = asset_power if asset_power else math.inf else: # DTK requires capacity as the maximum power reference and not based on # heatTransferCoefficient. Power could also be based on heatTransferCoefficient if we # use an option to select it. - max_power = asset.attributes["capacity"] if asset.attributes["capacity"] else math.inf + asset_capacity = self._get_asset_max_size_input(asset, "capacity") + max_power = asset_capacity if asset_capacity else math.inf + if max_power == math.inf: get_potential_errors().add_potential_issue( MesidoAssetIssueType.HEAT_EXCHANGER_POWER, @@ -1180,10 +1306,10 @@ def convert_heat_pump(self, asset: Asset) -> Tuple[ else: cop = asset.attributes["COP"] - if not asset.attributes["power"]: + power_secondary = self._get_asset_max_size_input(asset, "power") + if not power_secondary: raise _ESDLInputException(f"{asset.name} has no power specified") else: - power_secondary = asset.attributes["power"] power_electrical = power_secondary / cop params_t = self._supply_return_temperature_modifiers(asset) @@ -1267,10 +1393,22 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS "ElectricBoiler", } - max_supply = asset.attributes["power"] + aggregation_count = 0 + if asset.asset_type == "GeothermalSource": + max_supply = asset.attributes["power"] + aggregation_count = self._get_asset_max_size_input(asset, "aggregationCount") + + if not aggregation_count: + logger.error( + f"{asset.asset_type} '{asset.name}' has no aggregation count specified." + ) + assert int(aggregation_count) == aggregation_count and aggregation_count > 0 + + else: + max_supply = self._get_asset_max_size_input(asset, "power") if not max_supply: - logger.error(f"{asset.asset_type} '{asset.name}' has no max power specified. ") + logger.error(f"{asset.asset_type} '{asset.name}' has no max power specified.") assert max_supply > 0.0 # get price per unit of energy, @@ -1290,15 +1428,15 @@ def convert_heat_source(self, asset: Asset) -> Tuple[Type[HeatSource], MODIFIERS ) if asset.asset_type == "GeothermalSource": - modifiers["nr_of_doublets"] = asset.attributes["aggregationCount"] + modifiers["nr_of_doublets"] = aggregation_count modifiers["Heat_source"] = dict( min=0.0, - max=max_supply * asset.attributes["aggregationCount"], + max=max_supply * aggregation_count, nominal=max_supply / 2.0, ) modifiers["Heat_flow"] = dict( min=0.0, - max=max_supply * asset.attributes["aggregationCount"], + max=max_supply * aggregation_count, nominal=max_supply / 2.0, ) try: @@ -1398,6 +1536,7 @@ def convert_ates(self, asset: Asset) -> Tuple[Type[ATES], MODIFIERS]: hfr_charge_max = asset.attributes.get("maxChargeRate", math.inf) hfr_discharge_max = asset.attributes.get("maxDischargeRate", math.inf) + single_doublet_power = hfr_discharge_max # We assume the efficiency is realized over a period of 100 days @@ -1412,28 +1551,31 @@ def convert_ates(self, asset: Asset) -> Tuple[Type[ATES], MODIFIERS]: cp = self.cp q_max_ates = hfr_discharge_max / (cp * rho * dt) - q_nominal = min( - self._get_connected_q_nominal(asset), q_max_ates * asset.attributes["aggregationCount"] - ) + aggregation_count = self._get_asset_max_size_input(asset, "aggregationCount") + if not aggregation_count: + logger.error(f"{asset.asset_type} '{asset.name}' has no aggregation count specified.") + assert int(aggregation_count) == aggregation_count and aggregation_count > 0 + + q_nominal = min(self._get_connected_q_nominal(asset), q_max_ates * aggregation_count) modifiers = dict( Q=dict( - min=-q_max_ates * asset.attributes["aggregationCount"], - max=q_max_ates * asset.attributes["aggregationCount"], + min=-q_max_ates * aggregation_count, + max=q_max_ates * aggregation_count, nominal=q_nominal, ), single_doublet_power=single_doublet_power, heat_loss_coeff=(1.0 - efficiency ** (1.0 / 100.0)) / (3600.0 * 24.0), - nr_of_doublets=asset.attributes["aggregationCount"], + nr_of_doublets=aggregation_count, Stored_heat=dict( min=0.0, - max=hfr_charge_max * asset.attributes["aggregationCount"] * 180.0 * 24 * 3600.0, - nominal=hfr_charge_max * asset.attributes["aggregationCount"] * 30.0 * 24 * 3600.0, + max=hfr_charge_max * aggregation_count * 180.0 * 24 * 3600.0, + nominal=hfr_charge_max * aggregation_count * 30.0 * 24 * 3600.0, ), **self._generic_modifiers(asset), **self._generic_heat_modifiers( - -hfr_discharge_max * asset.attributes["aggregationCount"], - hfr_charge_max * asset.attributes["aggregationCount"], + -hfr_discharge_max * aggregation_count, + hfr_charge_max * aggregation_count, q_nominal, ), **self._supply_return_temperature_modifiers(asset), @@ -1475,8 +1617,8 @@ def convert_ates(self, asset: Asset) -> Tuple[Type[ATES], MODIFIERS]: modifiers.update( dict( Heat_ates=dict( - min=-hfr_charge_max * asset.attributes["aggregationCount"], - max=hfr_discharge_max * asset.attributes["aggregationCount"], + min=-hfr_charge_max * aggregation_count, + max=hfr_discharge_max * aggregation_count, nominal=hfr_discharge_max / 2.0, ), T_amb=asset.attributes["aquiferMidTemperature"], @@ -2473,7 +2615,16 @@ def convert_heat_source_gas( """ assert asset.asset_type in {"GasHeater"} - max_supply = asset.attributes["power"] + max_supply = None + is_one_in_port = True if len(asset.in_ports) == 1 else False + if is_one_in_port: + max_supply = self._get_asset_max_size_input(asset, "power") + else: + # TODO: range constraint to be added instead of using this value for OPTIONAL asset + # The implementation of range constraints has not been tested for assets with more + # than one port, this should still be done. We need to be sure to check the proper + # quantity is used. + max_supply = asset.attributes["power"] if not max_supply: logger.error(f"{asset.asset_type} '{asset.name}' has no max power specified. ") @@ -2508,7 +2659,7 @@ def convert_heat_source_gas( Q_nominal_gas=q_nominal_gas_no_carrier, ) ) - if len(asset.in_ports) == 1: + if is_one_in_port: return HeatSourceGas, modifiers id_mapping = asset.global_properties["carriers"][asset.in_ports[0].carrier.id][ @@ -2581,7 +2732,12 @@ def convert_heat_source_elec( """ assert asset.asset_type in {"ElectricBoiler"} - max_supply = asset.attributes["power"] + max_supply = None + is_one_in_port = True if len(asset.in_ports) == 1 else False + if is_one_in_port: + max_supply = self._get_asset_max_size_input(asset, "power") + else: # TODO: range constraint to be added instead of using this value for OPTIONAL asset + max_supply = asset.attributes["power"] if not max_supply: logger.error(f"{asset.asset_type} '{asset.name}' has no max power specified. ") @@ -2589,7 +2745,7 @@ def convert_heat_source_elec( _, modifiers = self.convert_heat_source(asset) modifiers["elec_power_nominal"] = max_supply - if len(asset.in_ports) == 1: + if is_one_in_port: return HeatSourceElec, modifiers id_mapping = asset.global_properties["carriers"][asset.in_ports[0].carrier.id][ @@ -2717,6 +2873,8 @@ def __init__( assets: Dict[str, Asset], name_to_id_map: Dict[str, str], converter_class=AssetToHeatComponent, + esdl_version=str, + esdl_ranged_constraint_usage=bool, **kwargs, ): super().__init__(None) @@ -2727,6 +2885,8 @@ def __init__( **{ "primary_port_name_convention": self.primary_port_name_convention, "secondary_port_name_convention": self.secondary_port_name_convention, + "energy_system_esdl_version": esdl_version, + "esdl_ranged_constraint_usage": esdl_ranged_constraint_usage, }, } ) diff --git a/src/mesido/esdl/esdl_mixin.py b/src/mesido/esdl/esdl_mixin.py index 52cc77ad..be2623db 100644 --- a/src/mesido/esdl/esdl_mixin.py +++ b/src/mesido/esdl/esdl_mixin.py @@ -16,6 +16,7 @@ from mesido.esdl.asset_to_component_base import _AssetToComponentBase from mesido.esdl.common import Asset from mesido.esdl.edr_pipe_class import EDRGasPipeClass, EDRPipeClass +from mesido.esdl.esdl_additional_vars_mixin import get_asset_contraints from mesido.esdl.esdl_heat_model import ESDLHeatModel from mesido.esdl.esdl_model_base import _ESDLModelBase from mesido.esdl.esdl_parser import ESDLStringParser @@ -80,8 +81,11 @@ class ESDLMixin( __max_supply_temperature: Optional[float] = None + __esdl_ranged_constraint_usage: bool = False + # TODO: remove this once ESDL allows specifying a minimum pipe size for an optional pipe. __minimum_pipe_size_name: str = "DN150" + __use_user_defined_minimum_pipe_size: bool = False def __init__(self, *args, **kwargs) -> None: """ @@ -104,6 +108,8 @@ def __init__(self, *args, **kwargs) -> None: kwargs : esdl_string or esdl_file_name must be provided """ + self.__use_esdl_ranged_constraint: bool = kwargs.get("use_esdl_ranged_constraint", False) + self.esdl_parser_class: type = kwargs.get("esdl_parser", ESDLStringParser) esdl_string = kwargs.get("esdl_string", None) model_folder = kwargs.get("model_folder") @@ -182,11 +188,13 @@ def __init__(self, *args, **kwargs) -> None: energy_system=self.__energy_system_handler.energy_system, file_path=input_file_path, database_credentials=read_only_dbase_credentials, + use_esdl_ranged_contraint=self._ESDLMixin__use_esdl_ranged_constraint, ) else: # read from a file, no database credentials needed self.__profile_reader: BaseProfileReader = profile_reader_class( energy_system=self.__energy_system_handler.energy_system, file_path=input_file_path, + use_esdl_ranged_contraint=self._ESDLMixin__use_esdl_ranged_constraint, ) # This way we allow users to adjust the parsed ESDL assets @@ -197,7 +205,13 @@ def __init__(self, *args, **kwargs) -> None: name_to_id_map = {a.name: a.id for a in assets.values()} if isinstance(self, PhysicsMixin): - self.__model = ESDLHeatModel(assets, name_to_id_map, **self.esdl_heat_model_options()) + self.__model = ESDLHeatModel( + assets=assets, + name_to_id_map=name_to_id_map, + esdl_version=self._ESDLMixin__energy_system_handler.energy_system.esdlVersion, + use_esdl_ranged_constraint=self._ESDLMixin__use_esdl_ranged_constraint, + **self.esdl_heat_model_options(), + ) else: assert isinstance(self, QTHMixin) @@ -258,6 +272,83 @@ def pre(self) -> None: ) self.name_to_esdl_id_map[esdl_asset.name] = esdl_id + def find_index_of_pipe_or_next_up(self, pipes, target_dn): + """ + Find the index of the first pipe in the ordered pipes that meets the specified size + """ + for ii, p in enumerate(pipes): + if float(p.name.replace("DN", "")) >= target_dn: + return ii + return None + + def _get_pipe_max_size_input( + self, + asset: Asset, + ) -> str: + """ + This function gets the upper limit for the pipe DN. Either from the asset attribute input + DN size or from the PipeDiameterConstraint, depending on what is applicable. + Args: + asset: mesido common asset with all attributes + max_size_attribute: type of attribute e.g. powrr, volume etc. + max_value_attribute: value that of the attribute, + constraint_attribute: Does the atrribute value originate from a constraint + """ + # Backward compatibility: + # PipeDiameterConstraint vs pipe diameter attribute, as the value to be used for the + # maximum size of the pipe. We check the esdl version to determine which attribute to use + # for the maximum size of the pipe. + # While the front end is being updated to accomodate easy use of PipeDiameterConstraint: + # - If the esdl version caters for PipeDiameterConstraint: + # - check if a constraint exists, then use it. + # - if the contraint does not exist then use the diameter specific in the asset attribute + # - If PipeDiameterConstraint does not + esdl_version = self._ESDLMixin__energy_system_handler.energy_system.esdlVersion + # "v2602" contains the items needed for the ranged constraint implementation in MESIDO + if esdl_version is not None and esdl_version >= "v2602": + pipe_constraint, qty_pipe_constraint = get_asset_contraints( + self, asset, esdl.PipeDiameterConstraint + ) + + if qty_pipe_constraint > 1: + logger.exit( + f"More than 1 pipe diameter constraint has been specified to " + f"pipe named {asset.name}, currenlty only the 1st constraint is being used" + ) + exit(1) + elif ( + qty_pipe_constraint == 0 + or pipe_constraint[0].maximum == esdl.PipeDiameterEnum.VALUE_SPECIFIED + ): + logger.warning( + "Expected a pipe diameter contraint (upper size limit) for pipe named " + f"{asset.name}, but none has been specified. Therefore, the pipe diameter " + "specified in the asset attribute is being used." + ) + # Do not delete the code below. This will be enforced late, once the front end has + # been updated for this + + # get_potential_errors().add_potential_issue( + # MesidoAssetIssueType.ASSET_UPPER_LIMIT, + # asset.id, + # f"Pipe named {asset.name}: The upper limit of the pipe size has to be + # specified via a maximum value in a pipe diameter constraint." + # ) + # Raise the potential error here if applicable, with feedback to user + # Else a normal error exit might occer which will not give feedback to the user + # potential_error_to_error(self._error_type_check) + + return asset.attributes["diameter"].name + else: + logger.warning( + f"For pipe named {asset.name}, the pipe diameter constraint max value is " + "used for the pipe's diameter upper limit, instead of the pipe diameter " + "specified (if any) in the asset attribute." + ) + return pipe_constraint[0].maximum.name + else: + return asset.attributes["diameter"].name + def __override_pipe_classes_dicts( self, asset: Asset, @@ -275,17 +366,29 @@ def __override_pipe_classes_dicts( c = override_classes[p] = [] c.append(no_pipe_class) - min_size = self.__minimum_pipe_size_name - min_size_idx = [idx for idx, pipe in enumerate(pipe_classes) if pipe.name == min_size] - assert len(min_size_idx) == 1 - min_size_idx = min_size_idx[0] - - max_size = asset.attributes["diameter"].name + max_size_name = self._get_pipe_max_size_input(asset) - max_size_idx = [idx for idx, pipe in enumerate(pipe_classes) if pipe.name == max_size] + max_size_idx = [ + idx for idx, pipe in enumerate(pipe_classes) if pipe.name == max_size_name + ] assert len(max_size_idx) == 1 max_size_idx = max_size_idx[0] + # Update the minimum pipe DN size if user specified limit is allowed + # TODO: in the future the lower limit will make use of PipeDiameterConstriant + if self._ESDLMixin__use_user_defined_minimum_pipe_size: + user_defined_lower_limit_dn_mm = 40.0 + pipe_dn_max_vs_min_ratio = 10.0 # Relates to area ratio of 100.0 + max_size_dn_mm = float(pipe_classes[max_size_idx].name.replace("DN", "")) + min_dn_by_factor_mm = max_size_dn_mm / pipe_dn_max_vs_min_ratio + min_dn_mm = max(min_dn_by_factor_mm, user_defined_lower_limit_dn_mm) + min_size_idx = self.find_index_of_pipe_or_next_up(pipe_classes, min_dn_mm) + else: # use default minimum pipe DN size + min_size_idx = self.find_index_of_pipe_or_next_up( + pipe_classes, float(self.__minimum_pipe_size_name.replace("DN", "")) + ) + assert min_size_idx + if max_size_idx < min_size_idx: logger.warning( f"{p} has an upper DN size smaller than the used minimum size " @@ -341,6 +444,11 @@ def override_pipe_classes(self) -> None: pipe_classes[i], investment_costs=pipe_diameter_cost_map[pipe_class.name], ) + if ( + not self._ESDLMixin__use_user_defined_minimum_pipe_size + and float(pipe_classes[i].name.replace("DN", "")) != 20.0 + ): + self._ESDLMixin__use_user_defined_minimum_pipe_size = True else: del pipe_classes[i] diff --git a/src/mesido/esdl/profile_parser.py b/src/mesido/esdl/profile_parser.py index 0c738508..8639e327 100644 --- a/src/mesido/esdl/profile_parser.py +++ b/src/mesido/esdl/profile_parser.py @@ -11,6 +11,7 @@ from esdl.units.conversion import ENERGY_IN_J, POWER_IN_W, convert_to_unit from mesido.esdl.common import Asset +from mesido.esdl.esdl_additional_vars_mixin import get_asset_contraints from mesido.potential_errors import MesidoAssetIssueType, get_potential_errors import numpy as np @@ -44,12 +45,15 @@ def __init__( self, energy_system: esdl.EnergySystem, file_path: Optional[Path], + use_esdl_ranged_contraint: bool = False, ): self._profiles: Dict[int, Dict[str, np.ndarray]] = defaultdict(dict) self._energy_system: esdl.EnergySystem = energy_system self._file_path: Optional[Path] = file_path self._reference_datetimes: Optional[pd.DatetimeIndex] = None + self._use_esdl_ranged_constraint = use_esdl_ranged_contraint + def read_profiles( self, io: DataStore, @@ -110,9 +114,41 @@ def read_profiles( for component_type, var_name in self.component_type_to_var_name_map.items(): for component in energy_system_components.get(component_type, []): profile = self._profiles[ensemble_member].get(component + var_name, None) - asset_power = esdl_assets[esdl_asset_names_to_ids[component]].attributes[ - "power" - ] + asset = esdl_assets[esdl_asset_names_to_ids[component]] + asset_state = asset.attributes["state"] + + asset_power = None + if asset_state == esdl.AssetStateEnum.ENABLED: + asset_power = asset.attributes["power"] + elif asset_state == esdl.AssetStateEnum.OPTIONAL: + if ( + self._energy_system.esdlVersion is not None + and self._energy_system.esdlVersion >= "v2602" + and self._use_esdl_ranged_constraint + # "v2602" contains the items needed for the ranged constraint + # implementation in MESIDO + ): + range_constraints, qty_range_constraints = get_asset_contraints( + self, asset, esdl.RangedConstraint + ) + if qty_range_constraints == 1: + asset_power = range_constraints[0].range.maxValue + else: + logger.error( + f"Asset named {asset.name}: The code currently only " + f"caters for 1 RangedConstraint but {qty_range_constraints} " + "have been specified" + ) + sys.exit(1) + else: + asset_power = asset.attributes["power"] + else: + logger.warning( + f"Read profiles: asset {component} has a state {asset_state.name} " + "and currently the code only caters for asset states ENABLED or " + "OPTIONAL" + ) + if profile is not None: values = profile else: @@ -220,17 +256,20 @@ class InfluxDBProfileReader(BaseProfileReader): esdl.esdl.GasProducer: ".maximum_gas_source", esdl.esdl.ResidualHeatSource: ".maximum_heat_source", esdl.esdl.GeothermalSource: ".maximum_heat_source", + esdl.esdl.PVInstallation: ".maximum_electricity_source", } def __init__( self, energy_system: esdl.EnergySystem, file_path: Optional[Path], + use_esdl_ranged_contraint: bool = False, database_credentials: Optional[Dict[str, Tuple[str, str]]] = None, ): super().__init__( energy_system=energy_system, file_path=file_path, + use_esdl_ranged_contraint=use_esdl_ranged_contraint, ) self._df = pd.DataFrame() self._database_credentials = ( @@ -330,7 +369,6 @@ def _load_profiles_from_source( var_base_name = asset.name if variable_suffix in [ self.asset_type_to_variable_name_conversion[esdl.esdl.GasProducer], - self.asset_type_to_variable_name_conversion[esdl.esdl.ElectricityProducer], ]: logger.error( f"Profiles for {var_base_name} from esdl has not been tested yet but only" @@ -346,7 +384,6 @@ def _load_profiles_from_source( var_base_name = asset.name if var_base_name in [ self.asset_type_to_variable_name_conversion[esdl.esdl.GasProducer], - self.asset_type_to_variable_name_conversion[esdl.esdl.ElectricityProducer], ]: logger.error(f"Profiles for {var_base_name} from esdl has not been tested yet") sys.exit(1) @@ -608,10 +645,12 @@ def __init__( self, energy_system: esdl.EnergySystem, file_path: Path, + use_esdl_ranged_contraint: bool = False, ): super().__init__( energy_system=energy_system, file_path=file_path, + use_esdl_ranged_contraint=use_esdl_ranged_contraint, ) def _load_profiles_from_source( diff --git a/src/mesido/potential_errors.py b/src/mesido/potential_errors.py index 5eb1aea0..4f992e4e 100644 --- a/src/mesido/potential_errors.py +++ b/src/mesido/potential_errors.py @@ -46,6 +46,7 @@ class MesidoAssetIssueType(Enum): ASSET_PROFILE_MULTIPLIER = "asset_profile.multiplier" COLD_ASSET_TYPE_NOT_SUPPORTED = "cold_asset.type" ELECT_ASSET_TYPE_NOT_SUPPORTED = "elect_asset.type" + ASSET_UPPER_LIMIT = "asset.upper_limit_size" class PotentialErrors: diff --git a/src/mesido/workflows/gas_elect_workflow.py b/src/mesido/workflows/gas_elect_workflow.py index f05ad27a..946aa1ca 100644 --- a/src/mesido/workflows/gas_elect_workflow.py +++ b/src/mesido/workflows/gas_elect_workflow.py @@ -2,8 +2,7 @@ import os from mesido.esdl.esdl_additional_vars_mixin import ESDLAdditionalVarsMixin -from mesido.esdl.esdl_mixin import DBAccessType -from mesido.esdl.esdl_mixin import ESDLMixin +from mesido.esdl.esdl_mixin import DBAccessType, ESDLMixin from mesido.head_loss_class import HeadLossOption from mesido.network_common import NetworkSettings from mesido.techno_economic_mixin import TechnoEconomicMixin diff --git a/src/mesido/workflows/grow_workflow.py b/src/mesido/workflows/grow_workflow.py index 23017677..38a8ed0d 100644 --- a/src/mesido/workflows/grow_workflow.py +++ b/src/mesido/workflows/grow_workflow.py @@ -1,3 +1,4 @@ +import gc import locale import logging import os @@ -6,8 +7,7 @@ from typing import Dict from mesido.esdl.esdl_additional_vars_mixin import ESDLAdditionalVarsMixin -from mesido.esdl.esdl_mixin import DBAccessType -from mesido.esdl.esdl_mixin import ESDLMixin +from mesido.esdl.esdl_mixin import DBAccessType, ESDLMixin from mesido.head_loss_class import HeadLossOption from mesido.potential_errors import reset_potential_errors from mesido.techno_economic_mixin import TechnoEconomicMixin @@ -847,6 +847,9 @@ def run_end_scenario_sizing( priorities_output = solution._priorities_output + del solution + gc.collect() + solution = run_optimization_problem_solver( end_scenario_problem_class, solver_class=solver_class, diff --git a/src/mesido/workflows/utils/error_types.py b/src/mesido/workflows/utils/error_types.py index 931cefa3..61099a93 100644 --- a/src/mesido/workflows/utils/error_types.py +++ b/src/mesido/workflows/utils/error_types.py @@ -49,6 +49,7 @@ def mesido_issue_type_gen_message(issue_type: MesidoAssetIssueType) -> str: MesidoAssetIssueType.ASSET_PROFILE_MULTIPLIER: "Incorrect asset profile multiplier", MesidoAssetIssueType.COLD_ASSET_TYPE_NOT_SUPPORTED: "Unsupported assets are being used.", MesidoAssetIssueType.ELECT_ASSET_TYPE_NOT_SUPPORTED: "Unsupported assets are being used.", + MesidoAssetIssueType.ASSET_UPPER_LIMIT: " Asset upper limit input incompatibility.", } return type_and_general_meassage[issue_type] @@ -77,6 +78,7 @@ def potential_error_to_error(network_check_type: Enum) -> None: MesidoAssetIssueType.ASSET_PROFILE_MULTIPLIER, MesidoAssetIssueType.COLD_ASSET_TYPE_NOT_SUPPORTED, MesidoAssetIssueType.ELECT_ASSET_TYPE_NOT_SUPPORTED, + MesidoAssetIssueType.ASSET_UPPER_LIMIT, ], NetworkErrors.HEAT_AND_COOL_NETWORK_ERRORS: [ MesidoAssetIssueType.HEAT_DEMAND_POWER, diff --git a/tests/models/esdl_constraints/model/testing_network_1_all_enabled.esdl b/tests/models/esdl_constraints/model/testing_network_1_all_enabled.esdl new file mode 100644 index 00000000..924822c4 --- /dev/null +++ b/tests/models/esdl_constraints/model/testing_network_1_all_enabled.esdl @@ -0,0 +1,329 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/esdl_constraints/model/testing_network_1_all_optional.esdl b/tests/models/esdl_constraints/model/testing_network_1_all_optional.esdl new file mode 100644 index 00000000..d405cc21 --- /dev/null +++ b/tests/models/esdl_constraints/model/testing_network_1_all_optional.esdl @@ -0,0 +1,333 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/esdl_constraints/model/testing_network_1_all_optional_excluding_pipes.esdl b/tests/models/esdl_constraints/model/testing_network_1_all_optional_excluding_pipes.esdl new file mode 100644 index 00000000..91cc01e8 --- /dev/null +++ b/tests/models/esdl_constraints/model/testing_network_1_all_optional_excluding_pipes.esdl @@ -0,0 +1,333 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/esdl_constraints/model/testing_network_1_supply_only.esdl b/tests/models/esdl_constraints/model/testing_network_1_supply_only.esdl new file mode 100644 index 00000000..0b02daad --- /dev/null +++ b/tests/models/esdl_constraints/model/testing_network_1_supply_only.esdl @@ -0,0 +1,210 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/esdl_constraints/src/run_1.py b/tests/models/esdl_constraints/src/run_1.py new file mode 100644 index 00000000..d864ee8f --- /dev/null +++ b/tests/models/esdl_constraints/src/run_1.py @@ -0,0 +1,31 @@ +from pathlib import Path + +from mesido.esdl.esdl_parser import ESDLFileParser +from mesido.workflows import EndScenarioSizingStaged + +# from rtctools.util import run_optimization_problem + + +class NetworkProblem(EndScenarioSizingStaged): + ... + # Nothing to be added for now + + +if __name__ == "__main__": + base_folder = Path(__file__).resolve().parent.parent + model_folder = base_folder / "model" + + problem = EndScenarioSizingStaged( + # base_folder=base_folder, + # model_folder=model_folder, + esdl_parser=ESDLFileParser, + esdl_file_name="testing_network_1.esdl", + ) + problem.pre() + + # For potential future + # This network has not been setup to be able to optimize yet + # solution = run_optimization_problem( + # problem, + # ) + # results = solution.extract_results() diff --git a/tests/models/test_case_small_network_ates_buffer_optional_assets/model/test_case_small_network_all_optional_pipe_catalog_pipe_DN.esdl b/tests/models/test_case_small_network_ates_buffer_optional_assets/model/test_case_small_network_all_optional_pipe_catalog_pipe_DN.esdl new file mode 100644 index 00000000..7491031a --- /dev/null +++ b/tests/models/test_case_small_network_ates_buffer_optional_assets/model/test_case_small_network_all_optional_pipe_catalog_pipe_DN.esdl @@ -0,0 +1,574 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/unit_cases_electricity/source_sink_cable/input/timeseries_with_pv.csv b/tests/models/unit_cases_electricity/source_sink_cable/input/timeseries_with_pv.csv new file mode 100644 index 00000000..7e5571c4 --- /dev/null +++ b/tests/models/unit_cases_electricity/source_sink_cable/input/timeseries_with_pv.csv @@ -0,0 +1,4 @@ +DateTime,ElectricityDemand_2af6,PV +1-1-2019 00:00,100,150 +1-1-2019 01:00,300,250 +1-1-2019 02:00,100,150 diff --git a/tests/models/unit_cases_electricity/source_sink_cable/model/pv_with_csv_profile.esdl b/tests/models/unit_cases_electricity/source_sink_cable/model/pv_with_csv_profile.esdl new file mode 100644 index 00000000..811a2ff5 --- /dev/null +++ b/tests/models/unit_cases_electricity/source_sink_cable/model/pv_with_csv_profile.esdl @@ -0,0 +1,37 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/tests/models/unit_cases_electricity/source_sink_cable/src/example.py b/tests/models/unit_cases_electricity/source_sink_cable/src/example.py index 7883c213..eca3786e 100644 --- a/tests/models/unit_cases_electricity/source_sink_cable/src/example.py +++ b/tests/models/unit_cases_electricity/source_sink_cable/src/example.py @@ -37,15 +37,36 @@ def function( return optimization_problem.state(self.state) +class MinimizeElecProduction(Goal): + priority = 3 + + order = 1 + + def function(self, optimization_problem, ensemble_member): + sum_ = 0 + for source in optimization_problem.energy_system_components.get("electricity_source", []): + if source not in optimization_problem.energy_system_components.get("solar_pv", []): + sum_ += optimization_problem.state(f"{source}.Electricity_source") + return sum_ + + +class MinimizeElecProductionSize(Goal): + priority = 2 + + order = 1 + + def __init__(self, source): + self.source = source + self.target_max = 0.0 + self.function_range = (0.0, 2.0 * 1e6) + self.function_nominal = 1e3 + + def function(self, optimization_problem, ensemble_member): + return optimization_problem.extra_variable(f"{self.source}__max_size", ensemble_member) + + class _GoalsAndOptions: def path_goals(self): - """ - Add goal to meet the specified power demands in the electricity network. - - Returns - ------- - Extended goals list. - """ goals = super().path_goals().copy() for demand in self.energy_system_components["electricity_demand"]: @@ -63,6 +84,52 @@ def energy_system_options(self): return options +class ElectricityProblemPV( + _GoalsAndOptions, + TechnoEconomicMixin, + LinearizedOrderGoalProgrammingMixin, + GoalProgrammingMixin, + ESDLMixin, + CollocatedIntegratedOptimizationProblem, +): + """ + Problem to check the behaviour of a simple source, cable, demand network. + """ + + def path_goals(self): + goals = super().path_goals().copy() + goals.append(MinimizeElecProduction()) + return goals + + def goals(self): + """ + Add goal to minimize max_size of electricity producers while ensuring + that they are equal to each other. + + Returns + ------- + Extended goals list. + """ + goals = super().goals().copy() + for source in self.energy_system_components["electricity_source"]: + goals.append(MinimizeElecProductionSize(source=source)) + + return goals + + def constraints(self, ensemble_member): + constraints = super().constraints(ensemble_member) + elec_prod_size = self.extra_variable("ElectricityProducer_edde__max_size", ensemble_member) + pv_size = self.extra_variable("PV__max_size", ensemble_member) + nom = self.variable_nominal("PV__max_size") + constraints.append(((elec_prod_size - pv_size) / nom, 0.0, 0.0)) + return constraints + + def energy_system_options(self): + options = super().energy_system_options() + options["include_electric_cable_power_loss"] = False + return options + + class ElectricityProblem( _GoalsAndOptions, TechnoEconomicMixin, @@ -98,13 +165,6 @@ def read(self): self.set_timeseries(f"{d}.target_electricity_demand", new_timeseries) def path_goals(self): - """ - Modified targets for the demand matching goal to push up the current in the system. - - Returns - ------- - list with goals. - """ goals = super().path_goals().copy() for demand in self.energy_system_components["electricity_demand"]: diff --git a/tests/test_electric_source_sink.py b/tests/test_electric_source_sink.py index f2c0e1fa..d6c312c0 100644 --- a/tests/test_electric_source_sink.py +++ b/tests/test_electric_source_sink.py @@ -7,7 +7,7 @@ import numpy as np -from utils_tests import electric_power_conservation_test +from utils_tests import demand_matching_test, electric_power_conservation_test # TODO: still have to make test where elecitricity direction is switched: # e.g. 2 nodes, with at each node a producer and consumer, first one node medium demand, second @@ -15,6 +15,63 @@ class TestMILPElectricSourceSink(TestCase): + + def test_source_sink_pv_csv_profile(self): + """ + Tests for an electricity network that consist out of 1 PV and 1 Electricity Producer + as electricity sources, a cable and a sink. PV has profile that is read from input csv. + Electricity Producer has no profile constraint. + We check that scaled PV constraint profile is upper bound of PV production + and profile constraint scaling functionality in asset sizing works. + + Checks: + - Check demand matching and energy conservation + - Check profile constraint scaling functionality works. + """ + + import models.unit_cases_electricity.source_sink_cable.src.example as example + from models.unit_cases_electricity.source_sink_cable.src.example import ElectricityProblemPV + + base_folder = Path(example.__file__).resolve().parent.parent + + solution = run_esdl_mesido_optimization( + ElectricityProblemPV, + base_folder=base_folder, + esdl_file_name="pv_with_csv_profile.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="timeseries_with_pv.csv", + ) + results = solution.extract_results() + + tol = 1e-6 + + demand_matching_test(solution, results) + electric_power_conservation_test(solution, results) + + # Test that scaled PV constraint profile is upper bound of + # PV production and profile constraint scaling functionality + # in asset sizing works + profile_non_scaled = solution.get_timeseries("PV.maximum_electricity_source").values + max_profile_non_scaled = max(profile_non_scaled) + profile_scaled = profile_non_scaled / max_profile_non_scaled + + # Check that realized PV profile is scaled version of PV profile constraint + np.testing.assert_array_less( + results["PV.Electricity_source"], profile_scaled * results["PV__max_size"] + tol + ) + + # Check that maximum of realized PV profile smaller than maximum of PV profile constraint + np.testing.assert_array_less( + results["PV__max_size"], + max(solution.get_timeseries("PV.maximum_electricity_source").values) + tol, + ) + + # Check electricity producers max sizes are equal + np.testing.assert_allclose( + results["PV__max_size"], results["ElectricityProducer_edde__max_size"] + ) + def test_source_sink(self): """ Tests for an electricity network that consist out of a source, a cable and a sink. diff --git a/tests/test_end_scenario_sizing.py b/tests/test_end_scenario_sizing.py index 719e8d6a..538e8901 100644 --- a/tests/test_end_scenario_sizing.py +++ b/tests/test_end_scenario_sizing.py @@ -410,6 +410,62 @@ def __init__(self, *args, **kwargs): abs(results[f"{pipe}.dH"][ii]) + tol, ) + def test_end_scenario_sizing_pipe_catalog_lower_pipe_dn(self): + """ + Check that the temporary implementation of allowing the reduction of the lower pipe limit + via measures + """ + + import models.test_case_small_network_ates_buffer_optional_assets.src.run_ates as run_ates + + base_folder = Path(run_ates.__file__).resolve().parent.parent + + class ReducedDemandProblem(EndScenarioSizing): + def read(self): + super().read() + + # Reduce the heating demand to allwo usage of smaller pipes + for d in self.energy_system_components["heat_demand"]: + target = self.get_timeseries(f"{d}.target_heat_demand") + target.values[:] = target.values[:] / 1.0e2 + + self.io.set_timeseries( + f"{d}.target_heat_demand", + self.io._DataStore__timeseries_datetimes, + target.values, + 0, + ) + + solution = run_end_scenario_sizing( + ReducedDemandProblem, + base_folder=base_folder, + esdl_file_name="test_case_small_network_all_optional_pipe_catalog_pipe_DN.esdl", + esdl_parser=ESDLFileParser, + profile_reader=ProfileReaderFromFile, + input_timeseries_file="Warmte_test.csv", + error_type_check=NetworkErrors.NO_POTENTIAL_ERRORS_CHECK, + ) + + results = solution.extract_results() + parameters = solution.parameters(0) + + # Check min DN in the available pipe classes and the resulting pipe size + pipes_to_check = ["Pipe2", "Pipe2_ret", "Pipe4", "Pipe4_ret"] # user input minimum DN + for pipe in solution.energy_system_components.get("heat_pipe"): + available_pipe_classes = solution.pipe_classes(pipe) + + if pipe in pipes_to_check: + np.testing.assert_equal(available_pipe_classes[0].name == "DN40", True) + np.testing.assert_equal(parameters[f"{pipe}.diameter"], 0.0431) + np.testing.assert_equal(results[f"{pipe}__hn_diameter"], 0.0431) + elif pipe not in ["Pipe1", "Pipe1_ret"]: # this Pipe1 is not placed + if available_pipe_classes[0].name == "None": + np.testing.assert_equal(available_pipe_classes[1].name == "DN50", True) + else: + np.testing.assert_equal(available_pipe_classes[0].name == "DN50", True) + np.testing.assert_allclose(parameters[f"{pipe}.diameter"], 0.0545, atol=1e-6) + np.testing.assert_allclose(results[f"{pipe}__hn_diameter"], 0.0545, atol=1e-6) + def test_end_scenario_sizing_pipe_catalog(self): """ Checks that the problem uses different pipe costs when a pipe catalog is provided along @@ -657,6 +713,7 @@ def calculate_objective_value_end_scenario_sizing_all_optional(self, solution): a.test_end_scenario_sizing_staged() a.test_end_scenario_sizing_discounted() a.test_end_scenario_sizing_head_loss() + a.test_end_scenario_sizing_pipe_catalog_lower_pipe_dn() a.test_end_scenario_sizing_pipe_catalog() a.test_end_scenario_sizing_pipe_catalog_with_templates() print("Execution time: " + time.strftime("%M:%S", time.gmtime(time.time() - start_time))) diff --git a/tests/test_esdl_constraints.py b/tests/test_esdl_constraints.py new file mode 100644 index 00000000..08ce0e74 --- /dev/null +++ b/tests/test_esdl_constraints.py @@ -0,0 +1,295 @@ +from pathlib import Path +from unittest import TestCase + +from mesido.esdl.esdl_parser import ESDLFileParser + +import numpy as np + + +class TestRangedConstraints(TestCase): + """ + Check that RangedConstraint are correctly handled. + + """ + + def test_asset_max_capacity(self): + """ + Test that RangeConstraints or asset attribute input values are used as upper limits for + assets variables, depending if the asset is ENABLED or OPTIONAL. + OPTIONAL: use RangeConstraint specified at the asset, instead of the asset input attribute + ENABLED: use asset attribute input values + + In the esdl files the asset input value (power, volume or aggregation) count has been + specified to be smaller than the values specified in the RangedCosntraints + + Notes regarding the esdl files: + - testing_network_1_all_enabled: All the assets are enabled. + - testing_network_1_all_optional_excluding_pipes: + - All assets OPTIONAL and pipes are ENABLED + - Pipe 1 and 2 have PipeDiameterConstraint values + - All pipes (including pipe 1 & 2) have a Diameter specified that is smaller than + PipeDiameterConstraint value + - testing_network_1_all_optional: + - All the assets, including the pipes are OPTIONAL + - Pipe 1 & 2 have a PipeDiameterConstraint values. + - All the pipes (including pipe 1 & 2) have a Diameter attribute value specified. The + PipeDiameterConstraint value is larger than the pipe Diameter specified. + + Checks: + * Heat assets (excluding the pipes) + - Power and capacity: + - OPTIONAL assets: Check that the RangedContraint value is used for max size, and that + this values is larger than the asset attribute input value + - ENABLED assets: Check that the asset attribute input value is used for max size, and + that this values is smaller than the RangedContraint value + - Aggregation count assets: + - OPTIONAL assets: Check that the RangedContraint value is used for aggregation count + and that this values is larger than the asset attribute input value + - ENABLED assets: Check that the asset attribute input value is used for aggregation + count, and that this values is smaller than the RangedContraint value + - Volume assets: + - OPTIONAL assets: Check that the RangedContraint value is used for volume and the max + capacity of a buffer tank, and that this volume is larger than the asset attribute + input value + - ENABLED assets: # Check that the asset attribute input value is used for volume and + the max capacity of a buffer tank, and that this volume is smaler than the + RangedContraint value + - All OPTIONAL assets: Check the range contraint name + + * Pipes: + - Pipe DN size: + - OPTIONAL with PipeDiameterConstraint value and Diameter specified: Check that the + available pipe classes go up to PipeDiameterConstraint value and that this value is + larger than the Diameter specified + - OPTIONAL with only a Diameter specified: Check that the available pipe classes go up + to the Diameter specified + - ENABLED pipes: Check that the pipe has a diameter as specified in the attribute + Diameter value + + """ + import models.esdl_constraints.src.run_1 as run_1 + from models.esdl_constraints.src.run_1 import NetworkProblem + + base_folder = Path(run_1.__file__).resolve().parent.parent + model_folder = base_folder / "model" + + kwargs = {"use_esdl_ranged_constraint": True} + + # the esdl version must be >= "v2602" + esdl_files = { + "testing_network_1_all_enabled.esdl", + "testing_network_1_all_optional_excluding_pipes.esdl", + "testing_network_1_all_optional.esdl", + } + + for esdl_file in esdl_files: + problem = NetworkProblem( + base_folder=base_folder, + model_folder=model_folder, + esdl_parser=ESDLFileParser, + esdl_file_name=esdl_file, + **kwargs, + ) + + problem.pre() + + all_optional = False + if "all_optional" in esdl_file: + all_optional = True + + # ------------------------------------------------------------------------------------- + # Check the assets exlcuding the pipes + + asset_tested = { + # assets where power or capcicity is used for upper limit of size + "power_assets": [ + "HeatPump_6e86", + "HeatProducer_1c31", + "GasHeater_11d9", + "ElectricBoiler_ac4c", + "HeatExchange_1cf5", + ], + # assets where aggregation count is used for upper limit + "aggregation_assets": ["GeothermalSource_0a38", "ATES_49ac"], + # assets where volume is used for upper limit of size + "volume_assets": ["HeatStorage_50b6"], + } + + for key, asset_names in asset_tested.items(): + + for asset_name in asset_names: + esdl_asset = problem._esdl_assets[ + problem.esdl_asset_name_to_id_map[f"{asset_name}"] + ] + range_contraint_max = None + asset_input_max = None + + if "power_assets" == key: + power_or_capacity = "" + if asset_name == "HeatExchange_1cf5": + power_or_capacity = "capacity" + else: + power_or_capacity = "power" + + range_contraint_max = esdl_asset.attributes["constraint"][0].range.maxValue + asset_input_max = esdl_asset.attributes[power_or_capacity] + + if all_optional: + # Check that the RangedContraint value is used for max size, and that + # this values is larger than the asset attribute input value + np.testing.assert_equal( + problem.bounds()[f"{asset_name}__max_size"][1], + range_contraint_max, + ) + np.testing.assert_array_less( + asset_input_max, + problem.bounds()[f"{asset_name}__max_size"][1], + ) + # Check the range contraint name + np.testing.assert_equal( + esdl_asset.attributes["constraint"][0].attributeReference, + power_or_capacity, + ) + else: + # Check that the asset attribute input value is used for max size, and + # that this values is smaller than the RangedContraint value + np.testing.assert_equal( + problem.bounds()[f"{asset_name}__max_size"][1], + asset_input_max, + ) + np.testing.assert_array_less( + problem.bounds()[f"{asset_name}__max_size"][1], + range_contraint_max, + ) + + asset_variable = "" + if asset_name == "HeatPump_6e86": + asset_variable = "Secondary_heat" + else: + asset_variable = "Heat_flow" + # Check that the correct value is used for heat flow variable bound + np.testing.assert_equal( + problem.bounds()[f"{asset_name}__max_size"][1], + problem.bounds()[f"{asset_name}.{asset_variable}"][1], + ) + + if "aggregation_assets" == key: + + range_contraint_max = esdl_asset.attributes["constraint"][0].range.maxValue + asset_input_max = esdl_asset.attributes["aggregationCount"] + + if all_optional: + # Check that the RangedContraint value is used for aggregation count, + # and that this values is larger than the asset attribute input value + np.testing.assert_equal( + problem.bounds()[f"{asset_name}_aggregation_count"][1], + range_contraint_max, + ) + np.testing.assert_array_less( + asset_input_max, + problem.bounds()[f"{asset_name}_aggregation_count"][1], + ) + else: + # Check that the asset attribute input value is used for aggregation + # count, and that this values is smaller than the RangedContraint value + np.testing.assert_equal( + problem.bounds()[f"{asset_name}_aggregation_count"][1], + asset_input_max, + ) + np.testing.assert_array_less( + problem.bounds()[f"{asset_name}_aggregation_count"][1], + range_contraint_max, + ) + # Check the range contraint name + np.testing.assert_equal( + esdl_asset.attributes["constraint"][0].attributeReference, + "aggregationCount", + ) + + if "volume_assets" == key: + range_contraint_max = esdl_asset.attributes["constraint"][0].range.maxValue + asset_input_max = esdl_asset.attributes["volume"] + + if all_optional: + # Check that the RangedContraint value is used for volume and the max + # capacity of a buffer tank, and that this volume is larger than the + # asset attribute input value + asset_max_capacity = ( + problem.parameters(0)[f"{asset_name}.rho"] + * range_contraint_max + * problem.parameters(0)[f"{asset_name}.cp"] + * problem.parameters(0)[f"{asset_name}.dT"] + ) + np.testing.assert_equal( + asset_max_capacity, + problem.bounds()[f"{asset_name}__max_size"][1], + ) + np.testing.assert_array_less( + asset_input_max, + problem.parameters(0)[f"{asset_name}.volume"], + ) + np.testing.assert_almost_equal( + problem.parameters(0)[f"{asset_name}.volume"], + range_contraint_max, + ) + # Check the range contraint name + np.testing.assert_equal( + esdl_asset.attributes["constraint"][0].attributeReference, + "volume", + ) + else: + # Check that the asset attribute input value is used for volume and the + # max capacity of a buffer tank, and that this volume is smaler than the + # RangedContraint + asset_max_capacity = ( + problem.parameters(0)[f"{asset_name}.rho"] + * asset_input_max + * problem.parameters(0)[f"{asset_name}.cp"] + * problem.parameters(0)[f"{asset_name}.dT"] + ) + np.testing.assert_equal( + asset_max_capacity, + problem.bounds()[f"{asset_name}__max_size"][1], + ) + np.testing.assert_array_less( + problem.parameters(0)[f"{asset_name}.volume"], + range_contraint_max, + ) + np.testing.assert_almost_equal( + problem.parameters(0)[f"{asset_name}.volume"], + asset_input_max, + ) + + # ------------------------------------------------------------------------------------- + # Check the pipes with and without PipeDiameterConstraint + pipes_with_ranged_constraint = [""] + + if esdl_file == "testing_network_1_all_optional.esdl": + pipes_with_ranged_constraint = ["Pipe1", "Pipe1_ret", "Pipe2", "Pipe2_ret"] + + for pipe_name in problem.energy_system_components.get("heat_pipe", []): + esdl_asset = problem._esdl_assets[ + problem.esdl_asset_name_to_id_map[f"{pipe_name}"] + ] + + pipe_classes = list(problem._heat_pipe_topo_pipe_class_map[f"{pipe_name}"]) + pipe_input_size = esdl_asset.attributes["diameter"].name + + if pipe_name in pipes_with_ranged_constraint: + range_contraint_max = esdl_asset.attributes["constraint"][0].maximum.name + np.testing.assert_equal(pipe_classes[-1].name, range_contraint_max) + np.testing.assert_array_less( + int(pipe_input_size.replace("DN", "")), + int(range_contraint_max.replace("DN", "")), + ) + else: + np.testing.assert_equal(pipe_classes[-1].name, pipe_input_size) + elif esdl_file == "testing_network_1_all_optional_excluding_pipes.esdl": + for pipe_name in problem.energy_system_components.get("heat_pipe", []): + (problem.parameters(0)[f"{pipe_name}.diameter"], 0.695) # DN700 + + +if __name__ == "__main__": + + a = TestRangedConstraints() + a.test_asset_max_capacity() diff --git a/tests/test_potential_errors.py b/tests/test_potential_errors.py index de9795e2..7dcf3bb7 100644 --- a/tests/test_potential_errors.py +++ b/tests/test_potential_errors.py @@ -25,10 +25,12 @@ def __init__( self, energy_system: esdl.EnergySystem, file_path: Optional[Path], + use_esdl_ranged_contraint: bool = False, ): super().__init__( energy_system=energy_system, file_path=file_path, + use_esdl_ranged_contraint=use_esdl_ranged_contraint, ) self._loaded_profiles = pd.read_csv( file_path, diff --git a/tests/test_producer_profiles.py b/tests/test_producer_profiles.py index ea71b688..62e44d9a 100644 --- a/tests/test_producer_profiles.py +++ b/tests/test_producer_profiles.py @@ -217,6 +217,6 @@ def test_max_producer_esdl_scaled_profile(self): if __name__ == "__main__": a = TestProducerMaxProfile() - # a.test_max_producer_scaled_profile() + a.test_max_producer_scaled_profile() a.test_max_producer_esdl_unscaled_profile() - # a.test_max_producer_esdl_scaled_profile() + a.test_max_producer_esdl_scaled_profile() diff --git a/tests/test_profile_parsing.py b/tests/test_profile_parsing.py index 0184657c..b5479a5e 100644 --- a/tests/test_profile_parsing.py +++ b/tests/test_profile_parsing.py @@ -26,9 +26,15 @@ def __init__( self, energy_system: esdl.EnergySystem, file_path: Optional[Path], + use_esdl_ranged_contraint: bool, database_credentials: Optional[Dict[str, Tuple[str, str]]] = None, ): - super().__init__(energy_system, file_path, database_credentials=database_credentials) + super().__init__( + energy_system, + file_path, + use_esdl_ranged_contraint=use_esdl_ranged_contraint, + database_credentials=database_credentials, + ) self._loaded_profiles = pd.read_csv( file_path, index_col="DateTime", From a25ff019ea66c6eba01b8e1cb07342f6fa2d3972 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 15:09:37 +0100 Subject: [PATCH 21/23] Linting, formatting. --- src/mesido/esdl/esdl_model_base.py | 3 ++- src/mesido/financial_mixin.py | 6 ++++-- tests/test_multicommodity.py | 23 +++++++++++------------ 3 files changed, 17 insertions(+), 15 deletions(-) diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index 84018874..918e5e78 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -227,7 +227,8 @@ def __set_primary_secondary_heat_ports(): f"Heat out_ports " ) elif ( - asset.asset_type == "ElectricBoiler" or asset.asset_type == "GeothermalSource" + asset.asset_type == "ElectricBoiler" + or asset.asset_type == "GeothermalSource" and len(asset.out_ports) == 1 and len(asset.in_ports) == 2 ): diff --git a/src/mesido/financial_mixin.py b/src/mesido/financial_mixin.py index 09783d79..ac91ef31 100644 --- a/src/mesido/financial_mixin.py +++ b/src/mesido/financial_mixin.py @@ -1549,9 +1549,11 @@ def __annualized_capex_constraints(self, ensemble_member): # Input is assumed as as annual percentage discount_percentage = parameters[f"{asset_name}.discount_rate"] if np.isnan(asset_life_years) or np.isnan(discount_percentage): - logger.warning(f"Annualized cost cannot be computed for \ + logger.warning( + f"Annualized cost cannot be computed for \ {asset_name} since technical_life \ - or discount_rate are not set.") + or discount_rate are not set." + ) continue symbol_name = self._annualized_capex_var_map[asset_name] diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index 79cb114c..38b3b2d9 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -401,8 +401,10 @@ def test_geothermal_source_elec(self): # Equations check np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) - np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], - results["GeothermalSource_a77b.ElectricityIn.Power"]) + np.testing.assert_allclose( + results["ElectricityProducer_4dde.ElectricityOut.Power"], + results["GeothermalSource_a77b.ElectricityIn.Power"], + ) np.testing.assert_allclose( parameters["GeothermalSource_a77b.cop"] * results["GeothermalSource_a77b.Power_elec"], results["GeothermalSource_a77b.Heat_source"], @@ -422,7 +424,7 @@ def test_geothermal_source_elec(self): / parameters["GeothermalSource_a77b.cop"], results["GeothermalSource_a77b__variable_operational_cost"], ) - + def test_geothermal_source_elec_no_cop(self): """ This tests checks the electric geothermal producer asset when no cop is provided. @@ -453,10 +455,11 @@ def test_geothermal_source_elec_no_cop(self): # Equations check np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) - np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], - results["GeothermalSource_a77b.ElectricityIn.Power"]) - np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], - 0.0) + np.testing.assert_allclose( + results["ElectricityProducer_4dde.ElectricityOut.Power"], + results["GeothermalSource_a77b.ElectricityIn.Power"], + ) + np.testing.assert_allclose(results["ElectricityProducer_4dde.ElectricityOut.Power"], 0.0) # Test electricity port np.testing.assert_allclose( @@ -465,8 +468,4 @@ def test_geothermal_source_elec_no_cop(self): ) # Variable operational cost check. - np.testing.assert_allclose( - 0.0, - results["GeothermalSource_a77b__variable_operational_cost"] - ) - + np.testing.assert_allclose(0.0, results["GeothermalSource_a77b__variable_operational_cost"]) From 855c876789524bf0f7c45c11581443aa2ed75254 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 15:14:53 +0100 Subject: [PATCH 22/23] Formatting. --- src/mesido/financial_mixin.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/mesido/financial_mixin.py b/src/mesido/financial_mixin.py index ac91ef31..09783d79 100644 --- a/src/mesido/financial_mixin.py +++ b/src/mesido/financial_mixin.py @@ -1549,11 +1549,9 @@ def __annualized_capex_constraints(self, ensemble_member): # Input is assumed as as annual percentage discount_percentage = parameters[f"{asset_name}.discount_rate"] if np.isnan(asset_life_years) or np.isnan(discount_percentage): - logger.warning( - f"Annualized cost cannot be computed for \ + logger.warning(f"Annualized cost cannot be computed for \ {asset_name} since technical_life \ - or discount_rate are not set." - ) + or discount_rate are not set.") continue symbol_name = self._annualized_capex_var_map[asset_name] From 7b7c45a243a65634a764e7ee0c52955860a27f23 Mon Sep 17 00:00:00 2001 From: Santiago Patterson Date: Thu, 12 Mar 2026 17:06:27 +0100 Subject: [PATCH 23/23] PR comments. --- src/mesido/esdl/esdl_model_base.py | 15 +++++++-------- .../milp/heat/geothermal_source_elec.py | 1 - tests/test_multicommodity.py | 3 --- 3 files changed, 7 insertions(+), 12 deletions(-) diff --git a/src/mesido/esdl/esdl_model_base.py b/src/mesido/esdl/esdl_model_base.py index 918e5e78..83a96247 100644 --- a/src/mesido/esdl/esdl_model_base.py +++ b/src/mesido/esdl/esdl_model_base.py @@ -180,7 +180,7 @@ def __set_primary_secondary_heat_ports(): f"milp(4) and electricity (1) ports" ) elif ( - (asset.asset_type == "HeatPump" or asset.asset_type == "GeothermalSource") + (asset.asset_type == "HeatPump") and len(asset.out_ports) == 1 and len(asset.in_ports) in [1, 2] ): @@ -200,13 +200,12 @@ def __set_primary_secondary_heat_ports(): f"in port and 1 Heat out_ports " ) else: - if asset.asset_type == "HeatPump": - raise Exception( - f"{asset.name} has incorrect number of in/out ports. HeatPumps allow " - f"to have 1 in and 1 out port for air-water HP, 2 in ports and 2 out " - f"ports when modelling a water-water HP, or 3 in ports and 2 out ports " - f"when the electricity connection of the water-water HP is modelled." - ) + raise Exception( + f"{asset.name} has incorrect number of in/out ports. HeatPumps allow " + f"to have 1 in and 1 out port for air-water HP, 2 in ports and 2 out " + f"ports when modelling a water-water HP, or 3 in ports and 2 out ports " + f"when the electricity connection of the water-water HP is modelled." + ) elif ( asset.asset_type == "GasHeater" diff --git a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py index e5718333..8b9b6996 100644 --- a/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py +++ b/src/mesido/pycml/component_library/milp/heat/geothermal_source_elec.py @@ -28,7 +28,6 @@ def __init__(self, name, **modifiers): ) self.component_subtype = "geothermal_source_elec" - self.id_mapping_carrier = -1 self.min_voltage = nan self.add_variable(ElectricityPort, "ElectricityIn") diff --git a/tests/test_multicommodity.py b/tests/test_multicommodity.py index 38b3b2d9..bbb2013d 100644 --- a/tests/test_multicommodity.py +++ b/tests/test_multicommodity.py @@ -400,7 +400,6 @@ def test_geothermal_source_elec(self): electric_power_conservation_test(heat_problem, results) # Equations check - np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) np.testing.assert_allclose( results["ElectricityProducer_4dde.ElectricityOut.Power"], results["GeothermalSource_a77b.ElectricityIn.Power"], @@ -411,7 +410,6 @@ def test_geothermal_source_elec(self): ) # Test electricity port - np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.ElectricityIn.Power"]) np.testing.assert_allclose( results["GeothermalSource_a77b.ElectricityIn.Power"], results["GeothermalSource_a77b.Power_elec"], @@ -454,7 +452,6 @@ def test_geothermal_source_elec_no_cop(self): electric_power_conservation_test(heat_problem, results) # Equations check - np.testing.assert_array_less(0.0, results["GeothermalSource_a77b.Heat_source"]) np.testing.assert_allclose( results["ElectricityProducer_4dde.ElectricityOut.Power"], results["GeothermalSource_a77b.ElectricityIn.Power"],