fix: starting to test the workflow for finding the system sizing period

energy_system_modelling_package/energy_system_modelling_factories/hvac_dhw_systems_simulation_models/domestic_hot_water_heat_pump_with_tes.py

@@ -6,7 +6,7 @@ from hub.helpers.monthly_values import MonthlyValues
 
 class DomesticHotWaterHeatPumpTes:
   def __init__(self, hp, tes, hourly_dhw_demand_joules, upper_limit_tes,
-               outdoor_temperature, dt=900):
+               outdoor_temperature, dt=None):
     self.hp = hp
     self.tes = tes
     self.dhw_demand = [demand / cte.WATTS_HOUR_TO_JULES for demand in hourly_dhw_demand_joules]
@@ -30,8 +30,8 @@ class DomesticHotWaterHeatPumpTes:
     hp_delta_t = 8
     number_of_ts = int(cte.HOUR_TO_SECONDS / self.dt)
     source_temperature_hourly = self.hp_characteristics.hp_source_temperature()
-    source_temperature = [x for x in source_temperature_hourly for _ in range(number_of_ts)]
-    demand = [x for x in self.dhw_demand for _ in range(number_of_ts)]
+    source_temperature = [0] + [x for x in source_temperature_hourly for _ in range(number_of_ts)]
+    demand = [0] + [x for x in self.dhw_demand for _ in range(number_of_ts)]
     variable_names = ["t_sup_hp", "t_tank", "m_ch", "m_dis", "q_hp", "q_coil", "hp_cop",
                       "hp_electricity", "available hot water (m3)", "refill flow rate (kg/s)", "total_consumption"]
     num_hours = len(demand)

energy_system_modelling_package/energy_system_modelling_factories/hvac_dhw_systems_simulation_models/heat_pump_boiler_tes_heating.py

@@ -8,12 +8,15 @@ from energy_system_modelling_package.energy_system_modelling_factories.hvac_dhw_
 
 class HeatPumpBoilerTesHeating:
   def __init__(self, hp, boiler, tes, hourly_heating_demand_joules, heating_peak_load_watts, upper_limit_tes,
-               outdoor_temperature, dt=900):
+               outdoor_temperature, dt=None):
     self.hp = hp
     self.boiler = boiler
     self.tes = tes
     self.heating_demand = [demand / cte.WATTS_HOUR_TO_JULES for demand in hourly_heating_demand_joules]
-    self.heating_peak_load = heating_peak_load_watts
+    if heating_peak_load_watts is not None:
+      self.heating_peak_load = heating_peak_load_watts
+    else:
+      self.heating_peak_load = max(hourly_heating_demand_joules) / cte.HOUR_TO_SECONDS
     self.upper_limit_tes = upper_limit_tes
     self.hp_characteristics = HeatPump(self.hp, outdoor_temperature)
     self.t_out = outdoor_temperature
@@ -130,23 +133,25 @@ class HeatPumpBoilerTesHeating:
         coil_sum = 0
     hp.energy_consumption[cte.HEATING] = {}
     hp.energy_consumption[cte.HEATING][cte.HOUR] = hp_hourly
-    hp.energy_consumption[cte.HEATING][cte.MONTH] = MonthlyValues.get_total_month(
-      hp.energy_consumption[cte.HEATING][cte.HOUR])
-    hp.energy_consumption[cte.HEATING][cte.YEAR] = [
-      sum(hp.energy_consumption[cte.HEATING][cte.MONTH])]
     boiler.energy_consumption[cte.HEATING] = {}
     boiler.energy_consumption[cte.HEATING][cte.HOUR] = boiler_hourly
-    boiler.energy_consumption[cte.HEATING][cte.MONTH] = MonthlyValues.get_total_month(
-      boiler.energy_consumption[cte.HEATING][cte.HOUR])
-    boiler.energy_consumption[cte.HEATING][cte.YEAR] = [
-      sum(boiler.energy_consumption[cte.HEATING][cte.MONTH])]
+    if len(self.heating_demand) == 8760:
+      hp.energy_consumption[cte.HEATING][cte.MONTH] = MonthlyValues.get_total_month(
+        hp.energy_consumption[cte.HEATING][cte.HOUR])
+      hp.energy_consumption[cte.HEATING][cte.YEAR] = [
+        sum(hp.energy_consumption[cte.HEATING][cte.MONTH])]
+      boiler.energy_consumption[cte.HEATING][cte.MONTH] = MonthlyValues.get_total_month(
+        boiler.energy_consumption[cte.HEATING][cte.HOUR])
+      boiler.energy_consumption[cte.HEATING][cte.YEAR] = [
+        sum(boiler.energy_consumption[cte.HEATING][cte.MONTH])]
     if tes.heating_coil_capacity is not None:
       tes.heating_coil_energy_consumption[cte.HEATING] = {}
