Solar ice example with two storages

data/examples/run_example_config.py

@@ -1,6 +1,6 @@
 import pandas as pd
 import os
-import pathlib as pl
+import pathlib as _pl
 
 from optihood.IO.writers import ScenarioFileWriterExcel
 from optihood.energy_network import EnergyNetworkGroup as EnergyNetwork
@@ -9,26 +9,27 @@
 
     optimization_type = "cost"
     scenario = "HP_GS_PV"
-    main_result_folder = "..\\results"
     # set a time period for the optimization problem
     timePeriod = pd.date_range("2018-01-01 00:00:00", "2018-01-31 23:00:00", freq="60min")
 
     # define paths for input and result files
-    inputFilePath = pl.Path(f"..\\configs\\basic_example_config")
+    curDir = _pl.Path(__file__).resolve().parent
+    inputFilePath = curDir / ".." / "configs" / "basic_example_config"
     configFileName = f"scenario_{scenario}_group.ini"
     inputfileName = f"scenario_{scenario}.xls"  # excel file which would be created by createScenarioFile()
 
-    resultFilePath = f".\\{main_result_folder}\\{scenario}\\{optimization_type}\\group"
-    resultFileName = "scenario_HP_GS_PV.xls"
+    resultFilePath = curDir / ".." / "results"
+    resultFileName = "results_example_config.xlsx"
 
-    print("result file path: " + resultFilePath)
+    print("Scenario config file path: " + os.path.join(inputFilePath, configFileName))
+    print("Scenario excel file path: " + os.path.join(inputFilePath, inputfileName))
+    print("Result file path: " + os.path.join(resultFilePath, resultFileName))
 
     # initialize parameters
     numberOfBuildings = 4
     optimizationType = "costs"  # set as "env" for environmental optimization and "costs" for cost optimization
     mergeLinkBuses = True
     merge_buses = ["electricity", "space_heat", "domestic_hot_water"]
