Select Git revision
CMakeLists.txt
demand_generator.py 10.51 KiB
"""
MIT License
Copyright (c) 2023 RWTH Aachen University
Permission is hereby granted, free of charge, to any person obtaining a copy
of this software and associated documentation files (the "Software"), to deal
in the Software without restriction, including without limitation the rights
to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
copies of the Software, and to permit persons to whom the Software is
furnished to do so, subject to the following conditions:
The above copyright notice and this permission notice shall be included in
all copies or substantial portions of the Software.
THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
THE SOFTWARE.
"""
import calendar
import datetime
import numpy as np
import os
import pandas as pd
# based on system described in https://www.bdew.de/media/documents/2000131_Anwendung-repraesentativen_Lastprofile-Step-by-step.pdf
# profiles from https://www.bdew.de/energie/standardlastprofile-strom/
class ElectricityDemandGenerator:
def __init__(self):
# load profiles (Hardcoded H0: Haushalt)
# TODO: sort the csv by its index
profiles_df = pd.read_csv(
os.path.join(os.path.dirname(__file__), "selp.csv"), index_col=[0, 1]
)
self.profiles = dict()
profile_columns = profiles_df.columns
#loading profiles into self.profiles dictionary, indexed py period, weekday and profile
for profile in profile_columns:
for period in ["winter", "transition", "summer"]:
for weekday in ["working_day", "saturday", "sunday"]:
self.profiles[period, weekday, profile] = profiles_df.loc[
(period, weekday), profile
].values
# set standart period ranges
self.period_ranges = [
["winter", 1, 1, 3, 20],
["transition", 3, 21, 5, 14],
["summer", 5, 15, 9, 14],
["transition", 9, 15, 10, 31],
["winter", 11, 1, 12, 31],
]
# precompute dynamic factor
factor_function = (
lambda t: -3.92 * 10**-10 * t**4
+ 3.2 * 10**-7 * t**3
- 7.02 * 10**-5 * t**2
+ 2.1 * 10**-3 * t
+ 1.24
)
self.dynamic_factor = factor_function(np.arange(1, 366, dtype=float))
def generate(self, year, annual_demand, profile_type="H0"):
# year: year to generate the data for in astronomical numbering
# annual_demand: annual electricity demand in kWh
# days in the year, depending on wehter the given year a leap year
days_in_year = 365 if not calendar.isleap(year) else 366
# generate series of days in the year and to which period they belong
periods = pd.Series(
index=pd.date_range(datetime.datetime(year, 1, 1), periods=days_in_year),
dtype=str,
)
for period, month_start, day_start, month_end, day_end in self.period_ranges:
periods[
datetime.datetime(year, month_start, day_start) : datetime.datetime(
year, month_end, day_end
)
] = period
# construct profile and apply dynamic factor
profile = np.empty(days_in_year * 24 * 4, dtype=float)
for i, (date, period) in enumerate(periods.items()):
weekday = date.weekday()
if weekday < 5:
profile[i * 24 * 4 : (i + 1) * 24 * 4] = (
self.profiles[period, "working_day", profile_type] * self.dynamic_factor[i]
)
elif weekday == 5:
profile[i * 24 * 4 : (i + 1) * 24 * 4] = (
self.profiles[period, "saturday", profile_type] * self.dynamic_factor[i]
)
else:
profile[i * 24 * 4 : (i + 1) * 24 * 4] = (
self.profiles[period, "sunday", profile_type] * self.dynamic_factor[i]
)
# normalize profile and convert to kW
profile *= annual_demand / (1000.0 * 1000.0)
return profile
# based on system descibed in https://www.bdew.de/media/documents/Leitfaden_20160630_Abwicklung-Standardlastprofile-Gas.pdf (calculation of the reference value: page 108, calculation of the demand profile: page 47)
# profile sources depend of building class (see table 6 (page 163) in https://www.bdew.de/media/documents/Leitfaden_20160630_Abwicklung-Standardlastprofile-Gas.pdf)
# profile source TUM2002 https://mediatum.ub.tum.de/doc/601557/601557.pdf
# profile source TUM2006 http://www.gwb-netz.de/wa_files/05_bgw_leitfaden_lastprofile_56550.pdf
# profile source P2007/13 https://www.eichsfeldwerke.de/fileadmin/user_upload/Praxisinformation_P2007_13.pdf
class HeatingDemandGenerator:
def __init__(self):
# load factors (Hardcoded D13/DE_HEF03: Deutschland, bundesweit, Einfamilienhaushalt, Ausprägung o)
profiles_df = pd.read_csv(
os.path.join(os.path.dirname(__file__), "shlp.csv"), index_col=0
)
hour_factors_df = pd.read_csv(
os.path.join(os.path.dirname(__file__), "shlp_hour_factors.csv"),
index_col=[0, 1],
)
self.a = profiles_df["A"]
self.b = profiles_df["B"]
self.c = profiles_df["C"]
self.d = profiles_df["D"]
self.theta_0 = profiles_df["theta_0"]
self.weekday_factors = dict()
for profile in profiles_df.index:
self.weekday_factors[profile] = profiles_df.loc[
profile, ["F_Mon", "F_Tue", "F_Wed", "F_Thu", "F_Fri", "F_Sat", "F_Sun"]
].values
# self.weekday_factors = profiles_df.loc[
# "D13", ["F_Mon", "F_Tue", "F_Wed", "F_Thu", "F_Fri", "F_Sat", "F_Sun"]
# ].values
self.hour_factors = dict()
for profile in hour_factors_df.index.get_level_values(0).unique():
self.hour_factors[profile] = hour_factors_df.loc[profile].values
#self.hour_factors = hour_factors_df.loc["D13"].values
def generate(self, year, annual_demand, temperatures, seperate_warm_water=False, profile_type="D13"):
# year: year to generate the data for in astronomical numbering
# annual_demand: annual heating demand in kWh
# temperatures: daily temperatures in °C
# seperate_warm_water: if set to false, the function returns one heating demand profile, if set to true, the function returns two demand profiles, the first is the heating demand withoud the warm water demand and the second is the warm water demand
# days in the year, depending on wehter the given year a leap year
days_in_year = 365 if not calendar.isleap(year) else 366
if len(temperatures) != days_in_year:
raise ValueError(f"Lengh of daily tamperatures has to be {days_in_year}!")
# compute geometric mean of the temperatures
allocation_temperatures = (
temperatures
+ 0.5 * np.roll(temperatures, 1)
+ 0.25 * np.roll(temperatures, 2)
+ 0.125 * np.roll(temperatures, 3)
) / 1.875
# compute sigmoid function, weekday factors and hour factors
temp = self.b[profile_type] / (allocation_temperatures - 40.0)
h = (self.a[profile_type] / (1.0 + np.sign(temp) * (np.abs(temp) ** self.c[profile_type]))) + self.d[profile_type]
f_wd = np.array(
pd.date_range(datetime.datetime(year, 1, 1), periods=days_in_year).map(
lambda a: self.weekday_factors[profile_type][a.weekday()]
),
dtype=float,
)
hour_factors = self.hour_factors[profile_type] if profile_type in self.hour_factors.keys() else self.hour_factors["D13"]
f_h = hour_factors[
(
(allocation_temperatures <= -15.0) * 1
+ (
(-15.0 < allocation_temperatures)
& (allocation_temperatures <= -10.0)
)
* 2
+ (
(-10.0 < allocation_temperatures)
& (allocation_temperatures <= -5.0)
)
* 3
+ ((-5.0 < allocation_temperatures) & (allocation_temperatures <= 0.0))
* 4
+ ((0.0 < allocation_temperatures) & (allocation_temperatures <= 5.0))
* 5
+ ((5.0 < allocation_temperatures) & (allocation_temperatures <= 10.0))
* 6
+ ((10.0 < allocation_temperatures) & (allocation_temperatures <= 15.0))
* 7
+ ((15.0 < allocation_temperatures) & (allocation_temperatures <= 20.0))
* 8
+ ((20.0 < allocation_temperatures) & (allocation_temperatures <= 25.0))
* 9
+ (25.0 < allocation_temperatures) * 10
)
- 1
]
# compute reference daily demand (kundenwert)
reference_temp = self.b[profile_type] / (np.full(days_in_year, 8.0, dtype=float) - 40.0)
reference_h = (
self.a[profile_type]
/ (1.0 + np.sign(reference_temp) * (np.abs(reference_temp) ** self.c[profile_type]))
) + self.d[profile_type]
reference_daily_demand = annual_demand / np.sum(reference_h * f_wd)
# construct profiles
if seperate_warm_water:
profile_sum = np.empty(days_in_year * 24, dtype=float)
profile_warm_water = np.empty(days_in_year * 24, dtype=float)
for i in range(days_in_year):
profile_sum[i * 24 : (i + 1) * 24] = (
reference_daily_demand * h[i] * f_wd[i] * (f_h[i] / 100.0)
)
profile_warm_water[i * 24 : (i + 1) * 24] = (
reference_daily_demand * self.d[profile_type] * f_wd[i] * (f_h[i] / 100.0)
)
return profile_sum - profile_warm_water, profile_warm_water
else:
profile = np.empty(days_in_year * 24, dtype=float)
for i in range(days_in_year):
profile[i * 24 : (i + 1) * 24] = (
reference_daily_demand * h[i] * f_wd[i] * (f_h[i] / 100.0)
)
return profile