-      tes.heating_coil_energy_consumption[cte.HEATING][cte.HOUR] = coil_hourly
-      tes.heating_coil_energy_consumption[cte.HEATING][cte.MONTH] = MonthlyValues.get_total_month(
-        tes.heating_coil_energy_consumption[cte.HEATING][cte.HOUR])
-      tes.heating_coil_energy_consumption[cte.HEATING][cte.YEAR] = [
-        sum(tes.heating_coil_energy_consumption[cte.HEATING][cte.MONTH])]
+      if len(self.heating_demand) == 8760:
+        tes.heating_coil_energy_consumption[cte.HEATING][cte.HOUR] = coil_hourly
+        tes.heating_coil_energy_consumption[cte.HEATING][cte.MONTH] = MonthlyValues.get_total_month(
+          tes.heating_coil_energy_consumption[cte.HEATING][cte.HOUR])
+        tes.heating_coil_energy_consumption[cte.HEATING][cte.YEAR] = [
+          sum(tes.heating_coil_energy_consumption[cte.HEATING][cte.MONTH])]
     self.results['Heating Demand (W)'] = demand
     self.results['HP Heat Output (W)'] = q_hp
     self.results['HP Source Temperature'] = source_temperature

energy_system_modelling_package/energy_system_modelling_factories/system_sizing_methods/genetic_algorithm/genetic_algorithm.py

@@ -10,46 +10,178 @@ from energy_system_modelling_package.energy_system_modelling_factories.system_si
 
 
 class GeneticAlgorithm:
