(WIP) feature: add pipe sizing to dhn analysis

Reviewed-on: https://nextgenerations-cities.encs.concordia.ca/gitea/s_ranjbar/energy_system_modelling_workflow/pulls/17

energy_system_retrofit.py

@@ -92,9 +92,8 @@ for building in city.buildings:
                        fuel_tariffs=['Electricity-D', 'Gas-Energir']).life_cycle
   lcc_dataframe.to_csv(cost_analysis_output_path / f'{building.name}_retrofitted_lcc.csv')
   retrofitted_life_cycle_cost[f'{building.name}'] = cost_data(building, lcc_dataframe, cost_retrofit_scenario)
-for i in range(12):
-  dhw_consumption = 0
-  for building in city.buildings:
-    dhw_consumption += building.domestic_hot_water_consumption[cte.MONTH][i] / 3.6e6
+EnergySystemRetrofitReport(city, output_path, 'PV Implementation and System Retrofit',
+                           current_status_energy_consumption, retrofitted_energy_consumption,
+                           current_status_life_cycle_cost, retrofitted_life_cycle_cost).create_report()
 
 

input_files/test_geojson.geojson

@@ -1,55 +0,0 @@
-{
-  "type": "FeatureCollection",
-  "features": [
-    {
-      "type": "Feature",
-      "geometry": {
-        "type": "Polygon",
-        "coordinates": [
-          [
-            [
-              -73.58000127109773,
-              45.49613461675315
-            ],
-            [
-              -73.57962787855432,
-              45.496524875557746
-            ],
-            [
-              -73.57996357265695,
-              45.49668114195629
-            ],
-            [
-              -73.57996427397713,
-              45.496680342403664
-            ],
-            [
-              -73.58034707390021,
-              45.49625804233725
-            ],
-            [
-              -73.58034697395713,
-              45.496257942524835
-            ],
-            [
-              -73.58000127109773,
-              45.49613461675315
-            ]
-          ]
-        ]
-      },
-      "id": 179764,
-      "properties": {
-        "name": "01119274",
-        "address": "rue Guy (MTL) 2157",
-        "function": "Mixed use",
-        "mixed_type_1": "commercial",
-        "mixed_type_1_percentage": 50,
-        "mixed_type_2": "6000",
-        "mixed_type_2_percentage": 50,
-        "height": 62,
-        "year_of_construction": 1954
-      }
-    }
-  ]
-}

input_files/test_geojson1.geojson

@@ -1,52 +0,0 @@
-{
-  "type": "FeatureCollection",
-  "features": [
-    {
-      "type": "Feature",
-      "geometry": {
-        "type": "Polygon",
-        "coordinates": [
-          [
-            [
-              -73.58000127109773,
-              45.49613461675315
-            ],
-            [
-              -73.57962787855432,
-              45.496524875557746
-            ],
-            [
-              -73.57996357265695,
-              45.49668114195629
-            ],
-            [
-              -73.57996427397713,
-              45.496680342403664
-            ],
-            [
-              -73.58034707390021,
-              45.49625804233725
-            ],
-            [
-              -73.58034697395713,
-              45.496257942524835
-            ],
-            [
-              -73.58000127109773,
-              45.49613461675315
-            ]
-          ]
-        ]
-      },
-      "id": 179764,
-      "properties": {
-        "name": "01119274",
-        "address": "rue Guy (MTL) 2157",
-        "function": "Mixed use",
-        "usages": [{"usage": "commercial", "percentage": 50}, {"usage": "6000", "percentage": 50}],
-        "height": 62,
-        "year_of_construction": 1954
-      }
-    }
-  ]
-}

main.py

@@ -1,4 +1,5 @@
 from pathlib import Path
+from scripts.district_heating_network.directory_manager import DirectoryManager
 import subprocess
 from scripts.ep_run_enrich import energy_plus_workflow
 from hub.imports.geometry_factory import GeometryFactory
@@ -22,31 +23,39 @@ import hub.helpers.constants as cte
 from hub.exports.exports_factory import ExportsFactory
 from scripts.pv_feasibility import pv_feasibility
 import matplotlib.pyplot as plt
-import numpy as np
-# Specify the GeoJSON file path
-data = {}
-input_files_path = (Path(__file__).parent / 'input_files')
-input_files_path.mkdir(parents=True, exist_ok=True)
-# geojson_file = process_geojson(x=-73.58001358793511, y=45.496445294438715, diff=0.0001)
-geojson_file_path = input_files_path / 'test_geojson1.geojson'
-output_path = (Path(__file__).parent / 'out_files').resolve()
-output_path.mkdir(parents=True, exist_ok=True)
-energy_plus_output_path = output_path / 'energy_plus_outputs'
-energy_plus_output_path.mkdir(parents=True, exist_ok=True)
-simulation_results_path = (Path(__file__).parent / 'out_files' / 'simulation_results').resolve()
-simulation_results_path.mkdir(parents=True, exist_ok=True)
-sra_output_path = output_path / 'sra_outputs'
-sra_output_path.mkdir(parents=True, exist_ok=True)
-cost_analysis_output_path = output_path / 'cost_analysis'
-cost_analysis_output_path.mkdir(parents=True, exist_ok=True)
+from scripts.district_heating_network.district_heating_network_creator import DistrictHeatingNetworkCreator
+from scripts.district_heating_network.road_processor import road_processor
+from scripts.district_heating_network.district_heating_factory import DistrictHeatingFactory
+
+base_path = Path(__file__).parent
+dir_manager = DirectoryManager(base_path)
+
+# Input files directory
+input_files_path = dir_manager.create_directory('input_files')
+geojson_file_path = input_files_path / 'output_buildings.geojson'
+
+# Output files directory
+output_path = dir_manager.create_directory('out_files')
+
+# Subdirectories for output files
+energy_plus_output_path = dir_manager.create_directory('out_files/energy_plus_outputs')
+simulation_results_path = dir_manager.create_directory('out_files/simulation_results')
+sra_output_path = dir_manager.create_directory('out_files/sra_outputs')
+cost_analysis_output_path = dir_manager.create_directory('out_files/cost_analysis')
+
+# Select city area
+location = [45.53067276979674, -73.70234652694087]
+process_geojson(x=location[1], y=location[0], diff=0.001)
+
+# Create city object
 city = GeometryFactory(file_type='geojson',
                        path=geojson_file_path,
                        height_field='height',
                        year_of_construction_field='year_of_construction',
                        function_field='function',
                        function_to_hub=Dictionaries().montreal_function_to_hub_function).city
