diff --git a/src/mesido/heat_physics_mixin.py b/src/mesido/heat_physics_mixin.py
index 21b00cd93..8a625882b 100644
--- a/src/mesido/heat_physics_mixin.py
+++ b/src/mesido/heat_physics_mixin.py
@@ -445,14 +445,22 @@ def _get_min_bound(bound):
if len(temperature_regimes) == 0:
temperature = temperatures["temperature"]
self.__temperature_regime_var_bounds[temp_var_name] = (temperature, temperature)
- elif len(temperature_regimes) == 1:
- temperature = temperature_regimes[0]
- self.__temperature_regime_var_bounds[temp_var_name] = (temperature, temperature)
else:
- self.__temperature_regime_var_bounds[temp_var_name] = (
- min(temperature_regimes),
- max(temperature_regimes),
- )
+ if max(temperature_regimes) >= temperatures["temperature"]:
+ logger.error(
+ f"The temperature provided for carrier with name "
+ f"{temperatures['name']} and id {temperatures['id']} is "
+ f"smaller than the largest value in the temperature regime "
+ f"provided for this carrier, please update the esdl."
+ )
+ if len(temperature_regimes) == 1:
+ temperature = temperature_regimes[0]
+ self.__temperature_regime_var_bounds[temp_var_name] = (temperature, temperature)
+ else:
+ self.__temperature_regime_var_bounds[temp_var_name] = (
+ min(temperature_regimes),
+ max(temperature_regimes),
+ )
for temperature_regime in temperature_regimes:
carrier_selected_var = carrier_id_number_mapping + f"_{temperature_regime}"
@@ -1503,6 +1511,7 @@ def __source_heat_to_discharge_path_constraints(self, ensemble_member):
sup_carrier = parameters[f"{s}.T_supply_id"]
supply_temperatures = self.temperature_regimes(sup_carrier)
+
big_m = 2.0 * self.bounds()[f"{s}.HeatOut.Heat"][1]
big_m = (
big_m
@@ -1510,43 +1519,104 @@ def __source_heat_to_discharge_path_constraints(self, ensemble_member):
else 2.0 * self.bounds()[f"{s}.Heat_source"][1] * parameters[f"{s}.T_supply"] / dt
)
- if len(supply_temperatures) == 0:
- constraints.append(
- (
- (heat_out - discharge * cp * rho * parameters[f"{s}.T_supply"])
- / heat_nominal,
- 0.0,
- 0.0,
- )
- )
- else:
- for supply_temperature in supply_temperatures:
- sup_temperature_is_selected = self.state(f"{sup_carrier}_{supply_temperature}")
+ # Check to see if the out carrier has a temperature profile assigned to it.
+ temp_out_profile, _, _, _ = self.__get_out_port_carrier_temp_profile(
+ parameters, s, "heat_source"
+ )
+ if temp_out_profile is None:
+ if len(supply_temperatures) == 0:
constraints.append(
(
- (
- heat_out
- - discharge * cp * rho * supply_temperature
- + (1.0 - sup_temperature_is_selected) * big_m
- )
- / constraint_nominal,
+ (heat_out - discharge * cp * rho * parameters[f"{s}.T_supply"])
+ / heat_nominal,
+ 0.0,
0.0,
- np.inf,
)
)
- constraints.append(
- (
+ else:
+ for supply_temperature in supply_temperatures:
+ sup_temperature_is_selected = self.state(
+ f"{sup_carrier}_{supply_temperature}"
+ )
+
+ constraints.append(
(
- heat_out
- - discharge * cp * rho * supply_temperature
- - (1.0 - sup_temperature_is_selected) * big_m
+ (
+ heat_out
+ - discharge * cp * rho * supply_temperature
+ + (1.0 - sup_temperature_is_selected) * big_m
+ )
+ / constraint_nominal,
+ 0.0,
+ np.inf,
+ )
+ )
+ constraints.append(
+ (
+ (
+ heat_out
+ - discharge * cp * rho * supply_temperature
+ - (1.0 - sup_temperature_is_selected) * big_m
+ )
+ / constraint_nominal,
+ -np.inf,
+ 0.0,
)
- / constraint_nominal,
- -np.inf,
- 0.0,
)
+
+ return constraints
+
+ def __source_heat_to_discharge_variable_temp_constraints(self, ensemble_member):
+ """
+ Adds the same type of constraints to the source as
+ __source_heat_to_discharge_path_constraints for cases where the out carrier
+ has a prescribed temperature profile. An important difference is that these
+ are conventional constraints, since every timestep will have a specific value.
+ """
+
+ constraints = []
+ parameters = self.parameters(ensemble_member)
+
+ for s in self.energy_system_components.get("heat_source", []):
+ heat_nominal = parameters[f"{s}.Heat_nominal"]
+ cp = parameters[f"{s}.cp"]
+ rho = parameters[f"{s}.rho"]
+ dt = parameters[f"{s}.dT"]
+
+ big_m = 2.0 * self.bounds()[f"{s}.HeatOut.Heat"][1]
+ big_m = (
+ big_m
+ if big_m != np.inf
+ else 2.0 * self.bounds()[f"{s}.Heat_source"][1] * parameters[f"{s}.T_supply"] / dt
+ )
+
+ temp_out_profile, sup_carrier_name, temp_out_prof_start_idx, temp_out_prof_end_idx = (
+ self.__get_out_port_carrier_temp_profile(parameters, s, "heat_source")
+ )
+
+ if (
+ temp_out_profile is not None
+ ): # Case where the out carrier has a temp profile assigned to it.
+ heat_out_vector = self.__state_vector_scaled(f"{s}.HeatOut.Heat", ensemble_member)
+ discharge_vector = self.__state_vector_scaled(f"{s}.Q", ensemble_member)
+
+ constraints.append(
+ (
+ (
+ heat_out_vector
+ - discharge_vector
+ * cp
+ * rho
+ * temp_out_profile.values[
+ temp_out_prof_start_idx : temp_out_prof_end_idx + 1
+ ]
+ )
+ / heat_nominal,
+ 0.0,
+ 0.0,
)
+ )
return constraints
@@ -1617,6 +1687,56 @@ def __cold_demand_heat_to_discharge_path_constraints(self, ensemble_member):
return constraints
+ def __get_out_port_carrier_temp_profile(self, parameters, asset_name, asset_type):
+ """
+ This function finds the carrier lined to the asset's out port and grabs the
+ temperature profile assigned to it, if there is one assigned to it.
+ It returns the temperature as a timeseries, the name of the carrier, and the
+ start and end index of the temperature profile according to the problem's timeseries
+ (this last one only relevant for problems that are sliced).
+ """
+ # TODO: modify profile parser so the temperature profile is not called a price profile.
+ # TODO: modify this and the profile parser to incorporate the esdl option to
+ # use power instead of temperature.
+ # TODO: the start/end indices are needed for a very specific problem-times slicing case.
+ # Modify the way problems are sliced to remove this.
+
+ try:
+ carriers = self.esdl_carriers
+ carriers_ids = carriers.keys()
+ except AttributeError:
+ carriers = None
+ carriers_ids = []
+ sup_carrier_name = None
+ temp_out_profile = None
+ temp_out_prof_start_idx = None
+ temp_out_prof_end_idx = None
+ carrier_id_types = {"heat_source": ".T_supply_id", "heat_pipe": ".carrier_id"}
+ for carrier_id in carriers_ids:
+ if (
+ carriers[carrier_id]["id_number_mapping"]
+ == parameters[f"{asset_name}{carrier_id_types[asset_type]}"]
+ ):
+ sup_carrier_name = carriers[carrier_id]["name"]
+ try:
+ temp_out_profile = self.get_timeseries(f"{sup_carrier_name}.price_profile")
+ temp_out_prof_start_idx = int(
+ np.where(
+ self.get_timeseries(f"{sup_carrier_name}.price_profile").times
+ == self.times()[0]
+ )[0]
+ )
+ temp_out_prof_end_idx = int(
+ np.where(
+ self.get_timeseries(f"{sup_carrier_name}.price_profile").times
+ == self.times()[-1]
+ )[0]
+ )
+ except KeyError:
+ pass
+
+ return temp_out_profile, sup_carrier_name, temp_out_prof_start_idx, temp_out_prof_end_idx
+
def __pipe_heat_to_discharge_path_constraints(self, ensemble_member):
"""
This function adds constraints linking the flow to the thermal power at the pipe assets.
@@ -1668,7 +1788,13 @@ def __pipe_heat_to_discharge_path_constraints(self, ensemble_member):
temperatures = self.temperature_regimes(carrier)
for heat in [scaled_heat_in, scaled_heat_out]:
- if self.energy_system_options()["neglect_pipe_heat_losses"]:
+ temp_out_profile, _, _, _ = self.__get_out_port_carrier_temp_profile(
+ parameters, p, "heat_pipe"
+ )
+ if (
+ self.energy_system_options()["neglect_pipe_heat_losses"]
+ and temp_out_profile is None
+ ):
temp = parameters[f"{p}.temperature"]
if len(temperatures) == 0:
constraints.append(
@@ -1705,7 +1831,7 @@ def __pipe_heat_to_discharge_path_constraints(self, ensemble_member):
0.0,
)
)
- else:
+ elif not self.energy_system_options()["neglect_pipe_heat_losses"]:
# Note that during cold delivery the line can be colder than the ground
# temperature.
# In this case we have to bound the heat flowing in the line with the ground
@@ -3648,6 +3774,10 @@ def constraints(self, ensemble_member):
constraints.extend(self.__heat_matching_demand_insulation_constraints(ensemble_member))
constraints.extend(self.__ates_max_stored_heat_constriants(ensemble_member))
+ constraints.extend(
+ self.__source_heat_to_discharge_variable_temp_constraints(ensemble_member)
+ )
+
return constraints
def history(self, ensemble_member):
diff --git a/tests/models/source_pipe_sink/input/timeseries_prod_test.csv b/tests/models/source_pipe_sink/input/timeseries_prod_test.csv
new file mode 100644
index 000000000..ee57a729f
--- /dev/null
+++ b/tests/models/source_pipe_sink/input/timeseries_prod_test.csv
@@ -0,0 +1,46 @@
+DateTime,demand,elec,heat
+19-5-2013 22:00,150000,10,70
+19-5-2013 23:00,150000,10,72.0
+20-5-2013 00:00,150000,10,74
+20-5-2013 01:00,150000,10,76
+20-5-2013 02:00,150000,10,80
+20-5-2013 03:00,150000,10,82
+20-5-2013 04:00,150000,10,84
+20-5-2013 05:00,150000,10,86
+20-5-2013 06:00,150000,10,82
+20-5-2013 07:00,150000,10,80
+20-5-2013 08:00,150000,10,78
+20-5-2013 09:00,150000,10,76
+20-5-2013 10:00,150000,10,74
+20-5-2013 11:00,150000,10,70
+20-5-2013 12:00,150000,10,70
+20-5-2013 13:00,150000,10,70
+20-5-2013 14:00,100000,10,70
+20-5-2013 15:00,100000,10,70
+20-5-2013 16:00,100000,10,70
+20-5-2013 17:00,100000,10,70
+20-5-2013 18:00,100000,10,70
+20-5-2013 19:00,100000,10,70
+20-5-2013 20:00,100000,10,70
+20-5-2013 21:00,100000,10,70
+20-5-2013 22:00,100000,10,70
+20-5-2013 23:00,100000,10,70
+21-5-2013 00:00,100000,10,70
+21-5-2013 01:00,100000,10,70
+21-5-2013 02:00,100000,10,70
+21-5-2013 03:00,50000,10,70
+21-5-2013 04:00,50000,10,70
+21-5-2013 05:00,50000,10,70
+21-5-2013 06:00,50000,10,70
+21-5-2013 07:00,50000,10,70
+21-5-2013 08:00,50000,10,70
+21-5-2013 09:00,50000,10,70
+21-5-2013 10:00,50000,10,70
+21-5-2013 11:00,50000,10,70
+21-5-2013 12:00,50000,10,70
+21-5-2013 13:00,50000,10,70
+21-5-2013 14:00,50000,10,70
+21-5-2013 15:00,50000,10,70
+21-5-2013 16:00,50000,10,70
+21-5-2013 17:00,50000,10,70
+21-5-2013 18:00,50000,10,70
diff --git a/tests/models/source_pipe_sink/model/sourcesink_prof_test.esdl b/tests/models/source_pipe_sink/model/sourcesink_prof_test.esdl
new file mode 100644
index 000000000..cfdf8eb8e
--- /dev/null
+++ b/tests/models/source_pipe_sink/model/sourcesink_prof_test.esdl
@@ -0,0 +1,57 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/tests/models/source_pipe_sink/model/sourcesink_prof_test_2prod.esdl b/tests/models/source_pipe_sink/model/sourcesink_prof_test_2prod.esdl
new file mode 100644
index 000000000..74e38c648
--- /dev/null
+++ b/tests/models/source_pipe_sink/model/sourcesink_prof_test_2prod.esdl
@@ -0,0 +1,71 @@
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
diff --git a/tests/test_heat.py b/tests/test_heat.py
index 7f076bb81..316692a16 100644
--- a/tests/test_heat.py
+++ b/tests/test_heat.py
@@ -88,6 +88,78 @@ def energy_system_options(self):
energy_conservation_test(case, results)
heat_to_discharge_test(case, results)
+ def test_heat_prod_profile(self):
+ """
+ Test to ensure a carrier on the out port of a producer with a temperature profile
+ follows said profile. Two cases are run, the first one with just one producer
+ and one consumer. The second one has two producers connected in series. The main goal
+ of the second case is to ensure the system runs when one of the producers has a
+ prescribed out temperature equal to its in one (effectively bypassing it).
+
+ Checks:
+ - Checks that the temperatures coming out of the producer match the input profile.
+
+ """
+ import models.source_pipe_sink.src.double_pipe_heat as double_pipe_heat
+ from models.source_pipe_sink.src.double_pipe_heat import SourcePipeSink
+
+ class Model(SourcePipeSink):
+ def energy_system_options(self):
+ options = super().energy_system_options()
+ options["neglect_pipe_heat_losses"] = True
+
+ return options
+
+ base_folder = Path(double_pipe_heat.__file__).resolve().parent.parent
+
+ case = run_esdl_mesido_optimization(
+ Model,
+ base_folder=base_folder,
+ esdl_file_name="sourcesink_prof_test.esdl",
+ esdl_parser=ESDLFileParser,
+ profile_reader=ProfileReaderFromFile,
+ input_timeseries_file="timeseries_prod_test.csv",
+ )
+
+ results = case.extract_results()
+ parameters = case.parameters(0)
+
+ demand_matching_test(case, results)
+ energy_conservation_test(case, results)
+ heat_to_discharge_test(case, results)
+
+ heat_flow_out_pipe1 = results["Pipe1.HeatOut.Heat"]
+ vol_flow_pipe1 = results["Pipe1.HeatIn.Q"]
+ cp = parameters["Pipe1.cp"]
+ rho = parameters["Pipe1.rho"]
+ temp_pipe1 = heat_flow_out_pipe1 / (cp * rho * vol_flow_pipe1)
+ temp_input_prof = case.get_timeseries("heat.price_profile").values
+ np.testing.assert_array_almost_equal(temp_pipe1, temp_input_prof)
+
+ case = run_esdl_mesido_optimization(
+ Model,
+ base_folder=base_folder,
+ esdl_file_name="sourcesink_prof_test_2prod.esdl",
+ esdl_parser=ESDLFileParser,
+ profile_reader=ProfileReaderFromFile,
+ input_timeseries_file="timeseries_prod_test.csv",
+ )
+
+ results = case.extract_results()
+ parameters = case.parameters(0)
+
+ demand_matching_test(case, results)
+ energy_conservation_test(case, results)
+ heat_to_discharge_test(case, results)
+
+ heat_flow_out_pipe1 = results["Pipe1.HeatOut.Heat"]
+ vol_flow_pipe1 = results["Pipe1.HeatIn.Q"]
+ cp = parameters["Pipe1.cp"]
+ rho = parameters["Pipe1.rho"]
+ temp_pipe1 = heat_flow_out_pipe1 / (cp * rho * vol_flow_pipe1)
+ temp_input_prof = case.get_timeseries("heat.price_profile").values
+ np.testing.assert_array_almost_equal(temp_pipe1, temp_input_prof)
+
class TestMinMaxPressureOptions(TestCase):
import models.source_pipe_sink.src.double_pipe_heat as double_pipe_heat
@@ -374,6 +446,5 @@ def test_disconnected_pipe_darcy_weisbach(self):
if __name__ == "__main__":
-
a = TestHeat()
a.test_heat_loss()
diff --git a/tests/utils_tests.py b/tests/utils_tests.py
index 2577ae89e..86b6e946c 100644
--- a/tests/utils_tests.py
+++ b/tests/utils_tests.py
@@ -7,6 +7,36 @@
import numpy as np
+def __get_out_port_temp_profile(solution, asset_name, asset_type):
+ """
+ Returns the temperature profile specified (if any) for the carrier assigned to
+ the carrier on the out port.
+ """
+
+ parameters = solution.parameters(0)
+ try:
+ carriers = solution.esdl_carriers
+ carriers_ids = carriers.keys()
+ except AttributeError:
+ carriers = None
+ carriers_ids = []
+ sup_carrier_name = None
+ temp_out_profile = None
+ carrier_id_types = {"heat_source": ".T_supply_id", "heat_pipe": ".carrier_id"}
+ for carrier_id in carriers_ids:
+ if (
+ carriers[carrier_id]["id_number_mapping"]
+ == parameters[f"{asset_name}{carrier_id_types[asset_type]}"]
+ ):
+ sup_carrier_name = carriers[carrier_id]["name"]
+ try:
+ temp_out_profile = solution.get_timeseries(f"{sup_carrier_name}.price_profile")
+ except KeyError:
+ pass
+
+ return temp_out_profile
+
+
def feasibility_test(solution):
feasibility = solution.solver_stats["return_status"]
@@ -159,6 +189,7 @@ def heat_to_discharge_test(solution, results):
results[f"{d}.HeatOut.Heat"], results[f"{d}.Q"] * rho * cp * supply_t
)
+ supply_temp_profiles = []
for d in solution.energy_system_components.get("heat_source", []):
cp = solution.parameters(0)[f"{d}.cp"]
rho = solution.parameters(0)[f"{d}.rho"]
@@ -171,7 +202,11 @@ def heat_to_discharge_test(solution, results):
# supply_t = solution.parameters(0)[f"{d}.T_supply"]
# return_t = solution.parameters(0)[f"{d}.T_return"]
supply_t, return_t, dt = _get_component_temperatures(solution, results, d)
-
+ temp_profile = __get_out_port_temp_profile(solution, d, "heat_source")
+ if temp_profile is not None:
+ supply_t = temp_profile.values
+ supply_temp_profiles.append(temp_profile.values)
+ print(d, max(abs(results[f"{d}.HeatOut.Heat"] - results[f"{d}.Q"] * rho * cp * supply_t)))
np.testing.assert_allclose(
results[f"{d}.HeatOut.Heat"],
results[f"{d}.Q"] * rho * cp * supply_t,
@@ -285,7 +320,9 @@ def heat_to_discharge_test(solution, results):
rho = solution.parameters(0)[f"{p}.rho"]
carrier_id = solution.parameters(0)[f"{p}.carrier_id"]
indices = results[f"{p}.Q"] > 0
- if f"{carrier_id}_temperature" in results.keys():
+ if supply_temp_profiles:
+ temperature = max(temp for prof in supply_temp_profiles for temp in prof)
+ elif f"{carrier_id}_temperature" in results.keys():
temperature = np.clip(
results[f"{carrier_id}_temperature"][indices],
solution.parameters(0)[f"{p}.T_ground"],