-  def __init__(self, city, population_size=100, generations=50, crossover_rate=0.8, mutation_rate=0.1,
+  def __init__(self, population_size=50, generations=50, crossover_rate=0.8, mutation_rate=0.1,
                optimization_scenario='cost', output_path=None):
-    self.city = city
     self.population_size = population_size
+    self.population = []
     self.generations = generations
     self.crossover_rate = crossover_rate
     self.mutation_rate = mutation_rate
     self.optimization_scenario = optimization_scenario
     self.list_of_solutions = []
+    self.best_solution = None
+    self.best_solution_generation = None
     self.output_path = output_path
 
   def energy_system_archetype_optimal_sizing(self):
     pass
 
-# Initialize Population
-  def initialize_population(self, energy_system):
+  # Initialize Population
+  def initialize_population(self, building, energy_system):
     """
     The population to optimize the sizes of generation and storage components of any energy system are optimized.
     The initial population will be a list of individuals where each of them represent a possible solution to the
     problem.
+    :param building:
     :param energy_system:
     :return: List
     """
-    population = []
     for i in range(self.population_size):
-      population.append(Individual(energy_system, self.optimization_scenario))
-    return population
+      self.population.append(Individual(building, energy_system, self.optimization_scenario))
 
-  def solve_ga(self, energy_system):
+  def order_population(self):
     """
-    Solving GA for a single energy system. Here are the steps:
-    1- Population is initialized using the "initialize_population" method in this class.
-    :param energy_system:
+    ordering the population based on the fitness score in ascending order
     :return:
     """
-    energy_system = energy_system
-    best_individuals = []
-    best_solution = None
-    population = self.initialize_population(energy_system)
-
-
+    self.population = sorted(self.population, key=lambda x: x.individual['fitness_score'])
+
+  def tournament_selection(self):
+    selected = []
+    for _ in range(len(self.population)):
+      i, j = random.sample(range(self.population_size), 2)
+      if self.population[i].individual['fitness_score'] < self.population[j].individual['fitness_score']:
+        selected.append(copy.deepcopy(self.population[i]))
+      else:
+        selected.append(copy.deepcopy(self.population[j]))
+    return selected
+
+  def crossover(self, parent1, parent2):
+    """
+    Crossover between two parents to produce two children.
+
+    Swaps generation components and storage components between the two parents with a 50% chance.
+
+    :param parent1: First parent individual.
+    :param parent2: Second parent individual.
+    :return: Two child individuals (child1 and child2).
+    """
+    if random.random() < self.crossover_rate:
+      # Deep copy of the parents to create children
+      child1, child2 = copy.deepcopy(parent1), copy.deepcopy(parent2)
+
+      # Crossover for Generation Components
+      for i in range(len(parent1.individual['Generation Components'])):
+        if random.random() < 0.5:
+          # Swap the entire generation component
+          child1.individual['Generation Components'][i], child2.individual['Generation Components'][i] = (
+            child2.individual['Generation Components'][i],
+            child1.individual['Generation Components'][i]
+          )
+
+      # Crossover for Energy Storage Components
+      for i in range(len(parent1.individual['Energy Storage Components'])):
+        if random.random() < 0.5:
+          # Swap the entire storage component
+          child1.individual['Energy Storage Components'][i], child2.individual['Energy Storage Components'][i] = (
+            child2.individual['Energy Storage Components'][i],
+            child1.individual['Energy Storage Components'][i]
+          )
+
+      return child1, child2
+    else:
+      # If crossover does not happen, return copies of the original parents
+      return copy.deepcopy(parent1), copy.deepcopy(parent2)
+
+  def mutate(self, individual, building):
+    """
+    Mutates the individual's generation and storage components.
+
+    - `individual`: The individual to mutate (contains generation and storage components).
+    - `building`: Building data that contains constraints such as peak heating load and available space.
+
+    Returns the mutated individual.
+    """
+    # Mutate Generation Components
+    for generation_component in individual['Generation Components']:
+      if random.random() < self.mutation_rate:
+        if generation_component['nominal_heating_efficiency'] is not None:
+          # Mutate heating capacity
+          generation_component['heating_capacity'] = random.uniform(
+            0, building.heating_peak_load[cte.YEAR][0])
+        if generation_component['nominal_cooling_efficiency'] is not None:
+          # Mutate cooling capacity
+          generation_component['cooling_capacity'] = random.uniform(
+            0, building.cooling_peak_load[cte.YEAR][0])
+
+    # Mutate Storage Components
+    for storage_component in individual['Energy Storage Components']:
+      if random.random() < self.mutation_rate:
+        if storage_component['type'] == f'{cte.THERMAL}_storage':
+          # Mutate the volume of thermal storage
+          max_available_space = 0.01 * building.volume / building.storeys_above_ground
+          storage_component['volume'] = random.uniform(0, max_available_space)
+          if storage_component['heating_coil_capacity'] is not None:
+            storage_component['volume'] = random.uniform(0, building.heating_peak_load[cte.YEAR][0])
+
+    return individual
+
+  def solve_ga(self, building, energy_system):
+    """
+    Solving GA for a single energy system. Here are the steps:
+    1. Initialize population using the "initialize_population" method in this class.
+    2. Evaluate the initial population using the "score_evaluation" method in the Individual class.
+    3. Sort population based on fitness score.
+    4. Repeat selection, crossover, and mutation for a fixed number of generations.
+    5. Track the best solution found during the optimization process.
+
+    :param building: Building object for the energy system.
+    :param energy_system: Energy system to optimize.
+    :return: Best solution after running the GA.
+    """
+    # Step 1: Initialize the population
+    self.initialize_population(building, energy_system)
+    # Step 2: Evaluate the initial population
+    for individual in self.population:
+      individual.score_evaluation()
+    # Step 3: Order population based on fitness scores
+    self.order_population()
+    # Track the best solution
+    self.best_solution = self.population[0]
+    self.list_of_solutions.append(copy.deepcopy(self.best_solution.individual))
+    # Step 4: Run GA for a fixed number of generations
+    for generation in range(self.generations):
+      print(f"Generation {generation}")
+    # Selection (using tournament selection)
+      selected_population = self.tournament_selection()
+    # Create the next generation through crossover and mutation
+      next_population = []
+      for i in range(0, self.population_size, 2):
+        parent1 = selected_population[i]
+        parent2 = selected_population[i + 1] if (i + 1) < len(selected_population) else selected_population[0]
+    # Step 5: Apply crossover
+        child1, child2 = self.crossover(parent1, parent2)
+    # Step 6: Apply mutation
+        self.mutate(child1.individual, building)
+        self.mutate(child2.individual, building)
+    # Step 7: Evaluate the children
+        child1.score_evaluation()
+        child2.score_evaluation()
+        next_population.extend([child1, child2])
+    # Replace old population with the new one
+      self.population = next_population
+    # Step 8: Sort the new population based on fitness
+      self.order_population()
+    # Track the best solution found in this generation
+      if self.population[0].individual['fitness_score'] < self.best_solution.individual['fitness_score']:
+        self.best_solution = self.population[0]
+        self.best_solution_generation = generation
+    # Store the best solution in the list of solutions
+      self.list_of_solutions.append(copy.deepcopy(self.population[0]))
+    print(f"Best solution found in generation {self.best_solution_generation}")
+    print(f"Best solution: {self.best_solution.individual}")
+
+    # Return the best solution after running GA
+    return self.best_solution
 
 

energy_system_modelling_package/energy_system_modelling_factories/system_sizing_methods/genetic_algorithm/individual.py

@@ -1,43 +1,50 @@
+import math
 import random
-
 from hub.helpers.dictionaries import Dictionaries
-
 from hub.catalog_factories.costs_catalog_factory import CostsCatalogFactory
-
 import hub.helpers.constants as cte
 from energy_system_modelling_package.energy_system_modelling_factories.hvac_dhw_systems_simulation_models.domestic_hot_water_heat_pump_with_tes import \
   DomesticHotWaterHeatPumpTes
 from energy_system_modelling_package.energy_system_modelling_factories.hvac_dhw_systems_simulation_models.heat_pump_boiler_tes_heating import \
   HeatPumpBoilerTesHeating
+import numpy_financial as npf
 
 
 class Individual:
-  def __init__(self, building, energy_system, energy_demands, optimization_scenario, heating_design_load=None,
-               cooling_design_load=None, available_space=None, dt=3600, electricity_tariff_type='fixed'):
+  def __init__(self, building, energy_system, design_period_energy_demands, design_period_outdoor_temperature,
+               optimization_scenario,heating_design_load=None, cooling_design_load=None, dt=900, fuel_price_index=0.05,
+               electricity_tariff_type='fixed', consumer_price_index=0.04, interest_rate=0.04,
+               discount_rate=0.03, percentage_credit=0, credit_years=15):
     """
     :param building: building object
     :param energy_system: energy system to be optimized
-    :param energy_demands: energy demands that the system should supply. A dictionary in the format of
-    {'demand_type1':[hourly_values], 'demand_type2':[hourly_values]} must be provided. all values must be in Joules.
+    :param design_period_energy_demands: A dictionary demand types as its keys and hourly values of the day with the
+    maximum daily demand and the days before and after it.
     :param optimization_scenario: a string indicating the objective function from minimization of cost,
     energy consumption, and both together
     :param heating_design_load: heating design load in W
     :param cooling_design_load: cooling design load in W
-    :param available_space: available space for the system. currently it is only used to have a constraint on the
-    maximum size of the storage tank
     :param dt the time step size used for simulations
+    :param fuel_price_index the price increase index of all fuels. A single value is used for all fuels.
     :param electricity_tariff_type the electricity tariff type between 'fixed' and 'variable' for economic optimization
+    :param consumer_price_index
     """
     self.building = building
     self.energy_system = energy_system
+    self.design_period_energy_demands = design_period_energy_demands
     self.demand_types = energy_system.demand_types
-    self.energy_demands = energy_demands
     self.optimization_scenario = optimization_scenario
     self.heating_design_load = heating_design_load
     self.cooling_design_load = cooling_design_load
-    self.available_space = available_space
+    self.available_space = building.volume / building.storeys_above_ground
     self.dt = dt
+    self.fuel_price_index = fuel_price_index
     self.electricity_tariff_type = electricity_tariff_type
+    self.consumer_price_index = consumer_price_index
+    self.interest_rate = interest_rate
+    self.discount_rate = discount_rate
+    self.credit_years = credit_years
+    self.percentage_credit = percentage_credit
     self.costs = self.costs_archetype()
     self.feasibility = True
     self.fitness_score = 0
@@ -66,14 +73,7 @@ class Individual:
         cooling_capacity = 0
         heat_efficiency = generation_component.heat_efficiency
         cooling_efficiency = generation_component.cooling_efficiency
-        unit_fuel_cost = 0
-        if generation_component.fuel_type == cte.GAS:
-          gas_cost = self.fuel_cost(generation_component.fuel_type, 'fixed')
-          unit_fuel_cost = gas_cost.variable.values[0]
-        elif generation_component.fuel_type == cte.ELECTRICITY:
-          electricity_cost = self.fuel_cost(generation_component.fuel_type, self.electricity_tariff_type)
-          unit_fuel_cost = electricity_cost.variable.values
-
+        unit_fuel_cost = self.fuel_cost_per_kwh(generation_component.fuel_type, 'fixed')
       self.individual['Generation Components'].append({
         'type': generation_component.system_type,
         'heating_capacity': heating_capacity,
@@ -85,7 +85,7 @@ class Individual:
         'fuel_type': generation_component.fuel_type,
         'unit_investment_cost': investment_cost,
         'unit_reposition_cost': reposition_cost,
-        'unit_fuel_cost': unit_fuel_cost,
+        'unit_fuel_cost(CAD/kWh)': unit_fuel_cost,
         'unit_maintenance_cost': maintenance_cost,
         'lifetime': lifetime
       })
@@ -96,10 +96,12 @@ class Individual:
                                       f'{energy_storage_system.type_energy_stored}_storage'))
           if energy_storage_system.type_energy_stored == cte.THERMAL:
             heating_coil_capacity = energy_storage_system.heating_coil_capacity