-# ConstructionFactory('nrcan', city).enrich()
-# UsageFactory('nrcan', city).enrich()
+ConstructionFactory('nrcan', city).enrich()
+UsageFactory('nrcan', city).enrich()
 # WeatherFactory('epw', city).enrich()
 # energy_plus_workflow(city, energy_plus_output_path)
 # data[f'{city.buildings[0].function}'] = city.buildings[0].heating_demand[cte.YEAR][0] / 3.6e9
@@ -83,4 +92,4 @@ city = GeometryFactory(file_type='geojson',
 # # Save the plot
 # plt.savefig('plot_nrcan.png', bbox_inches='tight')
 # plt.close()
-print('test')
+print('test')

network_edges.csv

@@ -0,0 +1,103 @@
+Start Node,End Node,Flow Rate,Diameter
+1,2,7.036628334723999,100.31636393633916
+1,3,4.21120645112289,77.60554307494739
+1,55,3.5317773545013225,71.06997042924019
+2,4,7.480662304076292,103.43309043327245
+2,69,0.5550424616902951,28.174264035849028
+3,29,3.996190684370673,75.5983953887572
+3,88,0.2687697084402423,19.605577350674327
+4,5,12.739252570195731,134.97752174009295
+4,70,6.57323783264924,96.95699548330485
+5,6,1.440813373451662,45.3934656880935
+5,7,14.180065943647435,142.40610275642487
+6,25,0.6885299581772856,31.379852409201234
+6,66,0.9403542690929205,36.67205322096126
+7,8,18.919213781899757,164.4905748238601
+7,94,5.923934797815284,92.04381499868917
+8,9,20.2831963688024,170.31687247962364
+8,62,1.704978233628294,49.379750978468536
+9,10,23.795992013358696,184.47664139533458
+9,74,4.390994555695343,79.24482766451129
+10,11,24.08473660300742,185.59250196750875
+10,85,0.36093073706075773,22.71963814788559
+11,12,29.818111567083733,206.5045112717333
+11,97,7.166718705095327,101.23942160182479
+12,13,30.27814518525654,208.09138919592667
+12,87,0.5750420227159534,28.677366552586644
+13,14,42.1891802453921,245.63487040066076
+13,60,14.888793825169447,145.92148273668778
+14,15,42.46951106929258,246.44959372826426
+14,84,0.3504135298756097,22.386175758632653
+15,16,42.71334684180593,247.15606782002874
+15,100,0.30479471564177585,20.878206682814064
+16,17,43.164835926660224,248.45888000248186
+16,82,0.5643613560678524,28.40979566604211
+17,18,43.40956096248419,249.1622090841966
+17,67,0.3059062947800095,20.91624319844782
+18,19,44.36915698128853,251.90110032375995
+18,71,1.1994950235054551,41.417959700534844
+19,20,44.628206053508535,252.63539166978182
+19,90,0.3238113402749617,21.5196648389043
+20,21,46.27014639132988,257.24083624469847
+20,75,2.052425422276719,54.178024923945806
+21,22,73.0360691675874,323.1901912619741
+21,23,26.455032668422902,194.5107583389826
+21,24,0.31089010783469323,21.08593812173577
+22,49,84.50840788728556,347.6477603374658
+22,81,14.340423399622829,143.2090497302142
+23,31,25.527872344225216,191.07188392806023
+23,68,1.158950405247153,40.71194974740254
+24,56,0.3886126347933523,23.574795504777274
+25,26,0.48524655942646255,26.34333218820642
+25,95,0.2541042484385019,19.063183969533693
+26,27,0.2982895780206631,20.654206550796754
+26,65,0.23369622675732482,18.28164730458197
+27,28,0.17282514505654928,15.721462386442168
+27,59,0.15683054120527573,14.976309152909318
+28,57,0.2160314313208167,17.577129300900673
+29,30,3.8183263880744693,73.89686263587296
+29,58,0.22233037037019066,17.831540803086607
+30,33,3.699868945604976,72.74156600819312
+30,77,0.1480718030867907,14.552099595104522
+31,32,20.88880856250903,172.8408171354717
+31,61,5.79882972714514,91.06671153446011
+32,46,19.8454823967694,168.46911947690657
+32,83,1.304157707174466,43.18714870672264
+33,34,3.4711065659703353,70.45688678723134
+33,102,0.2859529745433158,20.222590387225335
+34,35,3.4711065659703735,70.45688678723174
+35,36,2.259511486129423,56.84558663220748
+35,86,1.5144938498012546,46.539662871966804
+36,37,1.7967791166489435,50.691695985759985
+36,76,0.5784154618507291,28.76136031308409
+37,38,1.3079368277973644,43.24967618306041
+37,98,0.6110528610646417,29.56166334911917
+38,39,1.161455220287694,40.75592094041735
+38,93,0.18310200938718202,16.182142867916955
+39,40,1.0509571713728947,38.76876664450215
+39,101,0.1381225611435788,14.054706710110416
+40,41,0.9248009014795987,36.367512616181635
+40,73,0.15769533736659458,15.017543626027152
+41,42,0.7870341012651642,33.54951537030909
+41,79,0.1722085002680287,15.693390021976926
+42,43,0.6096172477277307,29.526916703688457
+42,72,0.22177106692165793,17.809097799856982
+43,44,0.5017777614037179,26.788301656702675
+43,63,0.13479935790496445,13.884600479307904
+44,45,0.257385665027338,19.185876901497235
+44,80,0.3054901204704009,20.902010464636216
+45,64,0.32173208128418995,21.45046247984601
+46,47,18.899025568880468,164.4027895397298
+46,99,1.1830710348610967,41.133426323360354
+47,48,17.204865763811963,156.86105088056533
+47,78,2.1176997563356554,55.032807545697466
+48,50,15.889475871949896,150.74546712455864
+48,89,1.644237364827739,48.492182678681466
+49,53,4.137682980635671,76.92510229741136
+49,54,88.64609086792123,356.0567884918405
+50,52,15.88947587194994,150.74546712455884
+51,52,15.107521674531634,146.98942430449688
+51,96,18.884402093164404,164.33917235920677
+52,91,0.977442746772752,37.38825018026857
+53,92,5.172103725794528,86.004878956568
+54,103,88.64609086792132,356.0567884918407

pv_assessment.py

@@ -14,7 +14,7 @@ import hub.helpers.constants as cte
 from hub.exports.exports_factory import ExportsFactory
 from scripts.pv_sizing_and_simulation import PVSizingSimulation
 # Specify the GeoJSON file path
-geojson_file = process_geojson(x=-73.5681295982132, y=45.49218262677643, diff=0.0005)
+geojson_file = process_geojson(x=-73.5681295982132, y=45.49218262677643, diff=0.0001)
 file_path = (Path(__file__).parent / 'input_files' / 'output_buildings.geojson')
 # Specify the output path for the PDF file
 output_path = (Path(__file__).parent / 'out_files').resolve()
@@ -34,7 +34,7 @@ ExportsFactory('sra', city, output_path).export()
 sra_path = (output_path / f'{city.name}_sra.xml').resolve()
 subprocess.run(['sra', str(sra_path)])
 ResultFactory('sra', city, output_path).enrich()
-energy_plus_workflow(city)
+energy_plus_workflow(city, output_path=output_path)
 solar_angles = CitySolarAngles(city.name,
                                city.latitude,
                                city.longitude,

scripts/district_heating_network/directory_manager.py

@@ -0,0 +1,16 @@
+from pathlib import Path
+
+
+class DirectoryManager:
+    def __init__(self, base_path):
+        self.base_path = Path(base_path)
+        self.directories = {}
+
+    def create_directory(self, relative_path):
+        full_path = self.base_path / relative_path
+        full_path.mkdir(parents=True, exist_ok=True)
+        self.directories[relative_path] = full_path
+        return full_path
+
+    def get_directory(self, relative_path):
+        return self.directories.get(relative_path, None)

scripts/district_heating_network/district_heating_exporter.py

@@ -0,0 +1,86 @@
+import json
+import csv
+import logging
+
+
+class NetworkGraphExporter:
+  """
+  A class to export a network graph to various formats like CSV and GeoJSON.
+  """
+
+  def __init__(self, graph):
+    """
+    Initialize the exporter with a network graph.
+
+    :param graph: A NetworkX graph object.
+    """
+    self.graph = graph
+
+  def to_csv(self, file_path):
+    """
+    Save the graph nodes with their type, names, and positions to a CSV file.
+
+    :param file_path: The path to the output CSV file.
+    """
+    try:
+      with open(file_path, mode='w', newline='') as file:
+        writer = csv.writer(file)
+        writer.writerow(['Node ID', 'Name', 'Type', 'Position'])
+        for node, data in self.graph.nodes(data=True):
+          writer.writerow([node, data['name'], data['type'], data['pos']])
+      logging.info(f"Graph successfully saved to CSV file: {file_path}")
+    except Exception as e:
+      logging.error(f"Error saving graph to CSV file: {e}")
+
+  def to_geojson(self, file_path):
+    """
+    Save the graph to a GeoJSON file.
+
+    :param file_path: The path to the output GeoJSON file.
+    """
+    try:
+      features = []
+
+      for node, data in self.graph.nodes(data=True):
+        feature = {
+          "type": "Feature",
+          "geometry": {
+            "type": "Point",
+            "coordinates": data['pos']
+          },
+          "properties": {
+            "id": node,
+            "name": data['name'],
+            "type": data['type']
+          }
+        }
+        features.append(feature)
+
+      for u, v, data in self.graph.edges(data=True):
+        feature = {
+          "type": "Feature",
+          "geometry": {
+            "type": "LineString",
+            "coordinates": [self.graph.nodes[u]['pos'], self.graph.nodes[v]['pos']]
+          },
+          "properties": {
+            "length": data['length']
+          }
+        }
+        features.append(feature)
+
+      geojson = {
+        "type": "FeatureCollection",
+        "features": features
+      }
+
+      with open(file_path, 'w') as f:
+        json.dump(geojson, f, indent=2)
+      logging.info(f"Graph successfully saved to GeoJSON file: {file_path}")
+    except Exception as e:
+      logging.error(f"Error saving graph to GeoJSON file: {e}")
+
+
+# Configure logging
+logging.basicConfig(level=logging.INFO, format="%(asctime)s - %(levelname)s - %(message)s")
+logging.getLogger("numexpr").setLevel(logging.ERROR)

scripts/district_heating_network/district_heating_factory.py

@@ -0,0 +1,174 @@
+import networkx as nx
+import logging
+import CoolProp.CoolProp as CP  # Ensure correct import
+import math
+import numpy as np
+import matplotlib.pyplot as plt
+import csv
+
+class DistrictHeatingFactory:
+  """
+  DistrictHeatingFactory class
+  """
+
+  def __init__(self, city, graph, supply_temperature, return_temperature, simultaneity_factor):
+    self._city = city
+    self._network_graph = graph
+    self._supply_temperature = supply_temperature
+    self._return_temperature = return_temperature
+    self.simultaneity_factor = simultaneity_factor
+    self.fluid = "Water"
+
+  def enrich(self):
+    """
+    Enrich the network graph nodes with the whole building object from the city buildings.
+    """
+    for node in self._network_graph.nodes(data=True):
+      node_id, node_attrs = node
+      if node_attrs.get('type') == 'building':
+        building_name = node_attrs.get('name')
+        building_found = False
+        for building in self._city.buildings:
+          if building.name == building_name:
+            self._network_graph.nodes[node_id]['building_obj'] = building
+            building_found = True
+            break
+        if not building_found:
+          logging.error(msg=f"Building with name '{building_name}' not found in city.")
+
+  def calculate_flow_rates(self, A, Gext):
+    """
+    Solve the linear system to find the flow rates in each branch.
+    """
+    try:
+      G = np.linalg.lstsq(A, Gext, rcond=None)[0]
+      return G
+    except np.linalg.LinAlgError as e:
+      logging.error(f"Error solving the linear system: {e}")
+      return None
+
+  def switch_nodes(self, A, edge_index, node_index, edge):
+    """
+    Switch the in and out nodes for the given edge.
+    """
+    u, v = edge
+    i = node_index[u]
+    j = node_index[v]
+    k = edge_index[edge]
+    A[i, k], A[j, k] = -A[i, k], -A[j, k]
+
+  def sizing(self):
+    """
+    Calculate the diameter of the pipes in the district heating network.
+    """
+    # Calculate the peak flow rates for buildings
+    total_peak_load = 0
+    for node in self._network_graph.nodes(data=True):
+      node_id, node_attrs = node
+      if node_attrs.get('type') == 'building':
+        building = node_attrs.get('building_obj')
+        if building and "year" in building.heating_peak_load:
+          peak_load = building.heating_peak_load["year"][0]
+          if peak_load:
+            try:
+              specific_heat_capacity = CP.PropsSI('C', 'T', (self._supply_temperature + self._return_temperature) / 2, 'P', 101325, self.fluid)
+              peak_mass_flow_rate = peak_load / (specific_heat_capacity * (self._supply_temperature - self._return_temperature))
+              self._network_graph.nodes[node_id]['peak_flow'] = peak_mass_flow_rate
+              total_peak_load += peak_mass_flow_rate
+            except Exception as e:
+              logging.error(f"Error calculating peak mass flow rate for building {building.name}: {e}")
+
+    # Assign the total peak load to the central plant node
+    for node in self._network_graph.nodes(data=True):
+      node_id, node_attrs = node
+      if node_attrs.get('type') == 'generation':
+        self._network_graph.nodes[node_id]['peak_flow'] = -total_peak_load
+
+    # Create the incidence matrix A and Gext vector
+    num_nodes = self._network_graph.number_of_nodes()
+    num_edges = self._network_graph.number_of_edges()
+
+    A = np.zeros((num_nodes, num_edges))
+    Gext = np.zeros(num_nodes)
+
+    node_index = {node: i for i, node in enumerate(self._network_graph.nodes())}
+    edge_index = {edge: i for i, edge in enumerate(self._network_graph.edges())}
+
+    for node, data in self._network_graph.nodes(data=True):
+      i = node_index[node]
+      Gext[i] = data.get('peak_flow', 0)
+
+    for edge in self._network_graph.edges():
+      u, v = edge
+      i = node_index[u]
+      j = node_index[v]
+      k = edge_index[edge]
+      A[i, k] = 1
+      A[j, k] = -1
+
+    # Calculate initial flow rates
+    G = self.calculate_flow_rates(A, Gext)
+    if G is None:
+      return
+
+    # Check for negative flow rates and switch nodes if necessary
+    iterations = 0
+    max_iterations = num_edges * 2  # Arbitrary limit to avoid infinite loop
+
+    while any(flow_rate < 0 for flow_rate in G) and iterations < max_iterations:
+      for idx, flow_rate in enumerate(G):
+        if flow_rate < 0:
+          edge = list(self._network_graph.edges())[idx]
+          self.switch_nodes(A, edge_index, node_index, edge)
+      G = self.calculate_flow_rates(A, Gext)
+      if G is None:
+        return
+      iterations += 1
+
+    # Update the edge diameters based on the calculated flow rates
+    for edge, flow_rate in zip(self._network_graph.edges(), G):
+      u, v = edge
+      if not (self._network_graph.nodes[u].get('type') == 'building' or self._network_graph.nodes[v].get('type') == 'building'):
+        flow_rate *= self.simultaneity_factor
+
+      try:
+        density = CP.PropsSI('D', 'T', (self._supply_temperature + self._return_temperature) / 2, 'P', 101325, self.fluid)
+        velocity = 0.9  # m/s
+        diameter = math.sqrt((4 * abs(flow_rate)) / (density * velocity * math.pi)) * 1000  # mm
+        self._network_graph.edges[u, v]['diameter'] = diameter
+        self._network_graph.edges[u, v]['flow_rate'] = flow_rate
+      except Exception as e:
+        logging.error(f"Error calculating diameter for edge ({u}, {v}): {e}")
+
+  def plot_network(self):
+    """
+    Plot the network graph with flow rates shown as colors.
+    """
+    pos = {node: (data['pos'][0], data['pos'][1]) for node, data in self._network_graph.nodes(data=True) if 'pos' in data}
+    edge_colors = [self._network_graph[u][v]['flow_rate'] for u, v in self._network_graph.edges()]
+
+    plt.figure(figsize=(18, 12))
+    nodes = nx.draw_networkx_nodes(self._network_graph, pos, node_color='skyblue', node_size=700)
+    edges = nx.draw_networkx_edges(self._network_graph, pos, edge_color=edge_colors, edge_cmap=plt.cm.viridis, width=2)
+    labels = nx.draw_networkx_labels(self._network_graph, pos, font_size=10, font_color='black')
+
+    sm = plt.cm.ScalarMappable(cmap=plt.cm.viridis, norm=plt.Normalize(vmin=min(edge_colors), vmax=max(edge_colors)))
+    sm.set_array([])
+    cbar = plt.colorbar(sm, label='Flow rate')
+    cbar.ax.tick_params(labelsize=10)
+
+    plt.title('District Heating Network with Flow Rates', fontsize=16)
+    plt.axis('off')
+    plt.legend(handles=[nodes], labels=['Nodes'], loc='upper left')
+    plt.show()
+
+  def save_edges_to_csv(self, file_path):
+    """
+    Save the edges and their attributes to a CSV file.
+    """
+    with open(file_path, mode='w', newline='') as file:
+      writer = csv.writer(file)
+      writer.writerow(['Start Node', 'End Node', 'Flow Rate', 'Diameter'])
+      for u, v, data in self._network_graph.edges(data=True):
+        writer.writerow([u, v, data.get('flow_rate', ''), data.get('diameter', '')])
+    logging.info(f"Edges successfully saved to CSV file: {file_path}")

scripts/district_heating_network/district_heating_network_creator.py

@@ -0,0 +1,372 @@
+import json
+import math
+import logging
+import matplotlib.pyplot as plt
+import networkx as nx
+from shapely.geometry import Polygon, Point, LineString
+from typing import List, Tuple
+from rtree import index
+
+# Configure logging
+logging.basicConfig(level=logging.INFO, format="%(asctime)s - %(levelname)s - %(message)s")
+logging.getLogger("numexpr").setLevel(logging.ERROR)
+
+def haversine(lon1, lat1, lon2, lat2):
+  """
+  Calculate the great-circle distance between two points
+  on the Earth specified by their longitude and latitude.
+  """
+  R = 6371000  # Radius of the Earth in meters
+  phi1 = math.radians(lat1)
+  phi2 = math.radians(lat2)
+  delta_phi = math.radians(lat2 - lat1)
+  delta_lambda = math.radians(lon2 - lon1)
+
+  a = math.sin(delta_phi / 2.0) ** 2 + \
+      math.cos(phi1) * math.cos(phi2) * \
+      math.sin(delta_lambda / 2.0) ** 2
+
+  c = 2 * math.atan2(math.sqrt(a), math.sqrt(1 - a))
+  return R * c  # Output distance in meters
+
+class DistrictHeatingNetworkCreator:
+  def __init__(self, buildings_file: str, roads_file: str, central_plant_locations: List[Tuple[float, float]]):
+    """
+    Initialize the class with paths to the buildings and roads data files, and central plant locations.
+
+    :param buildings_file: Path to the GeoJSON file containing building data.
+    :param roads_file: Path to the GeoJSON file containing roads data.
+    :param central_plant_locations: List of tuples containing the coordinates of central plant locations.
+    """
+    if len(central_plant_locations) < 1:
+      raise ValueError("The list of central plant locations must have at least one member.")
+
+    self.buildings_file = buildings_file
+    self.roads_file = roads_file
+    self.central_plant_locations = central_plant_locations
+
+  def run(self) -> nx.Graph:
+    """
+    Main method to execute the district heating network creation process.
+    :return: NetworkX graph with nodes and edges representing the network.
+    """
+    try:
+      self._load_and_process_data()
+      self._find_nearest_roads()
+      self._find_nearest_points()
+      self._break_down_roads()
+      self._create_graph()
+      self._create_mst()
+      self._iteratively_remove_edges()
+      self._add_centroids_to_mst()
+      self._convert_edge_weights_to_meters()
+      self._create_final_network_graph()
+      return self.network_graph
+    except Exception as e:
+      logging.error(f"Error during network creation: {e}")
+      raise
+
+  def _load_and_process_data(self):
+    """
+    Load and process the building and road data.
+    """
+    try:
+      # Load building data
+      with open(self.buildings_file, 'r') as file:
+        city = json.load(file)
+
+      self.centroids = []
+      self.building_names = []
+      self.building_positions = []
+      buildings = city['features']
+      for building in buildings:
+        coordinates = building['geometry']['coordinates'][0]
+        building_polygon = Polygon(coordinates)
+        centroid = building_polygon.centroid
+        self.centroids.append(centroid)
+        self.building_names.append(str(building['id']))
+        self.building_positions.append((centroid.x, centroid.y))
+
+      # Add central plant locations as centroids
+      for i, loc in enumerate(self.central_plant_locations, start=1):
+        centroid = Point(loc)
+        self.centroids.append(centroid)
+        self.building_names.append(f'central_plant_{i}')
+        self.building_positions.append((centroid.x, centroid.y))
+
+      # Load road data
+      with open(self.roads_file, 'r') as file:
+        roads = json.load(file)
+
+      line_features = [feature for feature in roads['features'] if feature['geometry']['type'] == 'LineString']
+
+      self.lines = [LineString(feature['geometry']['coordinates']) for feature in line_features]
+      self.cleaned_lines = [LineString([line.coords[0], line.coords[-1]]) for line in self.lines]
+    except Exception as e:
+      logging.error(f"Error loading and processing data: {e}")
+      raise
+
+  def _find_nearest_roads(self):
+    """
+    Find the nearest road for each building centroid.
+    """
+    try:
+      self.closest_roads = []
+      unique_roads_set = set()
+
+      # Create spatial index for roads
+      idx = index.Index()
+      for pos, line in enumerate(self.cleaned_lines):
+        idx.insert(pos, line.bounds)
+
+      for centroid in self.centroids:
+        min_distance = float('inf')
+        closest_road = None
+        for pos in idx.nearest(centroid.bounds, 10):
+          road = self.cleaned_lines[pos]
+          distance = road.distance(centroid)
+          if distance < min_distance:
+            min_distance = distance
+            closest_road = road
+
+        if closest_road and closest_road.wkt not in unique_roads_set:
+          unique_roads_set.add(closest_road.wkt)
+          self.closest_roads.append(closest_road)
+    except Exception as e:
+      logging.error(f"Error finding nearest roads: {e}")
+      raise
+
+  def _find_nearest_points(self):
+    """
+    Find the nearest point on each closest road for each centroid.
+    """
+
+    def find_nearest_point_on_line(point: Point, line: LineString) -> Point:
+      return line.interpolate(line.project(point))
+
+    try:
+      self.nearest_points = []
+      for centroid in self.centroids:
+        min_distance = float('inf')
+        closest_road = None
+        for road in self.closest_roads:
+          distance = centroid.distance(road)
+          if distance < min_distance:
+            min_distance = distance
+            closest_road = road
+
+        if closest_road:
+          nearest_point = find_nearest_point_on_line(centroid, closest_road)
+          self.nearest_points.append(nearest_point)
+    except Exception as e:
+      logging.error(f"Error finding nearest points: {e}")
+      raise
+
+  def _break_down_roads(self):
+    """
+    Break down roads into segments connecting nearest points.
+    """
+
+    def break_down_roads(closest_roads: List[LineString], nearest_points_list: List[Point]) -> List[LineString]:
+      new_segments = []
+      for road in closest_roads:
+        coords = list(road.coords)
+        points_on_road = [point for point in nearest_points_list if road.distance(point) < 0.000000001]
+        sorted_points = sorted(points_on_road, key=lambda point: road.project(point))
+        sorted_points.insert(0, Point(coords[0]))
+        sorted_points.append(Point(coords[-1]))
+        for i in range(len(sorted_points) - 1):
+          segment = LineString([sorted_points[i], sorted_points[i + 1]])
+          new_segments.append(segment)
+      return new_segments
+
+    try:
+      self.new_segments = break_down_roads(self.closest_roads, self.nearest_points)
+      self.cleaned_lines = [line for line in self.cleaned_lines if line not in self.closest_roads]
+      self.cleaned_lines.extend(self.new_segments)
+    except Exception as e:
+      logging.error(f"Error breaking down roads: {e}")
+      raise
+
+  def _create_graph(self):
+    """
+    Create a NetworkX graph from the cleaned lines.
+    """
+    try:
+      self.G = nx.Graph()
+      for line in self.cleaned_lines:
+        coords = list(line.coords)
+        for i in range(len(coords) - 1):
+          u = coords[i]
+          v = coords[i + 1]
+          self.G.add_edge(u, v, weight=Point(coords[i]).distance(Point(coords[i + 1])))
+    except Exception as e:
+      logging.error(f"Error creating graph: {e}")
+      raise
+
+  def _create_mst(self):
+    """
+    Create a Minimum Spanning Tree (MST) from the graph.
+    """
+
+    def find_paths_between_nearest_points(g: nx.Graph, nearest_points: List[Point]) -> List[Tuple]:
+      edges = []
+      for i, start_point in enumerate(nearest_points):
+        start = (start_point.x, start_point.y)
+        for end_point in nearest_points[i + 1:]:
+          end = (end_point.x, end_point.y)
+          if nx.has_path(g, start, end):
+            path = nx.shortest_path(g, source=start, target=end, weight='weight')
+            path_edges = list(zip(path[:-1], path[1:]))
+            edges.extend((u, v, g[u][v]['weight']) for u, v in path_edges)
+      return edges
+
+    try:
+      edges = find_paths_between_nearest_points(self.G, self.nearest_points)
+      h = nx.Graph()
+      h.add_weighted_edges_from(edges)
+      mst = nx.minimum_spanning_tree(h, weight='weight')
+      final_edges = []
+      for u, v in mst.edges():
+        if nx.has_path(self.G, u, v):
+          path = nx.shortest_path(self.G, source=u, target=v, weight='weight')
+          path_edges = list(zip(path[:-1], path[1:]))
+          final_edges.extend((x, y, self.G[x][y]['weight']) for x, y in path_edges)
+      self.final_mst = nx.Graph()
+      self.final_mst.add_weighted_edges_from(final_edges)
+    except Exception as e:
+      logging.error(f"Error creating MST: {e}")
+      raise
+
+  def _iteratively_remove_edges(self):
+    """
+    Iteratively remove edges that do not have any nearest points and have one end with only one connection.
+    Also remove nodes that don't have any connections and street nodes with only one connection.
+    """
+    nearest_points_tuples = [(point.x, point.y) for point in self.nearest_points]
+
+    def find_edges_to_remove(graph: nx.Graph) -> List[Tuple]:
+      edges_to_remove = []
+      for u, v, d in graph.edges(data=True):
+        if u not in nearest_points_tuples and v not in nearest_points_tuples:
+          if graph.degree(u) == 1 or graph.degree(v) == 1:
+            edges_to_remove.append((u, v, d))
+      return edges_to_remove
+
+    def find_nodes_to_remove(graph: nx.Graph) -> List[Tuple]:
+      nodes_to_remove = []
+      for node in graph.nodes():
+        if graph.degree(node) == 0:
+          nodes_to_remove.append(node)
+      return nodes_to_remove
+
+    try:
+      edges_to_remove = find_edges_to_remove(self.final_mst)
+      self.final_mst_steps = [list(self.final_mst.edges(data=True))]
+
+      while edges_to_remove or find_nodes_to_remove(self.final_mst):
+        self.final_mst.remove_edges_from(edges_to_remove)
+        nodes_to_remove = find_nodes_to_remove(self.final_mst)
+        self.final_mst.remove_nodes_from(nodes_to_remove)
+        edges_to_remove = find_edges_to_remove(self.final_mst)
+        self.final_mst_steps.append(list(self.final_mst.edges(data=True)))
+
+      def find_single_connection_street_nodes(graph: nx.Graph) -> List[Tuple]:
+        single_connection_street_nodes = []
+        for node in graph.nodes():
+          if node not in nearest_points_tuples and graph.degree(node) == 1:
+            single_connection_street_nodes.append(node)
+        return single_connection_street_nodes
+
+      single_connection_street_nodes = find_single_connection_street_nodes(self.final_mst)
+
+      while single_connection_street_nodes:
+        for node in single_connection_street_nodes:
+          neighbors = list(self.final_mst.neighbors(node))
+          self.final_mst.remove_node(node)
+          for neighbor in neighbors:
+            if self.final_mst.degree(neighbor) == 0:
+              self.final_mst.remove_node(neighbor)
+        single_connection_street_nodes = find_single_connection_street_nodes(self.final_mst)
+        self.final_mst_steps.append(list(self.final_mst.edges(data=True)))
+    except Exception as e:
+      logging.error(f"Error iteratively removing edges: {e}")
+      raise
+
+  def _add_centroids_to_mst(self):
+    """
+    Add centroids to the final MST graph and connect them to their associated node on the graph.
+    """
+    try:
+      for i, centroid in enumerate(self.centroids):
+        building_name = self.building_names[i]
+        pos = (centroid.x, centroid.y)
+        node_type = 'building' if 'central_plant' not in building_name else 'generation'
+        self.final_mst.add_node(pos, type=node_type, name=building_name, pos=pos)
+
+        nearest_point = None
+        min_distance = float('inf')
+        for node in self.final_mst.nodes():
+          if self.final_mst.nodes[node].get('type') != 'building' and self.final_mst.nodes[node].get('type') != 'generation':
+            distance = centroid.distance(Point(node))
+            if distance < min_distance:
+              min_distance = distance
+              nearest_point = node
+
+        if nearest_point:
+          self.final_mst.add_edge(pos, nearest_point, weight=min_distance)
+    except Exception as e:
+      logging.error(f"Error adding centroids to MST: {e}")
+      raise
+
+  def _convert_edge_weights_to_meters(self):
+    """
+    Convert all edge weights in the final MST graph to meters using the Haversine formula.
+    """
+    try:
+      for u, v, data in self.final_mst.edges(data=True):
+        distance = haversine(u[0], u[1], v[0], v[1])
+        self.final_mst[u][v]['weight'] = distance
+    except Exception as e:
+      logging.error(f"Error converting edge weights to meters: {e}")
+      raise
+
+  def _create_final_network_graph(self):
+    """
+    Create the final network graph with the required attributes from the final MST.
+    """
+    self.network_graph = nx.Graph()
+    node_id = 1
+    node_mapping = {}
+    for node in self.final_mst.nodes:
+      pos = node
+      if 'type' in self.final_mst.nodes[node]:
+        if self.final_mst.nodes[node]['type'] == 'building':
+          name = self.final_mst.nodes[node]['name']
+          node_type = 'building'
+        elif self.final_mst.nodes[node]['type'] == 'generation':
+          name = self.final_mst.nodes[node]['name']
+          node_type = 'generation'
+      else:
+        name = f'junction_{node_id}'
+        node_type = 'junction'
+      self.network_graph.add_node(node_id, name=name, type=node_type, pos=pos)
+      node_mapping[node] = node_id
+      node_id += 1
+    for u, v, data in self.final_mst.edges(data=True):
+      u_new = node_mapping[u]
+      v_new = node_mapping[v]
+      length = data['weight']
+      self.network_graph.add_edge(u_new, v_new, length=length)
+
+  def plot_network_graph(self):
+    """
+    Plot the network graph using matplotlib and networkx.
+    """
+    plt.figure(figsize=(15, 10))
+    pos = {node: data['pos'] for node, data in self.network_graph.nodes(data=True)}
+    nx.draw_networkx_nodes(self.network_graph, pos, node_color='blue', node_size=50)
+    nx.draw_networkx_edges(self.network_graph, pos, edge_color='gray')
+    plt.title('District Heating Network Graph')
+    plt.axis('off')
+    plt.show()

scripts/district_heating_network/pipe_data.json

@@ -0,0 +1,191 @@
+[
+  {
+    "DN": 16,
+    "inner_diameter": 16.1,
+    "outer_diameter": 21.3,
+    "thickness": 2.6,
+    "cost_per_meter": 320
+  },
+  {
+    "DN": 20,
+    "inner_diameter": 21.7,
+    "outer_diameter": 26.9,
+    "thickness": 2.6,
+    "cost_per_meter": 320
+  },
+  {
+    "DN": 25,
+    "inner_diameter": 27.3,
+    "outer_diameter": 33.7,
+    "thickness": 3.2,
+    "cost_per_meter": 320
+  },
+  {
+    "DN": 32,
+    "inner_diameter": 37.2,
+    "outer_diameter": 42.4,
+    "thickness": 2.6,
+    "cost_per_meter": 350
+  },
+  {
+    "DN": 40,
+    "inner_diameter": 43.1,
+    "outer_diameter": 48.3,
+    "thickness": 2.6,
+    "cost_per_meter": 375
+  },
+  {
+    "DN": 50,
+    "inner_diameter": 54.5,
+    "outer_diameter": 60.3,
+    "thickness": 2.9,
+    "cost_per_meter": 400
+  },
+  {
+    "DN": 65,
+    "inner_diameter": 70.3,
+    "outer_diameter": 76.1,
+    "thickness": 2.9,
+    "cost_per_meter": 450
+  },
+  {
+    "DN": 80,
+    "inner_diameter": 82.5,
+    "outer_diameter": 88.9,
+    "thickness": 3.2,
+    "cost_per_meter": 480
+  },
+  {
+    "DN": 90,
+    "inner_diameter": 100.8,
+    "outer_diameter": 108,
+    "thickness": 3.6,
+    "cost_per_meter": 480
+  },
+  {
+    "DN": 100,
+    "inner_diameter": 107.1,
+    "outer_diameter": 114.3,
+    "thickness": 3.6,
+    "cost_per_meter": 550
+  },
+  {
+    "DN": 110,
+    "inner_diameter": 125.8,
+    "outer_diameter": 133,
+    "thickness": 3.6,
+    "cost_per_meter": 550
+  },
+  {
+    "DN": 125,
+    "inner_diameter": 132.5,
+    "outer_diameter": 139.7,
+    "thickness": 3.6,
+    "cost_per_meter": 630
+  },
+  {
+    "DN": 140,
+    "inner_diameter": 151,
+    "outer_diameter": 159,
+    "thickness": 4,
+    "cost_per_meter": 700
+  },
+  {
+    "DN": 150,
+    "inner_diameter": 160.3,
+    "outer_diameter": 168.3,
+    "thickness": 4,
+    "cost_per_meter": 700
+  },
+  {
+    "DN": 180,
+    "inner_diameter": 184.7,
+    "outer_diameter": 193.7,
+    "thickness": 4.5,
+    "cost_per_meter": 700
+  },
+  {
+    "DN": 200,
+    "inner_diameter": 210.1,
+    "outer_diameter": 219.1,
+    "thickness": 4.5,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 250,
+    "inner_diameter": 263,
+    "outer_diameter": 273,
+    "thickness": 5,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 300,
+    "inner_diameter": 312.7,
+    "outer_diameter": 323.9,
+    "thickness": 5.6,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 350,
+    "inner_diameter": 344.4,
+    "outer_diameter": 355.6,
+    "thickness": 5.6,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 400,
+    "inner_diameter": 393.8,
+    "outer_diameter": 406.4,
+    "thickness": 6.3,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 450,
+    "inner_diameter": 444.4,
+    "outer_diameter": 457,
+    "thickness": 6.3,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 500,
+    "inner_diameter": 495.4,
+    "outer_diameter": 508,
+    "thickness": 6.3,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 600,
+    "inner_diameter": 595.8,
+    "outer_diameter": 610,
+    "thickness": 7.1,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 700,
+    "inner_diameter": 696.8,
+    "outer_diameter": 711,
+    "thickness": 7.1,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 800,
+    "inner_diameter": 797,
+    "outer_diameter": 813,
+    "thickness": 8,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 900,
+    "inner_diameter": 894,
+    "outer_diameter": 914,
+    "thickness": 10,
+    "cost_per_meter": 860
+  },
+  {
+    "DN": 1000,
+    "inner_diameter": 996,
+    "outer_diameter": 1016,
+    "thickness": 10,
+    "cost_per_meter": 860
+  }
+]

scripts/district_heating_network/road_processor.py

@@ -0,0 +1,56 @@
+from pathlib import Path
+from shapely.geometry import Polygon, Point, shape
+import json
+
+
+def road_processor(x, y, diff):
+    """
+    Processes a .JSON file to find roads that have at least one node within a specified polygon.
+
+    Parameters:
+    x (float): The x-coordinate of the center of the selection box.
+    y (float): The y-coordinate of the center of the selection box.
+    diff (float): The half-width of the selection box.
+
+    Returns:
+    str: The file path of the output GeoJSON file containing the selected roads.
+    """
+    diff += 2 * diff
+    # Define the selection polygon
+    selection_box = Polygon([
+        [x + diff, y - diff],
+        [x - diff, y - diff],
+        [x - diff, y + diff],
+        [x + diff, y + diff]
+    ])
+
+    # Define input and output file paths
+    geojson_file = Path("./input_files/roads.json").resolve()
+    output_file = Path('./input_files/output_roads.geojson').resolve()
+
+    # Initialize a list to store the roads in the region
+    roads_in_region = []
+
+    # Read the GeoJSON file
+    with open(geojson_file, 'r') as file:
+        roads = json.load(file)
+        line_features = [feature for feature in roads['features'] if feature['geometry']['type'] == 'LineString']
+
+        # Check each road feature
+        for feature in line_features:
+            road_shape = shape(feature['geometry'])
+            # Check if any node of the road is inside the selection box
+            if any(selection_box.contains(Point(coord)) for coord in road_shape.coords):
+                roads_in_region.append(feature)
+
+    # Create a new GeoJSON structure with the selected roads
+    output_geojson = {
+        "type": "FeatureCollection",
+        "features": roads_in_region
+    }
+
+    # Write the selected roads to the output file
+    with open(output_file, 'w') as outfile:
+        json.dump(output_geojson, outfile)
+
+    return output_file

scripts/district_heating_network/simultinity_factor.py

@@ -0,0 +1,80 @@
+import pandas as pd
+import numpy as np
+
+class DemandShiftProcessor:
+  def __init__(self, city):
+    self.city = city
+
+  def random_shift(self, series):
+    shift_amount = np.random.randint(0, round(0.005 * len(series)))
+    return series.shift(shift_amount).fillna(series.shift(shift_amount - len(series)))
+
+  def process_demands(self):
+    heating_dfs = []
+    cooling_dfs = []
+
+    for building in self.city.buildings:
+      heating_df = self.convert_building_to_dataframe(building, 'heating')
+      cooling_df = self.convert_building_to_dataframe(building, 'cooling')
+      heating_df.set_index('Date/Time', inplace=True)
+      cooling_df.set_index('Date/Time', inplace=True)
+      shifted_heating_demands = heating_df.apply(self.random_shift, axis=0)
+      shifted_cooling_demands = cooling_df.apply(self.random_shift, axis=0)
+      self.update_building_demands(building, shifted_heating_demands, 'heating')
+      self.update_building_demands(building, shifted_cooling_demands, 'cooling')
+      heating_dfs.append(shifted_heating_demands)
+      cooling_dfs.append(shifted_cooling_demands)
+
+    combined_heating_df = pd.concat(heating_dfs, axis=1)
+    combined_cooling_df = pd.concat(cooling_dfs, axis=1)
+    self.calculate_and_set_simultaneity_factor(combined_heating_df, 'heating')
+    self.calculate_and_set_simultaneity_factor(combined_cooling_df, 'cooling')
+
+  def convert_building_to_dataframe(self, building, demand_type):
+    if demand_type == 'heating':
+      data = {
+        "Date/Time": self.generate_date_time_index(),
+        "Heating_Demand": building.heating_demand["hour"]
+      }
+    else:  # cooling
+      data = {
+        "Date/Time": self.generate_date_time_index(),
+        "Cooling_Demand": building.cooling_demand["hour"]
+      }
+    return pd.DataFrame(data)
+
+  def generate_date_time_index(self):
+    # Generate hourly date time index for a full year in 2024
+    date_range = pd.date_range(start="2013-01-01 00:00:00", end="2013-12-31 23:00:00", freq='H')
+    return date_range.strftime('%m/%d %H:%M:%S').tolist()
+
+  def update_building_demands(self, building, shifted_demands, demand_type):
+    if demand_type == 'heating':
+      shifted_series = shifted_demands["Heating_Demand"]
+      building.heating_demand = self.calculate_new_demands(shifted_series)
+    else:  # cooling
+      shifted_series = shifted_demands["Cooling_Demand"]
+      building.cooling_demand = self.calculate_new_demands(shifted_series)
+
+  def calculate_new_demands(self, shifted_series):
+    new_demand = {
+      "hour": shifted_series.tolist(),
+      "month": self.calculate_monthly_demand(shifted_series),
+      "year": [shifted_series.sum()]
+    }
+    return new_demand
+
+  def calculate_monthly_demand(self, series):
+    series.index = pd.to_datetime(series.index, format='%m/%d %H:%M:%S')
+    monthly_demand = series.resample('M').sum()
+    return monthly_demand.tolist()
+
+  def calculate_and_set_simultaneity_factor(self, combined_df, demand_type):
+    total_demand_original = combined_df.sum(axis=1)
+    peak_total_demand_original = total_demand_original.max()
+    individual_peak_demands = combined_df.max(axis=0)
+    sum_individual_peak_demands = individual_peak_demands.sum()
+    if demand_type == 'heating':
+      self.city.simultaneity_factor_heating = peak_total_demand_original / sum_individual_peak_demands
+    else:  # cooling
+      self.city.simultaneity_factor_cooling = peak_total_demand_original / sum_individual_peak_demands

scripts/optimization/energy_system_sizing_optimization.py

@@ -0,0 +1,26 @@
+import random
+import numpy as np
+import hub.helpers.constants as cte
+
+class EnergySystemOptimizer:
+  def __init__(self, building, objectives, population_size=20, number_of_generations=100, crossover_rate=0.8,
+               mutation_rate=0.1):
+    self.building = building
+    self.objectives = objectives
+    self.population_size = population_size
+    self.num_gen = number_of_generations
+    self.crossover_rate = crossover_rate
+    self.mutation_rate = mutation_rate
+
+# Initialize population
+  def initialize_population(self):
+    population = []
+    for _ in range(self.population_size):
+      individual = [
+        random.uniform(0.1, 10.0),  # hp_size
+        random.uniform(0.1, 10.0),  # boiler_size
+        random.uniform(0.1, 10.0),  # tes_size
+        random.uniform(40, 60)  # cutoff_temp
+      ]
+      population.append(individual)
+    return population

scripts/optimization/objectives.py

@@ -0,0 +1,16 @@
+# Objective functions
+MINIMUM_LCC = 1
+MINIMUM_EC = 2
+MINIMUM_EMISSIONS = 3
+MINIMUM_LCC_MINIMUM_EC = 4
+MINIMUM_LCC_MINIMUM_EMISSIONS = 5
+MINIMUM_EC_MINIMUM_EMISSIONS = 6
+MINIMUM_LCC_MINIMUM_EC_MINIMUM_EMISSIONS = 7
+OBJECTIVE_FUNCTIONS = [MINIMUM_LCC,
+                       MINIMUM_EC,
+                       MINIMUM_EMISSIONS,
+                       MINIMUM_LCC_MINIMUM_EC,
+                       MINIMUM_LCC_MINIMUM_EMISSIONS,
+                       MINIMUM_EC_MINIMUM_EMISSIONS,
+                       MINIMUM_LCC_MINIMUM_EC_MINIMUM_EMISSIONS]
+

scripts/random_assignation.py

@@ -28,8 +28,8 @@ residential_systems_percentage = {'system 1 gas': 15,
                                   'system 8 gas': 15,
                                   'system 8 electricity': 35}
 
-residential_new_systems_percentage = {'PV+ASHP+GasBoiler+TES': 100,
-                                      'PV+4Pipe+DHW': 0,
+residential_new_systems_percentage = {'PV+ASHP+GasBoiler+TES': 0,
+                                      'PV+4Pipe+DHW': 100,
                                       'PV+ASHP+ElectricBoiler+TES': 0,
                                       'PV+GSHP+GasBoiler+TES': 0,
                                       'PV+GSHP+ElectricBoiler+TES': 0,

scripts/system_simulation_models/archetype13.py

@@ -5,7 +5,7 @@ from hub.helpers.monthly_values import MonthlyValues
 
 
 class Archetype13:
-  def __init__(self, building, output_path):
+  def __init__(self, building, output_path, csv_output=True,):
     self._building = building
     self._name = building.name
     self._pv_system = building.energy_systems[0]
@@ -25,6 +25,7 @@ class Archetype13:
     self._t_out = building.external_temperature[cte.HOUR]
     self.results = {}
     self.dt = 900
+    self.csv_output = csv_output
 
   def hvac_sizing(self):
     storage_factor = 3
@@ -374,10 +375,11 @@ class Archetype13:
       MonthlyValues.get_total_month(self._building.domestic_hot_water_consumption[cte.HOUR]))
     self._building.domestic_hot_water_consumption[cte.YEAR] = [
       sum(self._building.domestic_hot_water_consumption[cte.MONTH])]
-    file_name = f'energy_system_simulation_results_{self._name}.csv'
-    with open(self._output_path / file_name, 'w', newline='') as csvfile:
-      output_file = csv.writer(csvfile)
-      # Write header
-      output_file.writerow(self.results.keys())
-      # Write data
-      output_file.writerows(zip(*self.results.values()))
+    if self.csv_output:
+      file_name = f'energy_system_simulation_results_{self._name}.csv'
+      with open(self._output_path / file_name, 'w', newline='') as csvfile:
+        output_file = csv.writer(csvfile)
+        # Write header
+        output_file.writerow(self.results.keys())
+        # Write data
+        output_file.writerows(zip(*self.results.values()))