-    dispatchMode = True  # Set to True to run the optimization in dispatch mode
 
     # solver specific command line options
     optimizationOptions = {
@@ -58,7 +59,7 @@
     # create an energy network and set the network parameters from an excel file
     network = EnergyNetwork(timePeriod)
     network.setFromExcel(os.path.join(inputFilePath, inputfileName), numberOfBuildings, opt=optimizationType,
-                         mergeLinkBuses=mergeLinkBuses, mergeBuses=merge_buses, dispatchMode=dispatchMode)
+                         mergeLinkBuses=mergeLinkBuses, mergeBuses=merge_buses)
 
     # optimize the energy network
     limit, capacitiesTransformers, capacitiesStorages = network.optimize(solver='gurobi',

data/examples/solar_ice_hp.py

@@ -16,6 +16,7 @@
 
     resultFilePath = curDir / ".." / "results"
     resultFileName = "results_solar_ice_hp.xlsx"
+    iceStoreFileName = "ice_storage_params.csv"
 
     print("Scenario file path: " + os.path.join(inputFilePath, inputfileName))
     print("Result file path: " + os.path.join(resultFilePath, resultFileName))
@@ -69,3 +70,4 @@
     if not os.path.exists(resultFilePath):
         os.makedirs(resultFilePath)
     network.exportToExcel(os.path.join(resultFilePath, resultFileName), mergeLinkBuses=mergeLinkBuses)
+    network.exportIceStorageModelParams(os.path.join(resultFilePath, iceStoreFileName))

data/examples/solar_ice_hp_2storages.py

@@ -0,0 +1,78 @@
+import pandas as pd
+import os
+import pathlib as _pl
+
+from optihood.IO.writers import ScenarioFileWriterExcel
+from optihood.energy_network import EnergyNetworkGroup as EnergyNetwork
+
+if __name__ == '__main__':
+    """
+    This example is a modification of the solar_ice_hp example.
+    Instead of a multi-layer thermal storage, this example has separate space heating and DHW storages
+    Therefore, this example has only 2 temperature levels --> SH and DHW temperatures
+    """
+    # set a time period for the optimization problem
+    timePeriod = pd.date_range("2018-06-01 00:00:00", "2018-12-31 23:00:00", freq="60min")
+
+    # define paths for input and result files
+    curDir = _pl.Path(__file__).resolve().parent
+    inputFilePath = curDir / ".." / "excels" / "solar_ice_hp_2storages"
+    inputfileName = "scenario.xls"
+
+    resultFilePath = curDir / ".." / "results"
+    resultFileName = "results_solar_ice_hp2.xlsx"
+    iceStoreFileName = "ice_storage_params2.csv"
+
+    print("Scenario file path: " + os.path.join(inputFilePath, inputfileName))
+    print("Result file path: " + os.path.join(resultFilePath, resultFileName))
+
+    # initialize parameters
+    numberOfBuildings = 4
+    optimizationType = "costs"  # set as "env" for environmental optimization and "costs" for cost optimization
+    mergeLinkBuses = False
+    mergeBuses = [] # "electricity", "heatPumpInputBus"
+    dispatchMode = True  # Set to True to run the optimization in dispatch mode
+
+    # solver specific command line options
+    optimizationOptions = {
+        "gurobi": {
+            "BarConvTol": 0.5,
+            # The barrier solver terminates when the relative difference between the primal and dual objective values
+            # is less than the specified tolerance (with a GRB_OPTIMAL status)
+            "OptimalityTol": 1e-4,
+            # Reduced costs must all be smaller than OptimalityTol in the improving direction in order for a model to
+            # be declared optimal
+            "MIPGap": 1e-2,
+            # Relative Tolerance between the best integer objective and the objective of the best node remaining
+            "MIPFocus": 2
+            # 1 feasible solution quickly. 2 proving optimality. 3 if the best objective bound is moving very
+            # slowly/focus on the bound "Cutoff": #Indicates that you aren't interested in solutions whose objective
+            # values are worse than the specified value., could be dynamically be used in moo
+        }
+        # add the command line options of the solver here, For example to use CBC add a new item to the dictionary
+        # "cbc": {"tee": False}
+    }
+
+    # create an energy network and set the network parameters from an excel file
+    network = EnergyNetwork(timePeriod)
+    network.setFromExcel(os.path.join(inputFilePath, inputfileName), numberOfBuildings, opt=optimizationType,
+                         mergeLinkBuses=mergeLinkBuses, mergeBuses=mergeBuses, dispatchMode=dispatchMode)
+
+    # optimize the energy network
+    limit, capacitiesTransformers, capacitiesStorages = network.optimize(solver='gurobi',
+                                                                         numberOfBuildings=numberOfBuildings,
+                                                                         options=optimizationOptions,
+                                                                         mergeLinkBuses=mergeLinkBuses)
+
+    # print optimization outputs i.e. costs, environmental impact and capacities selected for different components (
+    # with investment optimization)
+    network.printInvestedCapacities(capacitiesTransformers, capacitiesStorages)
+    network.printCosts()
+    network.printEnvImpacts()
+    network.printMetaresults()
+
+    # save results
+    if not os.path.exists(resultFilePath):
+        os.makedirs(resultFilePath)
+    network.exportToExcel(os.path.join(resultFilePath, resultFileName), mergeLinkBuses=mergeLinkBuses)
+    network.exportIceStorageModelParams(os.path.join(resultFilePath, iceStoreFileName))

optihood/energy_network.py

@@ -628,6 +628,33 @@ def printbuildingModelTemperatures(self, filename):
             df[f"epsilonIndoor_B{bNo}"] = epsilonIndoor
         df.to_csv(filename, sep=';', index=False)
 
+    def exportIceStorageModelParams(self, filename):
+        df = pd.DataFrame()
+        df["timestamp"] = self.timeindex[0:-1]
+        for i in range(self.__noOfBuildings):
+            bNo = i + 1
+            tStor = [v for k, v in self._optimizationModel.IceStorageBlock.tStor.get_values().items() if
+                       k[0].label.endswith(f"Building{bNo}")]
+            mIceStor = [v for k, v in self._optimizationModel.IceStorageBlock.mIceStor.get_values().items() if
+                             k[0].label.endswith(f"Building{bNo}")]
+            fIce = [v for k, v in self._optimizationModel.IceStorageBlock.fIce.get_values().items() if
+                     k[0].label.endswith(f"Building{bNo}")]
+            iceStatus = [v for k, v in self._optimizationModel.IceStorageBlock.iceStatus.get_values().items() if
+                             k[0].label.endswith(f"Building{bNo}")]
+            tStor_prev = [v for k, v in self._optimizationModel.IceStorageBlock.tStor_prev.get_values().items() if
+                            k[0].label.endswith(f"Building{bNo}")]
+            mIceStor_prev = [v for k, v in
+                                  self._optimizationModel.IceStorageBlock.mIceStor_prev.get_values().items() if
+                                  k[0].label.endswith(f"Building{bNo}")]
+            if tStor:
+                df[f"tStor_B{bNo}"] = tStor
+                df[f"mIceStor_B{bNo}"] = mIceStor
+                df[f"fIce_B{bNo}"] = fIce
+                df[f"iceStatus_prev_B{bNo}"] = iceStatus
+                df[f"tStor_prev_B{bNo}"] = tStor_prev
+                df[f"mIceStor_prev_B{bNo}"] = mIceStor_prev
+        df.to_csv(filename, sep=';', index=False)
+
     def saveUnprocessedResults(self, resultFile):
         with pd.ExcelWriter(resultFile) as writer:
             busLabelList = []

optihood/storages.py

@@ -377,8 +377,8 @@ def _create(self, group=None):
 
         #  ************* DECISION VARIABLES *****************************
         # temperature of ice storage
-        self.tStor = Var(self.icestorages, m.TIMESTEPS, within=NonNegativeReals, bounds=(0,200))
-        self.tStor_prev = Var(self.icestorages, m.TIMESTEPS, within=NonNegativeReals, bounds=(0,200))
+        self.tStor = Var(self.icestorages, m.TIMESTEPS, within=NonNegativeReals, bounds=(0,65))
+        self.tStor_prev = Var(self.icestorages, m.TIMESTEPS, within=NonNegativeReals, bounds=(0,65))
         # mass of ice in storage
         self.mIceStor = Var(self.icestorages, m.TIMESTEPS, within=NonNegativeReals, bounds=(0,100000))
         self.mIceStor_prev = Var(self.icestorages, m.TIMESTEPS, within=NonNegativeReals, bounds=(0,100000))
@@ -389,6 +389,7 @@ def _create(self, group=None):
 
         # for linearization of non-linear constraints with big M method
         M = 2000
+        epsilon = 0.000001
 
         #  ************* CONSTRAINTS *****************************
 
@@ -425,6 +426,17 @@ def _prev_temperature_rule(block):
         self.prev_temperature = Constraint(group, m.TIMESTEPS, noruleinit=True)
         self.prev_temperature_build = BuildAction(rule=_prev_temperature_rule)
 
+        def _prev_mass_ice_rule(block):
+            for g in group:
+                for t in m.TIMESTEPS:
+                    if t != 0:
+                        lhs = self.mIceStor_prev[g, t]
+                        rhs = self.mIceStor[g, t - 1]
+                        block.prev_mass_ice.add((g, t), (lhs == rhs))
+
+        self.prev_mass_ice = Constraint(group, m.TIMESTEPS, noruleinit=True)
+        self.prev_mass_ice_build = BuildAction(rule=_prev_mass_ice_rule)
+
         def _max_ice_fraction_rule(block):
             """set the value of ice fraction"""
             for g in group:
@@ -484,8 +496,8 @@ def _ice_state_rule_1(block):
             """rule for calculating the mass of ice in each timestep"""
             for g in group:
                 for t in m.TIMESTEPS:
-                    lhs = self.iceStatus[g,t]
-                    rhs = self.tStor[g,t]
+                    lhs = self.tStor[g,t]
+                    rhs = M*(1-self.iceStatus[g,t])
                     block.ice_state_1.add((g, t), (lhs <= rhs))
 
         self.ice_state_1 = Constraint(group, m.TIMESTEPS, noruleinit=True)
@@ -495,8 +507,8 @@ def _ice_state_rule_2(block):
             """rule for calculating the mass of ice in each timestep"""
             for g in group:
                 for t in m.TIMESTEPS:
-                    lhs = self.iceStatus[g,t]
-                    rhs = self.tStor[g,t] - M*(1-self.tStor[g,t])
+                    lhs = self.tStor[g,t]
+                    rhs = epsilon*(1-self.iceStatus[g,t])
                     block.ice_state_2.add((g, t), (lhs >= rhs))
 
         self.ice_state_2 = Constraint(group, m.TIMESTEPS, noruleinit=True)

tests/test_examples/test_config_examples.py

@@ -10,7 +10,7 @@
 cwd = _os.getcwd()
 packageDir = _pl.Path(_oh.__file__).resolve().parent
 scriptDir = packageDir / ".." / "data" / "examples"
-example_path = packageDir / ".." / "data" / "results" / "HP_GS_PV" / "cost" / "group"
+example_path = packageDir / ".." / "data" / "results"
 expected_data_dir = _pl.Path(__file__).resolve().parent / "expected_files"
 
 _SHEET_NAMES = ['gridBus__Building1', 'electricityProdBus__Building1', 'electricityInBus', 'shDemandBus',
@@ -46,7 +46,7 @@ def test_group_optimization_before_merge(self):
         _os.chdir(cwd)
         # =============================
 
-        excel_file_path = str(example_path / "scenario_HP_GS_PV.xls")
+        excel_file_path = str(example_path / "results_example_config.xlsx")
         expected_data_path = str(expected_data_dir / "test_run_example_config.xls")
 
         compare_xls_files(self, excel_file_path, expected_data_path, _SHEET_NAMES, manual_test=True)
@@ -66,7 +66,7 @@ def test_group_optimization_after_merge(self):
         _os.chdir(cwd)
         # =============================
 
-        excel_file_path = example_path / "scenario_HP_GS_PV.xls"
+        excel_file_path = example_path / "results_example_config.xlsx"
         assert excel_file_path.is_file()
         # expected_data_path = str(expected_data_dir / "test_run_example_config_after_merge.xls")