+            heating_coil_fuel_cost = self.fuel_cost_per_kwh(f'{cte.ELECTRICITY}', 'fixed')
             volume = 0
             capacity = None
           else:
             heating_coil_capacity = None
+            heating_coil_fuel_cost = None
             volume = None
             capacity = 0
           self.individual['Energy Storage Components'].append(
@@ -109,10 +111,10 @@ class Individual:
              'heating_coil_capacity': heating_coil_capacity,
              'unit_investment_cost': investment_cost,
              'unit_reposition_cost': reposition_cost,
+             'heating_coil_fuel_cost': heating_coil_fuel_cost,
              'unit_maintenance_cost': 0,
              'lifetime': lifetime
              })
-    print(self.individual)
 
   def initialization(self):
     """
@@ -124,29 +126,74 @@ class Individual:
     generation_components = self.individual['Generation Components']
     storage_components = self.individual['Energy Storage Components']
     for generation_component in generation_components:
-      if generation_component['nominal_heating_efficiency'] is not None and self.heating_design_load is not None:
-        generation_component['heating_capacity'] = random.uniform(0, self.heating_design_load)
+      if generation_component['nominal_heating_efficiency'] is not None:
+        if self.heating_design_load is not None:
+          generation_component['heating_capacity'] = random.uniform(0, self.heating_design_load)
+        else:
+          generation_component['heating_capacity'] = random.uniform(0, max(self.design_period_energy_demands[cte.HEATING]))
       else:
         generation_component['heating_capacity'] = None
-      if generation_component['nominal_cooling_efficiency'] is not None and self.cooling_design_load is not None:
-        generation_component['cooling_capacity'] = random.uniform(0, self.cooling_design_load)
+      if generation_component['nominal_cooling_efficiency'] is not None:
+        if self.cooling_design_load is not None:
+          generation_component['cooling_capacity'] = random.uniform(0, self.cooling_design_load)
+        else:
+          generation_component['cooling_capacity'] = random.uniform(0, self.design_period_energy_demands[cte.COOLING])
       else:
         generation_component['cooling_capacity'] = None
       if generation_component['nominal_electricity_efficiency'] is None:
         generation_component['electricity_capacity'] = None
     for storage_component in storage_components:
-      if storage_component['type'] == cte.THERMAL:
-        storage_component['volume'] = 0.05 * self.available_space
+      if storage_component['type'] == f'{cte.THERMAL}_storage':
+        storage_component['volume'] = random.uniform(0, 0.01 * self.available_space)
+        if storage_component['heating_coil_capacity'] is not None:
+          if self.heating_design_load is not None:
+            storage_component['heating_coil_capacity'] = random.uniform(0, self.heating_design_load)
+          else:
+            storage_component['heating_coil_capacity'] = random.uniform(0, self.design_period_energy_demands[cte.HEATING])
 
   def score_evaluation(self):
-    results = self.system_simulation()
+    self.initialization()
+    self.system_simulation()
+    self.individual['feasible'] = self.feasibility
     if self.feasibility:
-      if self.optimization_scenario == 'cost':
-        pass
-      elif self.optimization_scenario == 'energy_consumption':
-        pass
+      if 'cost' in self.optimization_scenario:
+        investment_cost = 0
+        operation_cost_year_0 = 0
+        maintenance_cost_year_0 = 0
+        for generation_system in self.individual['Generation Components']:
+          heating_capacity = 0 if generation_system['heating_capacity'] is None else generation_system[
+            'heating_capacity']
+          cooling_capacity = 0 if generation_system['cooling_capacity'] is None else generation_system[
+            'cooling_capacity']
+          capacity = max(heating_capacity, cooling_capacity)
+          investment_cost += capacity * generation_system['unit_investment_cost'] / 1000
+          maintenance_cost_year_0 += capacity * generation_system['unit_maintenance_cost'] / 1000
+          operation_cost_year_0 += (generation_system['yearly_energy_consumption(kWh)'] *
+                                    generation_system['unit_fuel_cost(CAD/kWh)'])
+        for storage_system in self.individual['Energy Storage Components']:
+          if cte.THERMAL in storage_system['type']:
+            investment_cost += storage_system['volume'] * storage_system['unit_investment_cost']
+            if storage_system['heating_coil_capacity'] is not None:
+              operation_cost_year_0 += (storage_system['yearly_energy_consumption(kWh)'] *
+                                        storage_system['heating_coil_fuel_cost'])
+        lcc = self.life_cycle_cost_calculation(investment_cost=investment_cost,
+                                               operation_cost_year_0=operation_cost_year_0,
+                                               maintenance_cost_year_0=maintenance_cost_year_0)
+        self.individual['lcc'] = lcc
+        self.individual['fitness_score'] = lcc
+      elif 'energy_consumption' in self.optimization_scenario:
+        total_energy_consumption = 0
+        for generation_system in self.individual['Generation Components']:
+          total_energy_consumption += generation_system['yearly_energy_consumption(kWh)']
+        for storage_system in self.individual['Energy Storage Components']:
+          if cte.THERMAL in storage_system['type'] and storage_system['heating_coil_capacity'] is not None:
+            total_energy_consumption += storage_system['yearly_energy_consumption(kWh)']
+        self.individual['total_energy_consumption'] = total_energy_consumption
+        self.individual['fitness_score'] = total_energy_consumption
       elif self.optimization_scenario == 'cost_energy_consumption':
         pass
+    else:
+      self.individual['fitness_score'] = float('inf')
 
   def system_simulation(self):
     """
