You cannot select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
Mooi-Kickstart/Simulations/notepad2.py

412 lines
14 KiB
Python

import warnings
from functools import partial, lru_cache
from numbers import Number
from itertools import count
import numpy as np
from numpy.polynomial import Polynomial
from scipy.optimize import minimize_scalar
# from .converters import *
class Asset:
"""Generic class for producing/consuming assets. Specific asset classes can
inherit from this class.
Parameters:
-----------
max_power : int/float
Maximum asset power in MW electric
min_power : int/float
Minimium asset load in MW electric
Usage:
------
Use the set_load and get_load methods to set and get asset status in MW.
Convention is negative values for inputs (consumption) and positive
values for outputs (production).
"""
_freq_to_multiplier = {"H": 1, "15T": (1 / 4), "1T": (1 / 60)}
_ids = count(0)
def __init__(self, name, max_power, min_power):
if min_power > max_power:
raise ValueError("'min_power' can not be larger than 'max_power'.")
self.name = name
self.id = next(self._ids)
self.max_power = max_power
self.min_power = min_power
self.modes = {"max": max_power, "min": min_power}
def __repr__(self):
return f"{self.__class__.__name__}(self, max_power={self.max_power}, min_power={self.min_power})"
def set_load(self, load):
"""Set Asset load in MW.
Convention is negative value for consumption and positive value
for production. Subclasses might use a different convention if
this seems more intiutive.
Returns the load that is set in MW.
"""
if load < self.min_power or load > self.max_power:
warnings.warn(
f"Chosen Asset load for {self.name} is out of range. "
f"Should be between {self.min_power} and {self.max_power}. "
f"Function will return boundary load level for now."
)
load = min(max(load, self.min_power), self.max_power)
return load
def set_mode(self, mode):
""" """
load = self.modes[mode]
return self.set_load(load)
def MW_to_MWh(self, MW):
"""Performs conversion from MW to MWh using the time_factor variable."""
return MW * self.time_factor
def MWh_to_MW(self, MWh):
"""Performs conversion from MWh to MW using the time_factor variable."""
return MWh / self.time_factor
def set_freq(self, freq):
"""
Function that aligns time frequency between Model and Asset.
Can be '1T', '15T' or 'H'
The time_factor variable is used in subclasses to perform MW to MWh conversions.
"""
self.freq = freq
self.time_factor = Asset._freq_to_multiplier[freq]
def set_financials(self, capex, opex, devex, lifetime=None, depreciate=True, salvage_value=0):
"""Set financial data of the asset."""
self.capex = capex
self.opex = opex
self.devex = devex
self.lifetime = lifetime
self.depreciate = depreciate
self.salvage_value = salvage_value
class WaterStorage(Asset):
def __init__(
self,
name,
max_power,
min_power,
roundtrip_eff,
capacity_per_volume,
volume,
lifetime,
temperature,
min_storagelevel,
initial_storagelevel=None
):
if min_power > max_power:
raise ValueError("'min_power' can not be larger than 'max_power'.")
super().__init__(name, max_power, min_power)
self.roundtrip_eff = roundtrip_eff
self.volume = volume
self.lifetime = lifetime
self.temperature = temperature
self.capacity_per_volume = capacity_per_volume #MWh/m3
self.max_storage_capacity = self.volume * self.capacity_per_volume #MWh
self.max_storagelevel = self.max_storage_capacity * 0.95
self.min_storagelevel = min_storagelevel
if initial_storagelevel:
self.storagelevel = initial_storagelevel
else:
self.storagelevel = self.min_storagelevel
def set_storagelevel(self, storagelevel):
"""Set the storagelevel in MWh."""
if storagelevel < self.min_storagelevel or storagelevel > self.max_storagelevel:
raise ValueError(
f"Tried to set Storage Level to {storagelevel}. "
f"Storage Level must be a value between "
f"{self.min_storagelevel} and {self.max_storagelevel} (in MWh)"
)
self.storagelevel = storagelevel
@property
def charging_power_limit(self):
max_charging_energy = self.max_storagelevel - self.storagelevel
max_charging_power = min(self.MWh_to_MW(max_charging_energy), -self.min_power)
return max_charging_power
@property
def discharging_power_limit(self):
max_discharging_energy = self.storagelevel - self.min_storagelevel
max_discharging_power = min(self.MWh_to_MW(max_discharging_energy), self.max_power)
return max_discharging_power
def charge(self, power):
energy = self.MW_to_MWh(power)
max_charging = self.max_storagelevel - self.storagelevel
bounded_energy = min (energy, max_charging)
# print('bounded_energy', bounded_energy)
new_cl = self.storagelevel + bounded_energy
# print('new_cl', new_cl)
self.set_storagelevel(new_cl)
power = self.MWh_to_MW(bounded_energy)
return power
def discharge(self, power):
energy = self.MW_to_MWh(power)
max_discharging = self.storagelevel - self.min_storagelevel
bounded_energy = min(energy, max_discharging)
# print('bounded_energy', bounded_energy)
new_cl = self.storagelevel - bounded_energy
# print('new_cl', new_cl)
self.set_storagelevel(new_cl)
power = self.MWh_to_MW(bounded_energy)
return power
def get_soc(self, storagelevel, max_storage_capacity):
"""Get the SoC in % (decimal value)"""
return self.storagelevel / self.max_storage_capacity
def set_financials(self,
capex_per_MW,
capex_per_MWh,
opex,
devex,
lifetime=None,
depreciate=True, salvage_value=0):
total_capex = (
capex_per_MW * self.max_power + capex_per_MWh * self.max_storage_capacity
)
super().set_financials(
total_capex,
opex,
devex,
lifetime=None,
depreciate=True,
salvage_value=0
)
def __repr__(self):
return (
f"{self.__class__.__name__}(name={self.name}, max_power={self.max_power}, "
f"min_power={self.min_power}, roundtrip_eff={self.roundtrip_eff}, capacity_per_volume={self.capacity_per_volume}, volume = {self.volume}),"
f"lifetime={self.lifetime}, temperature = {self.temperature}, min_storagelevel = {self.min_storagelevel})"
)
class Heatpump(Asset):
"""Subclass for a Heatpump.
Use cop parameter to set fixed COP (float/int) or COP curve (func).
COP curve should take load in MWhe and return COP.
Parameters:
-----------
max_th_power : numeric
Maximum thermal output in MW (positive value)
cop_curve : numeric or list or function
3 ways to set the COP of the Heatpump:
(1) Fixed COP based on [numeric] value.
(2) Polynomial with coefficients based on [list] input.
Input coeficients in format [c0, c1, c2, ..., c(n)],
will generate Polynomial p(x) = c0 + c1*x + c2*x^2 ... cn*x^n,
where x = % thermal load (in % of thermal capacity) as decimal value.
Example:
cop=[1, 2, 3, 4] will result in following COP curve:
p(x) = 1 + 2x + 3x**2 + 4x**3,
(3) [function] in format func(*args, **kwargs)
Function should return a Polynomial that takes 'load_perc' as parameter.
min_th_power : numeric
Minimum thermal output in MW (positive value)
Notes:
------
Sign convention:
Thermal power outputs have positive values
Electric power inputs have negative values
"""
def __init__(
self,
name,
max_th_power,
cop_curve,
min_th_power=0,
):
if max_th_power < 0 or min_th_power < 0:
raise ValueError("Thermal power can not have negative values.")
if min_th_power > max_th_power:
raise ValueError("'min_th_power' can not be larger than 'max_th_power'.")
self.name = name
self.max_th_power = max_th_power
self.min_th_power = min_th_power
self.cop_curve = self._set_cop_curve(cop_curve=cop_curve)
def __repr__(self):
return (
f"{self.__class__.__name__}(name='{self.name}', max_thermal_power={self.max_th_power}, "
f"cop_curve={self.cop_curve}, min_th_power={self.min_th_power})"
)
# Is turning everything into a Polynomial the best solution here?
def _set_cop_curve(self, cop_curve):
"""Generate COP curve function based on different inputtypes.
Returns a function that takes *args **kwargs and returns a Polynomial.
"""
if isinstance(cop_curve, list):
def func(*args, **kwargs):
return Polynomial(cop_curve)
return func
return cop_curve
@lru_cache(maxsize=None)
def get_cop(self, heat_output, Tsink=None, Tsource=None):
"""Get COP corresponding to certain load.
Parameters:
-----------
heat_output : numeric
Thermal load in MW
Tsink : numeric
Sink temperature in degrees celcius
Tsource : numeric
Source temperature in degrees celcius
Notes:
------
Sign convention:
Positive values for thermal load
Negative values for electric load
"""
load_perc = heat_output / self.max_th_power
cop_curve = self.cop_curve
if not callable(cop_curve):
return cop_curve
else:
return cop_curve(Tsink=Tsink, Tsource=Tsource)(load_perc)
def th_to_el_power(self, heat_output, Tsink=None, Tsource=None):
if not self.min_th_power <= heat_output <= self.max_th_power:
warnings.warn(
f"Chosen heat output is out of range [{self.min_th_power} - {self.max_th_power}]. "
"Heat output is being limited to the closest boundary."
)
heat_output = min(max(heat_output, self.min_th_power), self.max_th_power)
cop = self.get_cop(heat_output=heat_output, Tsink=Tsink, Tsource=Tsource)
return -heat_output / cop
def set_load(self, *args, **kwargs):
raise NotImplementedError(
"Directly setting the electric load of the heatpump is not possible (yet). "
"Functionality will be implemented if there is a specific usecase for it."
)
@lru_cache(maxsize=None)
def set_heat_output(self, heat_output, Tsink=None, Tsource=None):
"""Set heat output in MWth, returns load of heatpump as tuple (MWe, MWth)"""
if not self.min_th_power <= heat_output <= self.max_th_power:
warnings.warn(
f"Chosen heat output is out of range [{self.min_th_power} - {self.max_th_power}]. "
"Heat output is being limited to the closest boundary."
)
heat_output = min(max(heat_output, self.min_th_power), self.max_th_power)
# if Tsink is not None and Tsource is not None and Tsink <= Tsource:
# raise ValueError(f"Tsource '{Tsource}' can not be higher than '{Tsink}'.")
cop = self.get_cop(heat_output=heat_output, Tsink=Tsink, Tsource=Tsource)
e_load = -heat_output / cop
return e_load, heat_output
def _cost_function(self, x, c1, c2, c3, Tsink=None, Tsource=None):
"""Objective function for set_opt_load function.
x = heatpump thermal load in MW
c1 = electricity_cost
c2 = alt_heat_price
c3 = demand
"""
return (
x / self.get_cop(heat_output=x, Tsink=Tsink, Tsource=Tsource) * c1
+ (c3 - x) * c2
)
@lru_cache(maxsize=None)
def set_opt_load(
self,
electricity_cost,
alt_heat_price,
demand,
Tsink=None,
Tsource=None,
tolerance=0.01,
):
"""Set optimal load of Heatpump with minimal total heat costs.
Function uses np.minimize_scalar to minimize cost function.
Parameters:
-----------
electricity_cost:
Cost of input electricity in €/MWh(e)
alt_heat_price:
Price of heat from alternative source in €/MWh(th)
demand:
Heat demand in MW(th)
Returns:
--------
Optimal load of heatpump as tuple (MWe, MWth)
"""
c1 = electricity_cost
c2 = alt_heat_price
c3 = demand
cop_curve = self.cop_curve
if isinstance(cop_curve, Number):
if c1 / cop_curve <= c2:
return self.max_th_power
else:
return self.min_th_power
obj_func = partial(
self._cost_function, c1=c1, c2=c2, c3=c3, Tsink=Tsink, Tsource=Tsource
)
low_bound = 0
up_bound = min(c3, self.max_th_power)
opt_th_load = minimize_scalar(
obj_func,
bounds=(low_bound, up_bound),
method="bounded",
options={"xatol": tolerance},
).x
opt_e_load, opt_th_load = self.set_heat_output(
opt_th_load, Tsink=Tsink, Tsource=Tsource
)
return opt_e_load, opt_th_load