diff --git a/hub/city_model_structure/building.py b/hub/city_model_structure/building.py index b928535b..660f0f7e 100644 --- a/hub/city_model_structure/building.py +++ b/hub/city_model_structure/building.py @@ -896,13 +896,13 @@ class Building(CityObject): for energy_system in energy_systems: if cte.COOLING in energy_system.demand_types and cte.COOLING not in fuel_breakdown[key]: for generation_system in energy_system.generation_systems: - fuel_breakdown[generation_system.fuel_type][cte.COOLING] = self.cooling_consumption[cte.YEAR][0] / 3600 + fuel_breakdown[generation_system.fuel_type][cte.COOLING] = self.cooling_consumption[cte.YEAR][0] for fuel in heating_fuels: if cte.HEATING not in fuel_breakdown[fuel]: for energy_system in energy_systems: if cte.HEATING in energy_system.demand_types: for generation_system in energy_system.generation_systems: - fuel_breakdown[generation_system.fuel_type][cte.HEATING] = self.heating_consumption[cte.YEAR][0] / 3600 + fuel_breakdown[generation_system.fuel_type][cte.HEATING] = self.heating_consumption[cte.YEAR][0] for fuel in dhw_fuels: if cte.DOMESTIC_HOT_WATER not in fuel_breakdown[fuel]: for energy_system in energy_systems: diff --git a/hub/helpers/constants.py b/hub/helpers/constants.py index af22e94d..ad32c835 100644 --- a/hub/helpers/constants.py +++ b/hub/helpers/constants.py @@ -10,6 +10,7 @@ Project Coder Pilar Monsalvete Alvarez de Uribarri pilar.monsalvete@concordia.ca KELVIN = 273.15 WATER_DENSITY = 1000 # kg/m3 WATER_HEAT_CAPACITY = 4182 # J/kgK +WATER_THERMAL_CONDUCTIVITY = 0.65 # W/mK NATURAL_GAS_LHV = 36.6e6 # J/m3 AIR_DENSITY = 1.293 # kg/m3 AIR_HEAT_CAPACITY = 1005.2 # J/kgK diff --git a/hub/imports/results/ep_multiple_buildings.py b/hub/imports/results/ep_multiple_buildings.py index 9da12c7a..bd9a2dd3 100644 --- a/hub/imports/results/ep_multiple_buildings.py +++ b/hub/imports/results/ep_multiple_buildings.py @@ -60,9 +60,12 @@ class EnergyPlusMultipleBuildings: for building in self._city.buildings: building.heating_demand[cte.HOUR] = building_energy_demands[f'Building {building.name} Heating Demand (J)'] building.cooling_demand[cte.HOUR] = building_energy_demands[f'Building {building.name} Cooling Demand (J)'] - building.domestic_hot_water_heat_demand[cte.HOUR] = building_energy_demands[f'Building {building.name} DHW Demand (W)'] - building.appliances_electrical_demand[cte.HOUR] = building_energy_demands[f'Building {building.name} Appliances (W)'] - building.lighting_electrical_demand[cte.HOUR] = building_energy_demands[f'Building {building.name} Lighting (W)'] + building.domestic_hot_water_heat_demand[cte.HOUR] = \ + [x * cte.WATTS_HOUR_TO_JULES for x in building_energy_demands[f'Building {building.name} DHW Demand (W)']] + building.appliances_electrical_demand[cte.HOUR] = \ + [x * cte.WATTS_HOUR_TO_JULES for x in building_energy_demands[f'Building {building.name} Appliances (W)']] + building.lighting_electrical_demand[cte.HOUR] = \ + [x * cte.WATTS_HOUR_TO_JULES for x in building_energy_demands[f'Building {building.name} Lighting (W)']] building.heating_demand[cte.MONTH] = MonthlyValues.get_total_month(building.heating_demand[cte.HOUR]) building.cooling_demand[cte.MONTH] = MonthlyValues.get_total_month(building.cooling_demand[cte.HOUR]) building.domestic_hot_water_heat_demand[cte.MONTH] = ( diff --git a/hub/imports/weather/epw_weather_parameters.py b/hub/imports/weather/epw_weather_parameters.py index d8f3a002..0c362261 100644 --- a/hub/imports/weather/epw_weather_parameters.py +++ b/hub/imports/weather/epw_weather_parameters.py @@ -126,10 +126,10 @@ class EpwWeatherParameters: for building in self._city.buildings: building.external_temperature[cte.MONTH] = \ MonthlyValues().get_mean_values(building.external_temperature[cte.HOUR]) - building.external_temperature[cte.YEAR] = [sum(building.external_temperature[cte.HOUR]) / 9870] + building.external_temperature[cte.YEAR] = [sum(building.external_temperature[cte.HOUR]) / 8760] building.cold_water_temperature[cte.MONTH] = \ MonthlyValues().get_mean_values(building.cold_water_temperature[cte.HOUR]) - building.cold_water_temperature[cte.YEAR] = [sum(building.cold_water_temperature[cte.HOUR]) / 9870] + building.cold_water_temperature[cte.YEAR] = [sum(building.cold_water_temperature[cte.HOUR]) / 8760] # If the usage has already being imported, the domestic hot water missing values must be calculated here that # the cold water temperature is finally known diff --git a/hub/imports/weather/helpers/weather.py b/hub/imports/weather/helpers/weather.py index 7603cb5b..755f9ad3 100644 --- a/hub/imports/weather/helpers/weather.py +++ b/hub/imports/weather/helpers/weather.py @@ -8,7 +8,7 @@ Project Coder Pilar Monsalvete Alvarez de Uribarri pilar.monsalvete@concordia.ca import logging import math import hub.helpers.constants as cte - +from datetime import datetime, timedelta class Weather: """ @@ -55,25 +55,19 @@ class Weather: # and Craig Christensen, National Renewable Energy Laboratory # ambient temperatures( in °C) # cold water temperatures( in °C) - ambient_temperature_fahrenheit = [] - average_temperature = 0 - maximum_temperature = -1000 - minimum_temperature = 1000 - for temperature in ambient_temperature: - value = temperature * 9 / 5 + 32 - ambient_temperature_fahrenheit.append(value) - average_temperature += value / 8760 - if value > maximum_temperature: - maximum_temperature = value - if value < minimum_temperature: - minimum_temperature = value - delta_temperature = maximum_temperature - minimum_temperature - ratio = 0.4 + 0.01 * (average_temperature - 44) - lag = 35 - 1 * (average_temperature - 44) + t_out_fahrenheit = [1.8 * t_out + 32 for t_out in ambient_temperature] + t_out_average = sum(t_out_fahrenheit) / len(t_out_fahrenheit) + max_difference = max(t_out_fahrenheit) - min(t_out_fahrenheit) + ratio = 0.4 + 0.01 * (t_out_average - 44) + lag = 35 - (t_out_average - 35) + number_of_day = [a for a in range(1, 366)] + day_of_year = [day for day in number_of_day for _ in range(24)] + cold_temperature_fahrenheit = [] cold_temperature = [] - for temperature in ambient_temperature_fahrenheit: - radians = (0.986 * (temperature-15-lag) - 90) * math.pi / 180 - cold_temperature.append((average_temperature + 6 + ratio * (delta_temperature/2) * math.sin(radians) - 32) * 5/9) + for i in range(len(ambient_temperature)): + cold_temperature_fahrenheit.append(t_out_average + 6 + ratio * (max_difference / 2) * + math.sin(math.radians(0.986 * (day_of_year[i] - 15 - lag) - 90))) + cold_temperature.append((cold_temperature_fahrenheit[i] - 32) / 1.8) return cold_temperature def epw_file(self, region_code): diff --git a/main.py b/main.py index e69de29b..798b5d4c 100644 --- a/main.py +++ b/main.py @@ -0,0 +1,45 @@ +from scripts.geojson_creator import process_geojson +from pathlib import Path +import subprocess +from scripts.ep_run_enrich import energy_plus_workflow +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 +from hub.imports.results_factory import ResultFactory +from scripts.energy_system_analysis_report import EnergySystemAnalysisReport +from scripts import random_assignation +from hub.imports.energy_systems_factory import EnergySystemsFactory +from scripts.energy_system_sizing import SystemSizing +from scripts.energy_system_retrofit_results import system_results, new_system_results +from scripts.energy_system_sizing_and_simulation_factory import EnergySystemsSimulationFactory +from scripts.costs.cost import Cost +from scripts.costs.constants import SKIN_RETROFIT_AND_SYSTEM_RETROFIT_AND_PV, SYSTEM_RETROFIT_AND_PV +import hub.helpers.constants as cte +from hub.exports.exports_factory import ExportsFactory +# Specify the GeoJSON file path +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() +# Create city object from GeoJSON file +city = GeometryFactory('geojson', + path=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 +# Enrich city data +ConstructionFactory('nrcan', city).enrich() + +UsageFactory('nrcan', city).enrich() +WeatherFactory('epw', city).enrich() +energy_plus_workflow(city) +random_assignation.call_random(city.buildings, random_assignation.residential_new_systems_percentage) +EnergySystemsFactory('montreal_future', city).enrich() + +for building in city.buildings: + EnergySystemsSimulationFactory('archetype13', building=building, output_path=output_path).enrich() + +print('test') diff --git a/scripts/energy_system_sizing_and_simulation_factory.py b/scripts/energy_system_sizing_and_simulation_factory.py index 95a4d4e8..b0526e2c 100644 --- a/scripts/energy_system_sizing_and_simulation_factory.py +++ b/scripts/energy_system_sizing_and_simulation_factory.py @@ -6,6 +6,7 @@ Project Coder Saeed Ranjbar saeed.ranjbar@mail.concordia.ca """ from scripts.system_simulation_models.archetype13 import Archetype13 +from scripts.system_simulation_models.archetype13_stratified_tes import Archetype13Stratified from scripts.system_simulation_models.archetype1 import Archetype1 diff --git a/scripts/system_simulation_models/archetype13.py b/scripts/system_simulation_models/archetype13.py index 13904273..786115db 100644 --- a/scripts/system_simulation_models/archetype13.py +++ b/scripts/system_simulation_models/archetype13.py @@ -2,6 +2,8 @@ import math import hub.helpers.constants as cte import csv from hub.helpers.monthly_values import MonthlyValues + + class Archetype13: def __init__(self, building, output_path): self._building = building @@ -15,12 +17,13 @@ class Archetype13: self._heating_peak_load = building.heating_peak_load[cte.YEAR][0] self._cooling_peak_load = building.cooling_peak_load[cte.YEAR][0] self._domestic_hot_water_peak_load = building.domestic_hot_water_peak_load[cte.YEAR][0] - self._hourly_heating_demand = [0] + [demand / 3600 for demand in building.heating_demand[cte.HOUR]] - self._hourly_cooling_demand = [demand / 3600 for demand in building.cooling_demand[cte.HOUR]] - self._hourly_dhw_demand = [0] + building.domestic_hot_water_heat_demand[cte.HOUR] + self._hourly_heating_demand = [demand / cte.HOUR_TO_SECONDS for demand in building.heating_demand[cte.HOUR]] + self._hourly_cooling_demand = [demand / cte.HOUR_TO_SECONDS for demand in building.cooling_demand[cte.HOUR]] + self._hourly_dhw_demand = building.domestic_hot_water_heat_demand[cte.HOUR] self._output_path = output_path - self._t_out = [0] + building.external_temperature[cte.HOUR] + self._t_out = building.external_temperature[cte.HOUR] self.results = {} + self.dt = 900 def hvac_sizing(self): storage_factor = 3 @@ -47,7 +50,9 @@ class Archetype13: def heating_system_simulation(self): hp, boiler, tes = self.hvac_sizing() cop_curve_coefficients = [float(coefficient) for coefficient in hp.heat_efficiency_curve.coefficients] - demand = self._hourly_heating_demand + number_of_ts = int(cte.HOUR_TO_SECONDS / self.dt) + demand = [0] + [x for x in self._hourly_heating_demand for _ in range(number_of_ts)] + t_out = [0] + [x for x in self._t_out for _ in range(number_of_ts)] hp.source_temperature = self._t_out variable_names = ["t_sup_hp", "t_tank", "t_ret", "m_ch", "m_dis", "q_hp", "q_boiler", "hp_cop", "hp_electricity", "boiler_gas_consumption", "t_sup_boiler", "boiler_energy_consumption", @@ -56,9 +61,8 @@ class Archetype13: variables = {name: [0] * num_hours for name in variable_names} (t_sup_hp, t_tank, t_ret, m_ch, m_dis, q_hp, q_boiler, hp_cop, hp_electricity, boiler_gas_consumption, t_sup_boiler, boiler_energy_consumption, heating_consumption) = \ - [variables[name] for name in variable_names] + [variables[name] for name in variable_names] t_tank[0] = 55 - dt = 3600 hp_heating_cap = hp.nominal_heat_output boiler_heating_cap = boiler.nominal_heat_output hp_delta_t = 5 @@ -75,8 +79,8 @@ class Archetype13: for i in range(len(demand) - 1): t_tank[i + 1] = (t_tank[i] + (m_ch[i] * (t_sup_boiler[i] - t_tank[i]) + - (ua * (self._t_out[i] - t_tank[i])) / cte.WATER_HEAT_CAPACITY - - m_dis[i] * (t_tank[i] - t_ret[i])) * (dt / (cte.WATER_DENSITY * v))) + (ua * (t_out[i] - t_tank[i])) / cte.WATER_HEAT_CAPACITY - + m_dis[i] * (t_tank[i] - t_ret[i])) * (self.dt / (cte.WATER_DENSITY * v))) # hp operation if t_tank[i + 1] < 40: q_hp[i + 1] = hp_heating_cap @@ -93,7 +97,7 @@ class Archetype13: else: q_hp[i + 1], m_ch[i + 1], t_sup_hp[i + 1] = 0, 0, t_tank[i + 1] t_sup_hp_fahrenheit = 1.8 * t_sup_hp[i + 1] + 32 - t_out_fahrenheit = 1.8 * self._t_out[i + 1] + 32 + t_out_fahrenheit = 1.8 * t_out[i + 1] + 32 if q_hp[i + 1] > 0: hp_cop[i + 1] = (cop_curve_coefficients[0] + cop_curve_coefficients[1] * t_sup_hp_fahrenheit + @@ -109,39 +113,54 @@ class Archetype13: if q_hp[i + 1] > 0: if t_sup_hp[i + 1] < 45: q_boiler[i + 1] = boiler_heating_cap - elif demand[i + 1] > 0.5 * self._heating_peak_load / dt: + elif demand[i + 1] > 0.5 * self._heating_peak_load / self.dt: q_boiler[i + 1] = 0.5 * boiler_heating_cap boiler_energy_consumption[i + 1] = q_boiler[i + 1] / boiler_efficiency - boiler_gas_consumption[i + 1] = (q_boiler[i + 1] * dt) / (boiler_efficiency * cte.NATURAL_GAS_LHV) + boiler_gas_consumption[i + 1] = (q_boiler[i + 1] * self.dt) / (boiler_efficiency * cte.NATURAL_GAS_LHV) t_sup_boiler[i + 1] = t_sup_hp[i + 1] + (q_boiler[i + 1] / (m_ch[i + 1] * cte.WATER_HEAT_CAPACITY)) # storage discharging if demand[i + 1] == 0: m_dis[i + 1] = 0 t_ret[i + 1] = t_tank[i + 1] else: - if demand[i + 1] > 0.5 * self._heating_peak_load / dt: + if demand[i + 1] > 0.5 * self._heating_peak_load / cte.HOUR_TO_SECONDS: factor = 8 else: factor = 4 - m_dis[i + 1] = self._heating_peak_load / (cte.WATER_HEAT_CAPACITY * factor * dt) + m_dis[i + 1] = self._heating_peak_load / (cte.WATER_HEAT_CAPACITY * factor * cte.HOUR_TO_SECONDS) t_ret[i + 1] = t_tank[i + 1] - demand[i + 1] / (m_dis[i + 1] * cte.WATER_HEAT_CAPACITY) - + tes.temperature = [] + hp_electricity_j = [(x * cte.WATTS_HOUR_TO_JULES) / number_of_ts for x in hp_electricity] + boiler_consumption_j = [(x * cte.WATTS_HOUR_TO_JULES) / number_of_ts for x in boiler_energy_consumption] + hp_hourly = [] + boiler_hourly = [] + boiler_sum = 0 + hp_sum = 0 + for i in range(1, len(demand)): + hp_sum += hp_electricity_j[i] + boiler_sum += boiler_consumption_j[i] + if (i - 1) % number_of_ts == 0: + tes.temperature.append(t_tank[i]) + hp_hourly.append(hp_sum) + boiler_hourly.append(boiler_sum) + hp_sum = 0 + boiler_sum = 0 hp.energy_consumption[cte.HEATING] = {} - hp.energy_consumption[cte.HEATING][cte.HOUR] = hp_electricity + 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_energy_consumption + 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])] - tes.temperature = t_tank + self.results['Heating Demand (W)'] = demand self.results['HP Heat Output (W)'] = q_hp - self.results['HP Source Temperature'] = self._t_out + self.results['HP Source Temperature'] = t_out self.results['HP Supply Temperature'] = t_sup_hp self.results['HP COP'] = hp_cop self.results['HP Electricity Consumption (W)'] = hp_electricity @@ -152,42 +171,47 @@ class Archetype13: self.results['TES Charging Flow Rate (kg/s)'] = m_ch self.results['TES Discharge Flow Rate (kg/s)'] = m_dis self.results['Heating Loop Return Temperature'] = t_ret - return hp_electricity, boiler_energy_consumption + return hp_hourly, boiler_hourly def cooling_system_simulation(self): hp = self.hvac_sizing()[0] eer_curve_coefficients = [float(coefficient) for coefficient in hp.cooling_efficiency_curve.coefficients] cooling_efficiency = float(hp.cooling_efficiency) - demand = self._hourly_cooling_demand + number_of_ts = int(cte.HOUR_TO_SECONDS / self.dt) + demand = [0] + [x for x in self._hourly_cooling_demand for _ in range(number_of_ts)] + t_out = [0] + [x for x in self._t_out for _ in range(number_of_ts)] hp.source_temperature = self._t_out variable_names = ["t_sup_hp", "t_ret", "m", "q_hp", "hp_electricity", "hp_eer"] num_hours = len(demand) variables = {name: [0] * num_hours for name in variable_names} (t_sup_hp, t_ret, m, q_hp, hp_electricity, hp_eer) = [variables[name] for name in variable_names] t_ret[0] = 13 - dt = 3600 - for i in range(len(demand) - 1): + + for i in range(1, len(demand)): if demand[i] > 0: - m[i] = self._cooling_peak_load / (cte.WATER_HEAT_CAPACITY * 5 * dt) - if t_ret[i] > 13: - if demand[i] < 0.25 * self._cooling_peak_load / dt: + m[i] = self._cooling_peak_load / (cte.WATER_HEAT_CAPACITY * 5 * cte.HOUR_TO_SECONDS) + if t_ret[i - 1] >= 13: + if demand[i] < 0.25 * self._cooling_peak_load / cte.HOUR_TO_SECONDS: q_hp[i] = 0.25 * hp.nominal_cooling_output - elif demand[i] < 0.5 * self._cooling_peak_load / dt: + elif demand[i] < 0.5 * self._cooling_peak_load / cte.HOUR_TO_SECONDS: q_hp[i] = 0.5 * hp.nominal_cooling_output else: q_hp[i] = hp.nominal_cooling_output - t_sup_hp[i] = t_ret[i] - q_hp[i] / (m[i] * cte.WATER_HEAT_CAPACITY) + t_sup_hp[i] = t_ret[i - 1] - q_hp[i] / (m[i] * cte.WATER_HEAT_CAPACITY) else: q_hp[i] = 0 - t_sup_hp[i] = t_ret[i] - t_ret[i + 1] = t_sup_hp[i] + demand[i] / (m[i] * cte.WATER_HEAT_CAPACITY) + t_sup_hp[i] = t_ret[i - 1] + if m[i] == 0: + t_ret[i] = t_sup_hp[i] + else: + t_ret[i] = t_sup_hp[i] + demand[i] / (m[i] * cte.WATER_HEAT_CAPACITY) else: m[i] = 0 q_hp[i] = 0 - t_sup_hp[i] = t_ret[i] - t_ret[i + 1] = t_ret[i] + t_sup_hp[i] = t_ret[i -1] + t_ret[i] = t_ret[i - 1] t_sup_hp_fahrenheit = 1.8 * t_sup_hp[i] + 32 - t_out_fahrenheit = 1.8 * self._t_out[i] + 32 + t_out_fahrenheit = 1.8 * t_out[i] + 32 if q_hp[i] > 0: hp_eer[i] = (eer_curve_coefficients[0] + eer_curve_coefficients[1] * t_sup_hp_fahrenheit + @@ -195,12 +219,20 @@ class Archetype13: eer_curve_coefficients[3] * t_out_fahrenheit + eer_curve_coefficients[4] * t_out_fahrenheit ** 2 + eer_curve_coefficients[5] * t_sup_hp_fahrenheit * t_out_fahrenheit) - hp_electricity[i] = q_hp[i] / cooling_efficiency + hp_electricity[i] = q_hp[i] / hp_eer[i] else: hp_eer[i] = 0 hp_electricity[i] = 0 + hp_electricity_j = [(x * cte.WATTS_HOUR_TO_JULES) / number_of_ts for x in hp_electricity] + hp_hourly = [] + hp_sum = 0 + for i in range(1, len(demand)): + hp_sum += hp_electricity_j[i] + if (i - 1) % number_of_ts == 0: + hp_hourly.append(hp_sum) + hp_sum = 0 hp.energy_consumption[cte.COOLING] = {} - hp.energy_consumption[cte.COOLING][cte.HOUR] = hp_electricity + hp.energy_consumption[cte.COOLING][cte.HOUR] = hp_hourly hp.energy_consumption[cte.COOLING][cte.MONTH] = MonthlyValues.get_total_month( hp.energy_consumption[cte.COOLING][cte.HOUR]) hp.energy_consumption[cte.COOLING][cte.YEAR] = [ @@ -212,12 +244,14 @@ class Archetype13: self.results['HP Electricity Consumption'] = hp_electricity self.results['Cooling Loop Flow Rate (kg/s)'] = m self.results['Cooling Loop Return Temperature'] = t_ret - return hp_electricity + return hp_hourly def dhw_system_simulation(self): hp, tes = self.dhw_sizing() cop_curve_coefficients = [float(coefficient) for coefficient in hp.heat_efficiency_curve.coefficients] - demand = self._hourly_dhw_demand + number_of_ts = int(cte.HOUR_TO_SECONDS / self.dt) + demand = [0] + [x for x in self._hourly_dhw_demand for _ in range(number_of_ts)] + t_out = [0] + [x for x in self._t_out 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)"] num_hours = len(demand) @@ -226,7 +260,7 @@ class Archetype13: [variables[name] for name in variable_names] t_tank[0] = 70 v_dhw[0] = tes.volume - dt = 3600 + hp_heating_cap = hp.nominal_heat_output hp_delta_t = 8 v, h = float(tes.volume), float(tes.height) @@ -239,26 +273,24 @@ class Archetype13: ua = u_tot * (2 * a_top + a_side) freshwater_temperature = 18 for i in range(len(demand) - 1): - delta_t_demand = demand[i] * (dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) + delta_t_demand = demand[i] * (self.dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) if t_tank[i] < 65: q_hp[i] = hp_heating_cap - delta_t_hp = q_hp[i] * (dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) + delta_t_hp = q_hp[i] * (self.dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) if demand[i] > 0: dhw_needed = (demand[i] * cte.HOUR_TO_SECONDS) / (cte.WATER_HEAT_CAPACITY * t_tank[i] * cte.WATER_DENSITY) m_dis[i] = dhw_needed * cte.WATER_DENSITY / cte.HOUR_TO_SECONDS m_refill[i] = m_dis[i] - delta_t_freshwater = m_refill[i] * (t_tank[i] - freshwater_temperature) * (dt / (v * cte.WATER_DENSITY)) + delta_t_freshwater = m_refill[i] * (t_tank[i] - freshwater_temperature) * (self.dt / (v * cte.WATER_DENSITY)) diff = delta_t_freshwater + delta_t_demand - delta_t_hp if diff > 0: if diff > 0: - power = diff * (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v) / dt + power = diff * (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v) / self.dt if power <= float(tes.heating_coil_capacity): q_coil[i] = power else: q_coil[i] = float(tes.heating_coil_capacity) - elif t_tank[i] < 65: - q_coil[i] = float(tes.heating_coil_capacity) - delta_t_coil = q_coil[i] * (dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) + delta_t_coil = q_coil[i] * (self.dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) if q_hp[i] > 0: m_ch[i] = q_hp[i] / (cte.WATER_HEAT_CAPACITY * hp_delta_t) @@ -267,7 +299,7 @@ class Archetype13: m_ch[i] = 0 t_sup_hp[i] = t_tank[i] t_sup_hp_fahrenheit = 1.8 * t_sup_hp[i] + 32 - t_out_fahrenheit = 1.8 * self._t_out[i] + 32 + t_out_fahrenheit = 1.8 * t_out[i] + 32 if q_hp[i] > 0: hp_cop[i] = (cop_curve_coefficients[0] + cop_curve_coefficients[1] * t_sup_hp_fahrenheit + @@ -275,32 +307,47 @@ class Archetype13: cop_curve_coefficients[3] * t_out_fahrenheit + cop_curve_coefficients[4] * t_out_fahrenheit ** 2 + cop_curve_coefficients[5] * t_sup_hp_fahrenheit * t_out_fahrenheit) - hp_electricity[i] = q_hp[i] / 3.5 + hp_electricity[i] = q_hp[i] / hp_cop[i] else: hp_cop[i] = 0 hp_electricity[i] = 0 t_tank[i + 1] = t_tank[i] + (delta_t_hp - delta_t_freshwater - delta_t_demand + delta_t_coil) + tes.temperature = [] + hp_electricity_j = [(x * cte.WATTS_HOUR_TO_JULES) / number_of_ts for x in hp_electricity] + heating_coil_j = [(x * cte.WATTS_HOUR_TO_JULES) / number_of_ts for x in q_coil] + hp_hourly = [] + coil_hourly = [] + coil_sum = 0 + hp_sum = 0 + for i in range(1, len(demand)): + hp_sum += hp_electricity_j[i] + coil_sum += heating_coil_j[i] + if (i - 1) % number_of_ts == 0: + tes.temperature.append(t_tank[i]) + hp_hourly.append(hp_sum) + coil_hourly.append(coil_sum) + hp_sum = 0 + coil_sum = 0 hp.energy_consumption[cte.DOMESTIC_HOT_WATER] = {} - hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.HOUR] = hp_electricity + hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.HOUR] = hp_hourly hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.MONTH] = MonthlyValues.get_total_month( hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.HOUR]) hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.YEAR] = [ sum(hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.MONTH])] tes.heating_coil_energy_consumption = {} - tes.heating_coil_energy_consumption[cte.HOUR] = q_coil + tes.heating_coil_energy_consumption[cte.HOUR] = coil_hourly tes.heating_coil_energy_consumption[cte.MONTH] = MonthlyValues.get_total_month( tes.heating_coil_energy_consumption[cte.HOUR]) tes.heating_coil_energy_consumption[cte.YEAR] = [ sum(tes.heating_coil_energy_consumption[cte.MONTH])] tes.temperature = t_tank - self.results['DHW Demand (W)'] = demand self.results['DHW HP Heat Output (W)'] = q_hp self.results['DHW HP Electricity Consumption (W)'] = hp_electricity - self.results['DHW HP Source Temperature'] = self._t_out + self.results['DHW HP Source Temperature'] = t_out self.results['DHW HP Supply Temperature'] = t_sup_hp self.results['DHW HP COP'] = hp_cop self.results['DHW TES Heating Coil Heat Output (W)'] = q_coil @@ -309,7 +356,7 @@ class Archetype13: self.results['DHW Flow Rate (kg/s)'] = m_dis self.results['DHW TES Refill Flow Rate (kg/s)'] = m_refill self.results['Available Water in Tank (m3)'] = v_dhw - return hp_electricity, q_coil + return hp_hourly, coil_hourly def enrich_buildings(self): hp_heating, boiler_consumption = self.heating_system_simulation() @@ -320,11 +367,11 @@ class Archetype13: self._building.heating_consumption[cte.HOUR] = heating_consumption self._building.heating_consumption[cte.MONTH] = ( MonthlyValues.get_total_month(self._building.heating_consumption[cte.HOUR])) - self._building.heating_consumption[cte.YEAR] = sum(self._building.heating_consumption[cte.MONTH]) + self._building.heating_consumption[cte.YEAR] = [sum(self._building.heating_consumption[cte.MONTH])] self._building.cooling_consumption[cte.HOUR] = hp_cooling self._building.cooling_consumption[cte.MONTH] = ( MonthlyValues.get_total_month(self._building.cooling_consumption[cte.HOUR])) - self._building.cooling_consumption[cte.YEAR] = sum(self._building.cooling_consumption[cte.MONTH]) + self._building.cooling_consumption[cte.YEAR] = [sum(self._building.cooling_consumption[cte.MONTH])] self._building.domestic_hot_water_consumption[cte.HOUR] = dhw_consumption self._building.domestic_hot_water_consumption[cte.MONTH] = ( MonthlyValues.get_total_month(self._building.domestic_hot_water_consumption[cte.HOUR])) diff --git a/scripts/system_simulation_models/archetype13_stratified_tes.py b/scripts/system_simulation_models/archetype13_stratified_tes.py new file mode 100644 index 00000000..c5bdd1e7 --- /dev/null +++ b/scripts/system_simulation_models/archetype13_stratified_tes.py @@ -0,0 +1,392 @@ +import math +import hub.helpers.constants as cte +import csv +from hub.helpers.monthly_values import MonthlyValues +import numpy as np + + +class Archetype13Stratified: + def __init__(self, building, output_path): + self._building = building + self._name = building.name + self._pv_system = building.energy_systems[0] + self._hvac_system = building.energy_systems[1] + self._dhw_system = building.energy_systems[-1] + self._dhw_peak_flow_rate = (building.thermal_zones_from_internal_zones[0].total_floor_area * + building.thermal_zones_from_internal_zones[0].domestic_hot_water.peak_flow * + cte.WATER_DENSITY) + self._heating_peak_load = building.heating_peak_load[cte.YEAR][0] + self._cooling_peak_load = building.cooling_peak_load[cte.YEAR][0] + self._domestic_hot_water_peak_load = building.domestic_hot_water_peak_load[cte.YEAR][0] + self._hourly_heating_demand = [demand / 3600 for demand in building.heating_demand[cte.HOUR]] + self._hourly_cooling_demand = [demand / 3600 for demand in building.cooling_demand[cte.HOUR]] + self._hourly_dhw_demand = [0] + building.domestic_hot_water_heat_demand[cte.HOUR] + self._output_path = output_path + self._t_out = building.external_temperature[cte.HOUR] + self.results = {} + + def hvac_sizing(self): + storage_factor = 3 + heat_pump = self._hvac_system.generation_systems[1] + boiler = self._hvac_system.generation_systems[0] + thermal_storage = boiler.energy_storage_systems[0] + heat_pump.nominal_heat_output = round(0.5 * self._heating_peak_load / 3600) + heat_pump.nominal_cooling_output = round(self._cooling_peak_load / 3600) + boiler.nominal_heat_output = round(0.5 * self._heating_peak_load / 3600) + thermal_storage.volume = round( + (self._heating_peak_load * storage_factor) / (cte.WATER_HEAT_CAPACITY * cte.WATER_DENSITY * 25)) + return heat_pump, boiler, thermal_storage + + def dhw_sizing(self): + storage_factor = 3 + dhw_hp = self._dhw_system.generation_systems[0] + dhw_hp.nominal_heat_output = 0.7 * self._domestic_hot_water_peak_load + dhw_hp.source_temperature = self._t_out + dhw_tes = dhw_hp.energy_storage_systems[0] + dhw_tes.volume = round( + (self._domestic_hot_water_peak_load * storage_factor * 3600) / (cte.WATER_HEAT_CAPACITY * cte.WATER_DENSITY * 10)) + return dhw_hp, dhw_tes + + def heating_system_simulation_stratified(self): + hp, boiler, tes = self.hvac_sizing() + cop_curve_coefficients = [float(coefficient) for coefficient in hp.heat_efficiency_curve.coefficients] + demand = [0] + [x for x in self._hourly_heating_demand for _ in range(12)] + hp.source_temperature = self._t_out + t_out = [0] + [x for x in self._t_out for _ in range(12)] + variable_names = ["t_sup_hp", "t1", "t2", "t3", "t4", "t_tank", "t_ret", "m_ch", "m_dis", "q_hp", "q_boiler", + "hp_cop", "hp_electricity", "boiler_gas_consumption", "t_sup_boiler", "boiler_energy_consumption", + "heating_consumption"] + num_hours = len(demand) + variables = {name: [0] * num_hours for name in variable_names} + (t_sup_hp, t1, t2, t3, t4, t_tank, t_ret, m_ch, m_dis, q_hp, q_boiler, hp_cop, + hp_electricity, boiler_gas_consumption, t_sup_boiler, boiler_energy_consumption, heating_consumption) = \ + [variables[name] for name in variable_names] + t_tank[0] = 55 + t1[0] = 55 + t2[0] = 55 + t3[0] = 55 + t4[0] = 55 + dt = 300 + hp_heating_cap = hp.nominal_heat_output + boiler_heating_cap = boiler.nominal_heat_output + hp_delta_t = 5 + boiler_efficiency = float(boiler.heat_efficiency) + v, h = float(tes.volume) / 4, float(tes.height) / 4 + r_tot = sum(float(layer.thickness) / float(layer.material.conductivity) for layer in + tes.layers) + u_tot = 1 / r_tot + d = math.sqrt((4 * v) / (math.pi * h)) + a_side = math.pi * d * h + a_top = math.pi * d ** 2 / 4 + ua_side = u_tot * a_side + ua_top_bottom = u_tot * (a_top + a_side) + # storage temperature prediction + for i in range(len(demand) - 1): + t1[i + 1] = t1[i] + ((m_ch[i] * (t_sup_boiler[i] - t1[i])) + ( + np.heaviside((m_dis[i] - m_ch[i]), 0) * (m_ch[i] - m_dis[i]) * (t1[i] - t2[i])) + ( + ua_top_bottom * (t_out[i] - t1[i])) / cte.WATER_HEAT_CAPACITY - cte.WATER_THERMAL_CONDUCTIVITY * (a_top * (t1[i] - t2[i])) / ( + cte.WATER_HEAT_CAPACITY * h)) * (dt / (cte.WATER_DENSITY * v)) + t2[i + 1] = t2[i] + ((np.heaviside((m_dis[i] - m_ch[i]), 0) * (m_ch[i] - m_dis[i]) * (t2[i] - t3[i])) + ( + ua_side * (t_out[i] - t2[i])) / cte.WATER_HEAT_CAPACITY - (cte.WATER_THERMAL_CONDUCTIVITY * (a_top * (t2[i] - t1[i])) / (cte.WATER_HEAT_CAPACITY * h)) - ( + cte.WATER_THERMAL_CONDUCTIVITY * (a_top * (t2[i] - t3[i])) / (cte.WATER_HEAT_CAPACITY * h)) + ( + np.heaviside((m_ch[i] - m_dis[i]), 0) * (m_ch[i] - m_dis[i]) * ( + t1[i] - t2[i]))) * (dt / (cte.WATER_DENSITY * v)) + t3[i + 1] = t3[i] + ((np.heaviside((m_dis[i] - m_ch[i]), 0) * (m_ch[i] - m_dis[i]) * (t3[i] - t4[i])) + ( + ua_side * (t_out[i] - t3[i])) / cte.WATER_HEAT_CAPACITY - (cte.WATER_THERMAL_CONDUCTIVITY * (a_top * (t3[i] - t2[i])) / (cte.WATER_HEAT_CAPACITY * h)) - ( + cte.WATER_THERMAL_CONDUCTIVITY * (a_top * (t3[i] - t4[i])) / (cte.WATER_HEAT_CAPACITY * h)) + ( + np.heaviside((m_ch[i] - m_dis[i]), 0) * (m_ch[i] - m_dis[i]) * ( + t2[i] - t3[i]))) * (dt / (cte.WATER_DENSITY * v)) + t4[i + 1] = t4[i] + (np.heaviside((m_ch[i] - m_dis[i]), 0) * ((m_ch[i] - m_dis[i]) * (t3[i] - t4[i])) + ( + ua_top_bottom * (t_out[i] - t4[-1])) / cte.WATER_HEAT_CAPACITY - m_dis[i] * ((t4[i] - t_ret[i])) - ( + cte.WATER_THERMAL_CONDUCTIVITY * (a_top * (t4[i] - t3[i])) / (cte.WATER_HEAT_CAPACITY * h))) * (dt / (cte.WATER_DENSITY * v)) + # hp operation + if t1[i + 1] < 40: + q_hp[i + 1] = hp_heating_cap + m_ch[i + 1] = q_hp[i + 1] / (cte.WATER_HEAT_CAPACITY * hp_delta_t) + t_sup_hp[i + 1] = (q_hp[i + 1] / (m_ch[i + 1] * cte.WATER_HEAT_CAPACITY)) + t4[i + 1] + elif 40 <= t1[i + 1] < 55 and q_hp[i] == 0: + q_hp[i + 1] = 0 + m_ch[i + 1] = 0 + t_sup_hp[i + 1] = t4[i + 1] + elif 40 <= t1[i + 1] < 55 and q_hp[i] > 0: + q_hp[i + 1] = hp_heating_cap + m_ch[i + 1] = q_hp[i + 1] / (cte.WATER_HEAT_CAPACITY * hp_delta_t) + t_sup_hp[i + 1] = (q_hp[i + 1] / (m_ch[i + 1] * cte.WATER_HEAT_CAPACITY)) + t4[i + 1] + else: + q_hp[i + 1], m_ch[i + 1], t_sup_hp[i + 1] = 0, 0, t4[i + 1] + t_sup_hp_fahrenheit = 1.8 * t_sup_hp[i + 1] + 32 + t_out_fahrenheit = 1.8 * t_out[i + 1] + 32 + if q_hp[i + 1] > 0: + hp_cop[i + 1] = (cop_curve_coefficients[0] + + cop_curve_coefficients[1] * t_sup_hp_fahrenheit + + cop_curve_coefficients[2] * t_sup_hp_fahrenheit ** 2 + + cop_curve_coefficients[3] * t_out_fahrenheit + + cop_curve_coefficients[4] * t_out_fahrenheit ** 2 + + cop_curve_coefficients[5] * t_sup_hp_fahrenheit * t_out_fahrenheit) + hp_electricity[i + 1] = q_hp[i + 1] / hp_cop[i + 1] + else: + hp_cop[i + 1] = 0 + hp_electricity[i + 1] = 0 + # boiler operation + if q_hp[i + 1] > 0: + if t_sup_hp[i + 1] < 45: + q_boiler[i + 1] = boiler_heating_cap + elif demand[i + 1] > 0.5 * self._heating_peak_load / dt: + q_boiler[i + 1] = 0.5 * boiler_heating_cap + boiler_energy_consumption[i + 1] = q_boiler[i + 1] / boiler_efficiency + boiler_gas_consumption[i + 1] = (q_boiler[i + 1] * dt) / (boiler_efficiency * cte.NATURAL_GAS_LHV) + t_sup_boiler[i + 1] = t_sup_hp[i + 1] + (q_boiler[i + 1] / (m_ch[i + 1] * cte.WATER_HEAT_CAPACITY)) + # storage discharging + if demand[i + 1] == 0: + m_dis[i + 1] = 0 + t_ret[i + 1] = t1[i + 1] + else: + if demand[i + 1] > 0.5 * self._heating_peak_load / cte.HOUR_TO_SECONDS: + factor = 8 + else: + factor = 4 + m_dis[i + 1] = self._heating_peak_load / (cte.WATER_HEAT_CAPACITY * factor * cte.HOUR_TO_SECONDS) + t_ret[i + 1] = t1[i + 1] - demand[i + 1] / (m_dis[i + 1] * cte.WATER_HEAT_CAPACITY) + + hp_electricity_wh = [x / 12 for x in hp_electricity] + boiler_consumption_wh = [x / 12 for x in boiler_energy_consumption] + hp_hourly = [] + boiler_hourly = [] + tes.temperature = {} + tes.temperature['layer_1'] = [] + tes.temperature['layer_2'] = [] + tes.temperature['layer_3'] = [] + tes.temperature['layer_4'] = [] + for i in range(1, len(demand), 12): + tes.temperature['layer_1'].append(t1[i]) + tes.temperature['layer_2'].append(t2[i]) + tes.temperature['layer_3'].append(t3[i]) + tes.temperature['layer_4'].append(t4[i]) + demand_modified = demand[1:] + hp_hourly.append(hp_electricity[1]) + boiler_hourly.append(boiler_energy_consumption[1]) + boiler_sum = 0 + hp_sum = 0 + for i in range(1, len(demand_modified) + 1): + hp_sum += hp_electricity_wh[i] + boiler_sum += boiler_consumption_wh[i] + if i % 12 == 0: + hp_hourly.append(hp_sum) + boiler_hourly.append(boiler_sum) + hp_sum = 0 + boiler_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])] + + self.results['Heating Demand (W)'] = demand + self.results['HP Heat Output (W)'] = q_hp + self.results['HP Source Temperature'] = t_out + self.results['HP Supply Temperature'] = t_sup_hp + self.results['HP COP'] = hp_cop + self.results['HP Electricity Consumption (W)'] = hp_electricity + self.results['Boiler Heat Output (W)'] = q_boiler + self.results['Boiler Supply Temperature'] = t_sup_boiler + self.results['Boiler Gas Consumption'] = boiler_gas_consumption + self.results['TES Layer 1 Temperature'] = t1 + self.results['TES Layer 2 Temperature'] = t2 + self.results['TES Layer 3 Temperature'] = t3 + self.results['TES Layer 4 Temperature'] = t4 + self.results['TES Charging Flow Rate (kg/s)'] = m_ch + self.results['TES Discharge Flow Rate (kg/s)'] = m_dis + self.results['Heating Loop Return Temperature'] = t_ret + return hp_electricity, boiler_energy_consumption + + def cooling_system_simulation(self): + hp = self.hvac_sizing()[0] + eer_curve_coefficients = [float(coefficient) for coefficient in hp.cooling_efficiency_curve.coefficients] + cooling_efficiency = float(hp.cooling_efficiency) + demand = self._hourly_cooling_demand + hp.source_temperature = self._t_out + variable_names = ["t_sup_hp", "t_ret", "m", "q_hp", "hp_electricity", "hp_eer"] + num_hours = len(demand) + variables = {name: [0] * num_hours for name in variable_names} + (t_sup_hp, t_ret, m, q_hp, hp_electricity, hp_eer) = [variables[name] for name in variable_names] + t_ret[0] = 13 + dt = 3600 + for i in range(len(demand) - 1): + if demand[i] > 0: + m[i] = self._cooling_peak_load / (cte.WATER_HEAT_CAPACITY * 5 * dt) + if t_ret[i] > 13: + if demand[i] < 0.25 * self._cooling_peak_load / dt: + q_hp[i] = 0.25 * hp.nominal_cooling_output + elif demand[i] < 0.5 * self._cooling_peak_load / dt: + q_hp[i] = 0.5 * hp.nominal_cooling_output + else: + q_hp[i] = hp.nominal_cooling_output + t_sup_hp[i] = t_ret[i] - q_hp[i] / (m[i] * cte.WATER_HEAT_CAPACITY) + else: + q_hp[i] = 0 + t_sup_hp[i] = t_ret[i] + t_ret[i + 1] = t_sup_hp[i] + demand[i] / (m[i] * cte.WATER_HEAT_CAPACITY) + else: + m[i] = 0 + q_hp[i] = 0 + t_sup_hp[i] = t_ret[i] + t_ret[i + 1] = t_ret[i] + t_sup_hp_fahrenheit = 1.8 * t_sup_hp[i] + 32 + t_out_fahrenheit = 1.8 * self._t_out[i] + 32 + if q_hp[i] > 0: + hp_eer[i] = (eer_curve_coefficients[0] + + eer_curve_coefficients[1] * t_sup_hp_fahrenheit + + eer_curve_coefficients[2] * t_sup_hp_fahrenheit ** 2 + + eer_curve_coefficients[3] * t_out_fahrenheit + + eer_curve_coefficients[4] * t_out_fahrenheit ** 2 + + eer_curve_coefficients[5] * t_sup_hp_fahrenheit * t_out_fahrenheit) + hp_electricity[i] = q_hp[i] / cooling_efficiency + else: + hp_eer[i] = 0 + hp_electricity[i] = 0 + hp.energy_consumption[cte.COOLING] = {} + hp.energy_consumption[cte.COOLING][cte.HOUR] = hp_electricity + hp.energy_consumption[cte.COOLING][cte.MONTH] = MonthlyValues.get_total_month( + hp.energy_consumption[cte.COOLING][cte.HOUR]) + hp.energy_consumption[cte.COOLING][cte.YEAR] = [ + sum(hp.energy_consumption[cte.COOLING][cte.MONTH])] + # self.results['Cooling Demand (W)'] = demand + # self.results['HP Cooling Output (W)'] = q_hp + # self.results['HP Cooling Supply Temperature'] = t_sup_hp + # self.results['HP Cooling COP'] = hp_eer + # self.results['HP Electricity Consumption'] = hp_electricity + # self.results['Cooling Loop Flow Rate (kg/s)'] = m + # self.results['Cooling Loop Return Temperature'] = t_ret + return hp_electricity + + def dhw_system_simulation(self): + hp, tes = self.dhw_sizing() + cop_curve_coefficients = [float(coefficient) for coefficient in hp.heat_efficiency_curve.coefficients] + demand = self._hourly_dhw_demand + 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)"] + num_hours = len(demand) + variables = {name: [0] * num_hours for name in variable_names} + (t_sup_hp, t_tank, m_ch, m_dis, m_refill, q_hp, q_coil, hp_cop, hp_electricity, v_dhw) = \ + [variables[name] for name in variable_names] + t_tank[0] = 70 + v_dhw[0] = tes.volume + dt = 3600 + hp_heating_cap = hp.nominal_heat_output + hp_delta_t = 8 + v, h = float(tes.volume), float(tes.height) + r_tot = sum(float(layer.thickness) / float(layer.material.conductivity) for layer in + tes.layers) + u_tot = 1 / r_tot + d = math.sqrt((4 * v) / (math.pi * h)) + a_side = math.pi * d * h + a_top = math.pi * d ** 2 / 4 + ua = u_tot * (2 * a_top + a_side) + freshwater_temperature = 18 + for i in range(len(demand) - 1): + delta_t_demand = demand[i] * (dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) + if t_tank[i] < 65: + q_hp[i] = hp_heating_cap + delta_t_hp = q_hp[i] * (dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) + if demand[i] > 0: + dhw_needed = (demand[i] * cte.HOUR_TO_SECONDS) / (cte.WATER_HEAT_CAPACITY * t_tank[i] * cte.WATER_DENSITY) + m_dis[i] = dhw_needed * cte.WATER_DENSITY / cte.HOUR_TO_SECONDS + m_refill[i] = m_dis[i] + delta_t_freshwater = m_refill[i] * (t_tank[i] - freshwater_temperature) * (dt / (v * cte.WATER_DENSITY)) + diff = delta_t_freshwater + delta_t_demand - delta_t_hp + if diff > 0: + if diff > 0: + power = diff * (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v) / dt + if power <= float(tes.heating_coil_capacity): + q_coil[i] = power + else: + q_coil[i] = float(tes.heating_coil_capacity) + elif t_tank[i] < 65: + q_coil[i] = float(tes.heating_coil_capacity) + delta_t_coil = q_coil[i] * (dt / (cte.WATER_DENSITY * cte.WATER_HEAT_CAPACITY * v)) + + if q_hp[i] > 0: + m_ch[i] = q_hp[i] / (cte.WATER_HEAT_CAPACITY * hp_delta_t) + t_sup_hp[i] = (q_hp[i] / (m_ch[i] * cte.WATER_HEAT_CAPACITY)) + t_tank[i] + else: + m_ch[i] = 0 + t_sup_hp[i] = t_tank[i] + t_sup_hp_fahrenheit = 1.8 * t_sup_hp[i] + 32 + t_out_fahrenheit = 1.8 * self._t_out[i] + 32 + if q_hp[i] > 0: + hp_cop[i] = (cop_curve_coefficients[0] + + cop_curve_coefficients[1] * t_sup_hp_fahrenheit + + cop_curve_coefficients[2] * t_sup_hp_fahrenheit ** 2 + + cop_curve_coefficients[3] * t_out_fahrenheit + + cop_curve_coefficients[4] * t_out_fahrenheit ** 2 + + cop_curve_coefficients[5] * t_sup_hp_fahrenheit * t_out_fahrenheit) + hp_electricity[i] = q_hp[i] / 3.5 + else: + hp_cop[i] = 0 + hp_electricity[i] = 0 + + t_tank[i + 1] = t_tank[i] + (delta_t_hp - delta_t_freshwater - delta_t_demand + delta_t_coil) + + hp.energy_consumption[cte.DOMESTIC_HOT_WATER] = {} + hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.HOUR] = hp_electricity + hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.MONTH] = MonthlyValues.get_total_month( + hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.HOUR]) + hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.YEAR] = [ + sum(hp.energy_consumption[cte.DOMESTIC_HOT_WATER][cte.MONTH])] + tes.heating_coil_energy_consumption = {} + tes.heating_coil_energy_consumption[cte.HOUR] = q_coil + tes.heating_coil_energy_consumption[cte.MONTH] = MonthlyValues.get_total_month( + tes.heating_coil_energy_consumption[cte.HOUR]) + tes.heating_coil_energy_consumption[cte.YEAR] = [ + sum(tes.heating_coil_energy_consumption[cte.MONTH])] + tes.temperature = t_tank + + + # self.results['DHW Demand (W)'] = demand + # self.results['DHW HP Heat Output (W)'] = q_hp + # self.results['DHW HP Electricity Consumption (W)'] = hp_electricity + # self.results['DHW HP Source Temperature'] = self._t_out + # self.results['DHW HP Supply Temperature'] = t_sup_hp + # self.results['DHW HP COP'] = hp_cop + # self.results['DHW TES Heating Coil Heat Output (W)'] = q_coil + # self.results['DHW TES Temperature'] = t_tank + # self.results['DHW TES Charging Flow Rate (kg/s)'] = m_ch + # self.results['DHW Flow Rate (kg/s)'] = m_dis + # self.results['DHW TES Refill Flow Rate (kg/s)'] = m_refill + # self.results['Available Water in Tank (m3)'] = v_dhw + return hp_electricity, q_coil + + def enrich_buildings(self): + hp_heating, boiler_consumption = self.heating_system_simulation_stratified() + # hp_cooling = self.cooling_system_simulation() + # hp_dhw, heating_coil = self.dhw_system_simulation() + heating_consumption = [hp_heating[i] + boiler_consumption[i] for i in range(len(hp_heating))] + # dhw_consumption = [hp_dhw[i] + heating_coil[i] for i in range(len(hp_dhw))] + # self._building.heating_consumption[cte.HOUR] = heating_consumption + # self._building.heating_consumption[cte.MONTH] = ( + # MonthlyValues.get_total_month(self._building.heating_consumption[cte.HOUR])) + # self._building.heating_consumption[cte.YEAR] = sum(self._building.heating_consumption[cte.MONTH]) + # self._building.cooling_consumption[cte.HOUR] = hp_cooling + # self._building.cooling_consumption[cte.MONTH] = ( + # MonthlyValues.get_total_month(self._building.cooling_consumption[cte.HOUR])) + # self._building.cooling_consumption[cte.YEAR] = sum(self._building.cooling_consumption[cte.MONTH]) + # self._building.domestic_hot_water_consumption[cte.HOUR] = dhw_consumption + # self._building.domestic_hot_water_consumption[cte.MONTH] = ( + # 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()))