@@ -154,14 +201,20 @@ class Individual:
     Based on cluster id and demands, model is imported and run.
     :return: dictionary of results
     """
-    results = None
-    if self.building.energy_systems_archetype_cluster_id == 1:
+    if self.building.energy_systems_archetype_cluster_id == '1':
       if cte.HEATING in self.demand_types:
         boiler = self.energy_system.generation_systems[0]
+        boiler.nominal_heat_output = self.individual['Generation Components'][0]['heating_capacity']
         hp = self.energy_system.generation_systems[1]
+        hp.nominal_heat_output = self.individual['Generation Components'][1]['heating_capacity']
         tes = self.energy_system.generation_systems[0].energy_storage_systems[0]
-        heating_demand_joules = self.energy_demands[f'{cte.HEATING}']
-        heating_peak_load_watts = self.heating_design_load
+        tes.volume = self.individual['Energy Storage Components'][0]['volume']
+        tes.height = self.building.average_storey_height - 1
+        tes.heating_coil_capacity = self.individual['Energy Storage Components'][0]['heating_coil_capacity'] \
+            if self.individual['Energy Storage Components'][0]['heating_coil_capacity'] is not None else None
+        heating_demand_joules = self.design_period_energy_demands[cte.HEATING]
+        heating_peak_load_watts = max(self.design_period_energy_demands[cte.HEATING]) if (self.heating_design_load is not
+                                                               None) else self.building.heating_peak_load[cte.YEAR][0]
         upper_limit_tes_heating = 55
         outdoor_temperature = self.building.external_temperature[cte.HOUR]
         results = HeatPumpBoilerTesHeating(hp=hp,
@@ -176,8 +229,13 @@ class Individual:
           self.feasibility = False
       elif cte.DOMESTIC_HOT_WATER in self.demand_types:
         hp = self.energy_system.generation_systems[0]
+        hp.nominal_heat_output = self.individual['Generation Components'][0]['heating_capacity']
         tes = self.energy_system.generation_systems[0].energy_storage_systems[0]
-        dhw_demand_joules = self.energy_demands[f'{cte.DOMESTIC_HOT_WATER}']
+        tes.volume = self.individual['Energy Storage Components'][0]['volume']
+        tes.height = self.building.average_storey_height - 1
+        tes.heating_coil_capacity = self.individual['Energy Storage Components'][0]['heating_coil_capacity'] \
+            if self.individual['Energy Storage Components'][0]['heating_coil_capacity'] is not None else None
+        dhw_demand_joules = self.design_period_energy_demands[cte.DOMESTIC_HOT_WATER]
         upper_limit_tes = 65
         outdoor_temperature = self.building.external_temperature[cte.HOUR]
         results = DomesticHotWaterHeatPumpTes(hp=hp,
@@ -188,9 +246,30 @@ class Individual:
                                               dt=self.dt).simulation()
         if min(results['DHW TES Temperature']) < 55:
           self.feasibility = False
-    return results
+      if self.feasibility:
+        generation_system_types = [generation_system.system_type for generation_system in
+                                   self.energy_system.generation_systems]
+        for generation_component in self.individual['Generation Components']:
+          if generation_component['type'] in generation_system_types:
+            index = generation_system_types.index(generation_component['type'])
+            for demand_type in self.demand_types:
+              if demand_type in self.energy_system.generation_systems[index].energy_consumption:
+                generation_component['yearly_energy_consumption(kWh)'] = (
+                    self.energy_system.generation_systems[index].energy_consumption[demand_type][cte.YEAR][0] / 3.6e6)
+        for storage_component in self.individual['Energy Storage Components']:
+          if storage_component['type'] == f'{cte.THERMAL}_storage' and storage_component[
+            'heating_coil_capacity'] is not None:
+            for generation_system in self.energy_system.generation_systems:
+              if generation_system.energy_storage_systems is not None:
+                for storage_system in generation_system.energy_storage_systems:
+                  if storage_system.type_energy_stored == cte.THERMAL:
+                    for demand_type in self.demand_types:
+                      if demand_type in storage_system.heating_coil_energy_consumption:
+                        storage_component['yearly_energy_consumption(kWh)'] = (
+                            storage_system.heating_coil_energy_consumption[demand_type][cte.YEAR][0] / 3.6e6)
 
-  def life_cycle_cost_calculation(self, life_cycle_duration=41):
+  def life_cycle_cost_calculation(self, investment_cost, operation_cost_year_0, maintenance_cost_year_0,
+                                  life_cycle_duration=41):
     """
     Calculating the life cycle cost of the energy system. The original cost workflow in hub is not used to reduce
     computation time.Here are the steps:
@@ -198,10 +277,55 @@ class Individual:
         2- Capital costs (investment and reposition) are calculated and appended to a list to have the capital cash
            flow.
         3-
+    :param maintenance_cost_year_0:
+    :param operation_cost_year_0:
+    :param investment_cost:
     :param life_cycle_duration:
     :return:
     """
-    return 0
+    capital_costs_cash_flow = [investment_cost]
+    operational_costs_cash_flow = [0]
+    maintenance_costs_cash_flow = [0]
+    end_of_life_costs = [0] * (life_cycle_duration + 1)
+    for i in range(1, life_cycle_duration + 1):
+      yearly_operational_cost = math.pow(1 + self.fuel_price_index, i) * operation_cost_year_0
+      yearly_maintenance_cost = math.pow(1 + self.consumer_price_index, i) * maintenance_cost_year_0
+      yearly_capital_cost = 0
+      for generation_system in self.individual['Generation Components']:
+        if (i % generation_system['lifetime']) == 0 and i != (life_cycle_duration - 1):
+          cost_increase = math.pow(1 + self.consumer_price_index, i)
+          heating_capacity = 0 if generation_system['heating_capacity'] is None else generation_system[
+            'heating_capacity']
+          cooling_capacity = 0 if generation_system['cooling_capacity'] is None else generation_system[
+            'cooling_capacity']
+          capacity = max(heating_capacity, cooling_capacity)
+          yearly_capital_cost += generation_system['unit_reposition_cost'] * capacity * cost_increase / 1000
+          yearly_capital_cost += -npf.pmt(self.interest_rate, self.credit_years,
+                                          investment_cost * self.percentage_credit)
+      for storage_system in self.individual['Energy Storage Components']:
+        if (i % storage_system['lifetime']) == 0 and i != (life_cycle_duration - 1):
+          cost_increase = math.pow(1 + self.consumer_price_index, i)
+          capacity = storage_system['volume'] if cte.THERMAL in storage_system['type'] else storage_system['capacity']
+          yearly_capital_cost += storage_system['unit_reposition_cost'] * capacity * cost_increase
+          yearly_capital_cost += -npf.pmt(self.interest_rate, self.credit_years,
+                                          investment_cost * self.percentage_credit)
+      capital_costs_cash_flow.append(yearly_capital_cost)
+      operational_costs_cash_flow.append(yearly_operational_cost)
+      maintenance_costs_cash_flow.append(yearly_maintenance_cost)
+    for year in range(1, life_cycle_duration + 1):
+      price_increase = math.pow(1 + self.consumer_price_index, year)
+      if year == life_cycle_duration:
+        end_of_life_costs[year] = (
+            self.building.thermal_zones_from_internal_zones[0].total_floor_area *
+            self.individual['End of Life Cost'] * price_increase
+        )
+
+    life_cycle_capital_cost = npf.npv(self.discount_rate, capital_costs_cash_flow)
+    life_cycle_operational_cost = npf.npv(self.discount_rate, operational_costs_cash_flow)
+    life_cycle_maintenance_cost = npf.npv(self.discount_rate, maintenance_costs_cash_flow)
+    life_cycle_end_of_life_cost = npf.npv(self.discount_rate, end_of_life_costs)
+    total_life_cycle_cost = life_cycle_capital_cost + life_cycle_operational_cost + life_cycle_maintenance_cost + life_cycle_end_of_life_cost
+    return total_life_cycle_cost
 
   def costs_archetype(self):
     costs_catalogue = CostsCatalogFactory('montreal_new').catalog.entries().archetypes
@@ -255,7 +379,7 @@ class Individual:
         if generation_system.energy_storage_systems is not None:
           for energy_storage_system in generation_system.energy_storage_systems:
             if energy_storage_system.type_energy_stored == cte.THERMAL:
-              if energy_storage_system.heating_coil_capacity is None:
+              if energy_storage_system.heating_coil_capacity is not None:
                 investment_cost += float(capital_costs_chapter.item('D306010_storage_tank').initial_investment[0])
                 reposition_cost += float(capital_costs_chapter.item('D306010_storage_tank').reposition[0])
                 lifetime += float(capital_costs_chapter.item('D306010_storage_tank').lifetime)
@@ -294,10 +418,19 @@ class Individual:
         maintenance_cost += item.maintenance[0]
     return maintenance_cost
 
-  def fuel_cost(self, fuel_type, fuel_cost_tariff_type):
+  def fuel_cost_per_kwh(self, fuel_type, fuel_cost_tariff_type):
+    fuel_cost = 0
     catalog_fuels = self.costs.operational_cost.fuels
     for fuel in catalog_fuels:
       if fuel_type == fuel.type and fuel_cost_tariff_type == fuel.variable.rate_type:
-        return fuel
-
-
+        if fuel.type == cte.ELECTRICITY and fuel_cost_tariff_type == 'fixed':
+          fuel_cost = fuel.variable.values[0]
+        elif fuel.type == cte.ELECTRICITY and fuel_cost_tariff_type == 'variable':
+          fuel_cost = fuel.variable.values[0]
+        else:
+          if fuel.type == cte.BIOMASS:
+            conversion_factor = 1
+          else:
+            conversion_factor = fuel.density[0]
+          fuel_cost = fuel.variable.values[0] / (conversion_factor * fuel.lower_heating_value[0] * 0.277)
+    return fuel_cost

hub/data/energy_systems/montreal_future_systems.xml

@@ -1274,7 +1274,7 @@
 			<storage_type>sensible</storage_type>
 			<nominal_capacity/>
 			<losses_ratio/>
-            <heating_coil_capacity>5000</heating_coil_capacity>
+            <heating_coil_capacity>0</heating_coil_capacity>
         </templateStorages>
   </energy_storage_components>
   <materials>

main.py

@@ -3,8 +3,7 @@ import subprocess
 from building_modelling.ep_run_enrich import energy_plus_workflow
 from energy_system_modelling_package.energy_system_modelling_factories.montreal_energy_system_archetype_modelling_factory import \
   MontrealEnergySystemArchetypesSimulationFactory
-from energy_system_modelling_package.energy_system_modelling_factories.system_sizing_methods.genetic_algorithm.individual import \
-  Individual
+from energy_system_modelling_package.energy_system_modelling_factories.system_sizing_methods.genetic_algorithm.genetic_algorithm import GeneticAlgorithm
 from hub.imports.geometry_factory import GeometryFactory
 from hub.helpers.dictionaries import Dictionaries
 from hub.imports.construction_factory import ConstructionFactory
@@ -58,9 +57,7 @@ for building in city.buildings:
   heating_design_load = building.heating_peak_load[cte.YEAR][0]
   cooling_design_load = building.cooling_peak_load[cte.YEAR][0]
   available_space = building.volume / building.storeys_above_ground
-  Individual(building, building.energy_systems[2], energy_demands=None, optimization_scenario='cost',
-             heating_design_load=heating_design_load, cooling_design_load=cooling_design_load,
-             available_space=available_space, dt=3600, electricity_tariff_type='fixed').system_components()
+  GeneticAlgorithm(optimization_scenario='cost').solve_ga(building=building, energy_system=building.energy_systems[1])
 # EnergySystemsSizingFactory('pv_sizing', city).enrich()
 # EnergySystemsSizingFactory('peak_load_sizing', city).enrich()
 # EnergySystemsSizingFactory('heuristic_sizing', city).enrich()

test.py

@@ -0,0 +1,84 @@
+import csv
+from pathlib import Path
+from building_modelling.ep_run_enrich import energy_plus_workflow
+from energy_system_modelling_package.energy_system_modelling_factories.energy_system_sizing_factory import \
+  EnergySystemsSizingFactory
+from energy_system_modelling_package.energy_system_modelling_factories.hvac_dhw_systems_simulation_models.heat_pump_boiler_tes_heating import \
+  HeatPumpBoilerTesHeating
+from hub.imports.geometry_factory import GeometryFactory
+from hub.helpers.dictionaries import Dictionaries
+from hub.imports.construction_factory import ConstructionFactory
+from hub.imports.usage_factory import UsageFactory
+from hub.imports.weather_factory import WeatherFactory
+import hub.helpers.constants as cte
+from building_modelling.geojson_creator import process_geojson
+from energy_system_modelling_package import random_assignation
+from hub.imports.energy_systems_factory import EnergySystemsFactory
+
+# Specify the GeoJSON file path
+input_files_path = (Path(__file__).parent / 'input_files')
+input_files_path.mkdir(parents=True, exist_ok=True)
+geojson_file = process_geojson(x=-73.5681295982132, y=45.49218262677643, diff=0.00006)
+geojson_file_path = input_files_path / 'output_buildings.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)
+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()
+WeatherFactory('epw', city).enrich()
+energy_plus_workflow(city, energy_plus_output_path)
+random_assignation.call_random(city.buildings, random_assignation.residential_new_systems_percentage)
+EnergySystemsFactory('montreal_future', city).enrich()
+EnergySystemsSizingFactory('peak_load_sizing', city).enrich()
+hp = city.buildings[0].energy_systems[1].generation_systems[1]
+boiler = city.buildings[0].energy_systems[1].generation_systems[0]
+tes = city.buildings[0].energy_systems[1].generation_systems[0].energy_storage_systems[0]
+# Step 1: Calculate daily demands
+daily_demands = [sum(city.buildings[0].heating_demand[cte.HOUR][i:i + 24]) for i in
+                 range(0, len(city.buildings[0].heating_demand[cte.HOUR]), 24)]
+# Step 2: Find the day with maximum demand
+max_day_index = daily_demands.index(max(daily_demands))
+# Step 3: Extract the hourly demands for the day before, the day with max demand, and the day after
+# Ensure you don't go out of bounds if max_day_index is 0 or the last day
+if max_day_index > 0 and max_day_index < len(daily_demands) - 1:
+  start_index = (max_day_index - 1) * 24
+  end_index = (max_day_index + 2) * 24
+  three_day_demands = city.buildings[0].heating_demand[cte.HOUR][start_index:end_index]
+elif max_day_index == 0:  # If max is the first day, just include day 1 and day 2
+  start_index = max_day_index * 24
+  end_index = (max_day_index + 2) * 24
+  three_day_demands = city.buildings[0].heating_demand[cte.HOUR][start_index:end_index]
+else:  # If max is the last day, just include the last two days
+  start_index = (max_day_index - 1) * 24
+  end_index = len(city.buildings[0].heating_demand[cte.HOUR])
+  three_day_demands = city.buildings[0].heating_demand[cte.HOUR][start_index:end_index]
+results = HeatPumpBoilerTesHeating(hp=hp,
+                                   boiler=boiler,
+                                   tes=tes,
+                                   hourly_heating_demand_joules=three_day_demands,
+                                   heating_peak_load_watts=max(three_day_demands),
+                                   upper_limit_tes=55,
+                                   outdoor_temperature=city.buildings[0].external_temperature[cte.HOUR],
+                                   dt=900).simulation()
+heating_eui = city.buildings[0].heating_demand[cte.YEAR][0] / (
+    city.buildings[0].thermal_zones_from_internal_zones[0].total_floor_area * 3.6e6)
+file_name = f'energy_system_simulation_results_{city.buildings[0].name}.csv'
+with open(simulation_results_path / file_name, 'w', newline='') as csvfile:
+  output_file = csv.writer(csvfile)
+  # Write header
+  output_file.writerow(results.keys())
+  # Write data
+  output_file.writerows(zip(*results.values()))