feat(optimization): model battery LCOS and probabilistic bypass

This commit is contained in:
Andreas
2026-07-15 09:23:29 +02:00
parent 8ddb7ce754
commit bed1f0f275
27 changed files with 1362 additions and 890 deletions
+5 -1
View File
@@ -50,8 +50,12 @@ class BatteriesCommonSettings(DevicesBaseSettings):
levelized_cost_of_storage_kwh: float = Field(
default=0.0,
ge=0.0,
json_schema_extra={
"description": "Levelized cost of storage (LCOS), the average lifetime cost of delivering one kWh [€/kWh].",
"description": (
"Levelized cost of storage (LCOS), applied once to each kWh delivered "
"by the battery [€/kWh]."
),
"examples": [0.12],
},
)
@@ -34,6 +34,11 @@ class Battery:
self.initial_soc_percentage = self.parameters.initial_soc_percentage
self.charging_efficiency = self.parameters.charging_efficiency
self.discharging_efficiency = self.parameters.discharging_efficiency
self.levelized_cost_of_storage_kwh = (
self.parameters.levelized_cost_of_storage_kwh
if isinstance(self.parameters, SolarPanelBatteryParameters)
else 0.0
)
# Charge rates, in case of None use default
self.charge_rates = np.array(BATTERY_DEFAULT_CHARGE_RATES, dtype=float)
@@ -115,6 +120,10 @@ class Battery:
raw_soc_available_wh = max(self.soc_wh - self.min_soc_wh, 0.0)
return min(raw_power_remaining_wh, raw_soc_available_wh) * self.discharging_efficiency
def discharged_energy_wh(self, hour: int) -> float:
"""Return DC energy delivered by the battery in one optimization slot."""
return self._discharged_raw_wh_per_slot[hour] * self.discharging_efficiency
def set_discharge_per_hour(self, discharge_array: np.ndarray) -> None:
"""Sets the discharge values for each hour."""
if len(discharge_array) != self.prediction_hours:
+78 -126
View File
@@ -14,9 +14,6 @@ class Inverter:
battery: Optional[Battery] = None,
slot_duration_h: float = 1.0,
):
# slot_duration_h scales the per-slot energy cap (max_power_wh). It
# defaults to 1.0, which keeps the hourly behaviour for the default
# optimization interval of 3600 s.
self.parameters: InverterParameters = parameters
self.battery: Optional[Battery] = battery
self.slot_duration_h: float = slot_duration_h
@@ -28,16 +25,13 @@ class Inverter:
logger.error(error_msg)
raise ValueError(error_msg)
self.self_consumption_predictor = get_eos_load_interpolator()
# max_power_wh is supplied as a power [W] that the legacy hourly code
# treats as Wh-per-hour. Scale it to the actual slot length so a 15-min
# slot can move at most a quarter of that energy.
self.max_power_wh = (
self.parameters.max_power_wh * self.slot_duration_h
) # Maximum energy the inverter can move in one optimization slot
# max_power_wh is supplied as power [W] but used as the maximum energy
# the inverter can move during one optimization slot.
self.max_power_wh = self.parameters.max_power_wh * self.slot_duration_h
self.dc_to_ac_efficiency = self.parameters.dc_to_ac_efficiency
self.ac_to_dc_efficiency = self.parameters.ac_to_dc_efficiency
# max_ac_charge_power_w stays in Watts. It feeds a dimensionless,
# slot-agnostic power-ratio cap in genetic.py simulate().
# This value remains a power [W]. GeneticSimulation converts it into a
# slot-independent charge-factor limit.
self.max_ac_charge_power_w = self.parameters.max_ac_charge_power_w
def _discharge_battery_to_ac(self, requested_ac_wh: float, hour: int) -> tuple[float, float]:
@@ -58,128 +52,86 @@ class Inverter:
hour: int,
allow_battery_grid_export: bool = False,
) -> tuple[float, float, float, float]:
"""Process one slot using probabilistic direct PV-to-load overlap.
``generation`` and ``consumption`` are interval energies. The load
probability table is evaluated in watts and yields the expected direct
PV-to-load power. The remaining load and PV surplus are then handled
independently, because both can occur during different sub-intervals of
the same hourly or 15-minute slot.
"""
losses = 0.0
grid_export = 0.0
grid_import = 0.0
self_consumption = 0.0
generation = max(float(generation), 0.0)
consumption = max(float(consumption), 0.0)
# Cache inverter DC→AC efficiency for discharge path
dc_to_ac_eff = self.dc_to_ac_efficiency
if generation >= consumption:
if consumption > self.max_power_wh:
# If consumption exceeds maximum inverter power
losses += generation - self.max_power_wh
remaining_power = self.max_power_wh - consumption
grid_import = -remaining_power # Negative indicates feeding into the grid
self_consumption = self.max_power_wh
else:
# Calculate scr using cached results per energy management/optimization run.
# The interpolator expects power levels [W]; consumption/generation are
# energy per slot [Wh], so convert via the slot duration (identical at
# the hourly default, ×4 on the 15-minute grid).
scr = self.self_consumption_predictor.calculate_self_consumption(
consumption / self.slot_duration_h, generation / self.slot_duration_h
# Convert interval energy [Wh] to mean power [W] for the probability
# lookup, then convert its expected direct power back to slot energy.
if generation > 0.0 and consumption > 0.0:
expected_direct_power_w = (
self.self_consumption_predictor.calculate_expected_direct_consumption(
consumption / self.slot_duration_h,
generation / self.slot_duration_h,
)
# Remaining power after consumption
remaining_power = (generation - consumption) * scr # EVQ
# Remaining load Self Consumption not perfect
remaining_load_evq = (generation - consumption) * (1.0 - scr)
from_battery_dc = 0.0
if remaining_load_evq > 0:
# Akku muss den Restverbrauch decken
if self.battery:
# Request more DC from battery to account for DC→AC conversion loss
dc_request = remaining_load_evq / dc_to_ac_eff
from_battery_dc, discharge_losses = self.battery.discharge_energy(
dc_request, hour
)
# Convert DC output to AC
from_battery_ac = from_battery_dc * dc_to_ac_eff
inverter_discharge_losses = from_battery_dc - from_battery_ac
remaining_load_evq -= from_battery_ac
losses += discharge_losses + inverter_discharge_losses
else:
from_battery_ac = 0.0
# Wenn der Akku den Restverbrauch nicht vollständig decken kann, wird der Rest ins Netz gezogen
if remaining_load_evq > 0:
grid_import += remaining_load_evq
remaining_load_evq = 0
else:
from_battery_ac = 0.0
if remaining_power > 0:
# Load battery with excess energy (DC path, no inverter conversion needed)
charge_losses = 0.0
if self.battery:
charged_energie, charge_losses = self.battery.charge_energy(
remaining_power, hour
)
remaining_surplus = remaining_power - (charged_energie + charge_losses)
else:
remaining_surplus = remaining_power
# Feed-in to the grid based on remaining capacity
if remaining_surplus > self.max_power_wh - consumption:
grid_export = self.max_power_wh - consumption
losses += remaining_surplus - grid_export
else:
grid_export = remaining_surplus
losses += charge_losses
self_consumption = (
consumption + from_battery_ac
) # Self-consumption is equal to the load
if allow_battery_grid_export and self.battery:
export_capacity = max(self.max_power_wh - consumption - grid_export, 0.0)
remaining_battery_ac = (
self.battery.remaining_discharge_energy_wh(hour) * dc_to_ac_eff
)
export_capacity = min(export_capacity, remaining_battery_ac)
battery_export_ac, battery_export_losses = self._discharge_battery_to_ac(
export_capacity, hour
)
grid_export += battery_export_ac
losses += battery_export_losses
)
direct_pv_energy = expected_direct_power_w * self.slot_duration_h
else:
# Case 2: Insufficient generation, cover shortfall
shortfall = consumption - generation
available_ac_power = max(self.max_power_wh - generation, 0)
direct_pv_energy = 0.0
# Discharge battery to cover shortfall, if possible
if self.battery:
# Need shortfall in AC, request more DC from battery for DC→AC conversion
ac_needed = min(shortfall, available_ac_power)
dc_request = ac_needed / dc_to_ac_eff
battery_discharge_dc, discharge_losses = self.battery.discharge_energy(
dc_request, hour
)
# Convert DC output to AC
battery_discharge_ac = battery_discharge_dc * dc_to_ac_eff
inverter_discharge_losses = battery_discharge_dc - battery_discharge_ac
losses += discharge_losses + inverter_discharge_losses
else:
battery_discharge_ac = 0
# Direct PV is bounded by both input energies and by the AC energy the
# inverter can move during this slot.
direct_pv_energy = min(
max(direct_pv_energy, 0.0),
generation,
consumption,
self.max_power_wh,
)
remaining_load = max(consumption - direct_pv_energy, 0.0)
pv_surplus = max(generation - direct_pv_energy, 0.0)
remaining_inverter_ac_capacity = max(self.max_power_wh - direct_pv_energy, 0.0)
# Draw remaining required power from the grid (discharge_losses are already subtracted in the battery)
grid_import = shortfall - battery_discharge_ac
self_consumption = generation + battery_discharge_ac
# Load gaps and PV surplus may both occur within the same coarse slot.
# Cover the load gap first; this preserves the existing chronological
# approximation and can create headroom for later PV charging.
battery_discharge_ac = 0.0
if remaining_load > 0.0 and self.battery and remaining_inverter_ac_capacity > 0.0:
requested_ac_wh = min(remaining_load, remaining_inverter_ac_capacity)
battery_discharge_ac, battery_discharge_losses = self._discharge_battery_to_ac(
requested_ac_wh, hour
)
remaining_load = max(remaining_load - battery_discharge_ac, 0.0)
remaining_inverter_ac_capacity = max(
remaining_inverter_ac_capacity - battery_discharge_ac, 0.0
)
losses += battery_discharge_losses
if allow_battery_grid_export and self.battery and grid_import <= 0.0:
export_capacity = max(self.max_power_wh - consumption, 0.0)
remaining_battery_ac = (
self.battery.remaining_discharge_energy_wh(hour) * dc_to_ac_eff
)
export_capacity = min(export_capacity, remaining_battery_ac)
battery_export_ac, battery_export_losses = self._discharge_battery_to_ac(
export_capacity, hour
)
grid_export += battery_export_ac
losses += battery_export_losses
grid_import = remaining_load
# Charge from the probabilistic PV surplus on the DC path. Stored energy
# plus charge losses equals the PV energy accepted by the battery.
remaining_surplus = pv_surplus
if remaining_surplus > 0.0 and self.battery:
charged_energy, charge_losses = self.battery.charge_energy(remaining_surplus, hour)
remaining_surplus = max(remaining_surplus - charged_energy - charge_losses, 0.0)
losses += charge_losses
pv_grid_export = min(remaining_surplus, remaining_inverter_ac_capacity)
grid_export += pv_grid_export
remaining_inverter_ac_capacity = max(remaining_inverter_ac_capacity - pv_grid_export, 0.0)
# PV which can neither charge the battery nor pass through the inverter
# is curtailed and reported as a loss.
losses += max(remaining_surplus - pv_grid_export, 0.0)
if allow_battery_grid_export and self.battery and remaining_inverter_ac_capacity > 0.0:
remaining_battery_ac = (
self.battery.remaining_discharge_energy_wh(hour) * self.dc_to_ac_efficiency
)
export_capacity = min(remaining_inverter_ac_capacity, remaining_battery_ac)
battery_export_ac, battery_export_losses = self._discharge_battery_to_ac(
export_capacity, hour
)
grid_export += battery_export_ac
losses += battery_export_losses
self_consumption = direct_pv_energy + battery_discharge_ac
return grid_export, grid_import, losses, self_consumption
@@ -450,7 +450,18 @@ class GeneticSimulation(PydanticBaseModel):
feed_in_tariff_per_hour[hour_idx] = hourly_feed_in_tariff
# Financial calculations
costs_per_hour[hour_idx] = energy_consumption_grid_actual * hourly_electricity_price
grid_cost = energy_consumption_grid_actual * hourly_electricity_price
# LCOS is charged exactly once on battery-delivered DC energy. It is
# not charged on input energy, internal discharge losses, or the
# downstream DC-to-AC inverter loss.
battery_lcos_cost = 0.0
if battery_fast:
battery_lcos_cost = (
battery_fast.discharged_energy_wh(hour)
* battery_fast.levelized_cost_of_storage_kwh
/ 1000.0
)
costs_per_hour[hour_idx] = grid_cost + battery_lcos_cost
revenue_per_hour[hour_idx] = energy_feedin_grid_actual * hourly_feed_in_tariff
total_cost = np.nansum(costs_per_hour)
@@ -1244,8 +1255,14 @@ class GeneticOptimization(OptimizationBase):
if charge_price <= 0:
continue
# Price that a future discharge hour must reach to break even
break_even_price = charge_price / round_trip_eff
# Price that a future AC discharge hour must reach to break
# even. LCOS is defined per DC Wh delivered by the battery;
# dividing it by DC-to-AC efficiency converts it to the
# corresponding cost per useful/exported AC Wh.
lcos_per_wh_dc = getattr(bat, "levelized_cost_of_storage_kwh", 0.0) / 1000.0
break_even_price = (
charge_price / round_trip_eff + lcos_per_wh_dc / inv.dc_to_ac_efficiency
)
best_uncovered_price = best_prices[hour]
@@ -98,6 +98,17 @@ class BaseBatteryParameters(DeviceParameters):
class SolarPanelBatteryParameters(BaseBatteryParameters):
"""PV battery device simulation configuration."""
levelized_cost_of_storage_kwh: float = Field(
default=0.0,
ge=0.0,
json_schema_extra={
"description": (
"Levelized cost of storage applied once to each kWh delivered "
"by the battery [EUR/kWh]."
),
"examples": [0.12],
},
)
max_charge_power_w: Optional[float] = max_charging_power_field()
@@ -54,7 +54,10 @@ class GeneticEnergyManagementParameters(GeneticParametersBaseModel):
)
preis_euro_pro_wh_akku: float = Field(
json_schema_extra={
"description": "A float representing the cost of battery energy per watt-hour."
"description": (
"Terminal value of usable battery energy remaining at the end of the "
"optimization horizon [EUR/Wh]. This is not the battery LCOS."
)
}
)
gesamtlast: list[float] = Field(
@@ -416,7 +419,6 @@ class GeneticOptimizationParameters(
cls.config.devices.max_batteries = 1
if cls.config.devices.max_batteries == 0:
battery_params = None
battery_lcos_kwh = 0
else:
if cls.config.devices.batteries is None:
logger.info("No battery device data available - defaulting to demo data.")
@@ -428,6 +430,9 @@ class GeneticOptimizationParameters(
capacity_wh=battery_config.capacity_wh,
charging_efficiency=battery_config.charging_efficiency,
discharging_efficiency=battery_config.discharging_efficiency,
levelized_cost_of_storage_kwh=(
battery_config.levelized_cost_of_storage_kwh
),
max_charge_power_w=battery_config.max_charge_power_w,
min_soc_percentage=battery_config.min_soc_percentage,
max_soc_percentage=battery_config.max_soc_percentage,
@@ -441,14 +446,6 @@ class GeneticOptimizationParameters(
cls.config.devices.batteries = [{"device_id": "battery1", "capacity_wh": 8000}]
# Retry
continue
# Levelized cost of ownership
if battery_config.levelized_cost_of_storage_kwh is None:
logger.info(
"No battery device LCOS data available - defaulting to 0 €/kWh. Parameter preparation attempt {}.",
attempt,
)
battery_config.levelized_cost_of_storage_kwh = 0
battery_lcos_kwh = battery_config.levelized_cost_of_storage_kwh
# Initial SOC
try:
initial_soc_factor = cls.measurement.key_to_value(
@@ -654,7 +651,9 @@ class GeneticOptimizationParameters(
strompreis_euro_pro_wh=elecprice_marketprice_wh,
einspeiseverguetung_euro_pro_wh=feed_in_tariff_wh,
gesamtlast=loadforecast_power_w,
preis_euro_pro_wh_akku=battery_lcos_kwh / 1000,
preis_euro_pro_wh_akku=(
cls.config.optimization.terminal_value_euro_per_kwh / 1000
),
),
temperature_forecast=weather_temp_air,
pv_akku=battery_params,
@@ -101,6 +101,18 @@ class OptimizationCommonSettings(SettingsBaseModel):
},
)
terminal_value_euro_per_kwh: float = Field(
default=0.0,
json_schema_extra={
"description": (
"Value assigned to usable battery energy remaining at the end of the "
"optimization horizon [EUR/kWh]. This terminal value is independent "
"of the battery LCOS. Defaults to 0 EUR/kWh."
),
"examples": [0.0, 0.20],
},
)
genetic: GeneticCommonSettings = Field(
default_factory=GeneticCommonSettings,
json_schema_extra={
+70 -1
View File
@@ -17,8 +17,32 @@ class SelfConsumptionProbabilityInterpolator:
self.interpolator: RegularGridInterpolator = pickle.load(file) # noqa: S301
self.load_power_min_w = float(self.interpolator.grid[0][0])
self.load_power_max_w = float(self.interpolator.grid[0][-1])
self.minute_load_levels_w = np.asarray(self.interpolator.grid[1], dtype=float)
self.minute_load_max_w = float(self.interpolator.grid[1][-1])
def _load_distribution(self, mean_load_power_w: float) -> tuple[np.ndarray, np.ndarray]:
"""Return the conditional minute-load distribution for a mean load.
The table stores one probability mass for each 50 W minute-load bin.
Linear interpolation between its mean-load rows can introduce very small
numerical deviations, so negative masses are removed and the result is
normalized explicitly.
"""
bounded_mean_load_w = float(
np.clip(mean_load_power_w, self.load_power_min_w, self.load_power_max_w)
)
points = np.column_stack(
(
np.full(self.minute_load_levels_w.shape, bounded_mean_load_w),
self.minute_load_levels_w,
)
)
probabilities = np.maximum(np.asarray(self.interpolator(points), dtype=float), 0.0)
probability_sum = float(probabilities.sum())
if probability_sum <= 0.0:
return self.minute_load_levels_w, probabilities
return self.minute_load_levels_w, probabilities / probability_sum
def _generate_points(
self, mean_load_power_w: float, pv_power_w: float
) -> tuple[np.ndarray, np.ndarray]:
@@ -39,7 +63,12 @@ class SelfConsumptionProbabilityInterpolator:
@cache_energy_management
def calculate_self_consumption(self, mean_load_power_w: float, pv_power_w: float) -> float:
"""Calculate the PV self-consumption rate using RegularGridInterpolator.
"""Return the legacy cumulative minute-load probability.
This method is retained for API compatibility. Its result is the
probability that the minute load is no greater than ``pv_power_w``;
it is not an energy self-consumption ratio. New energy-flow code must
use :meth:`calculate_expected_direct_consumption`.
The results are cached until the start of the next energy management run/ optimization.
@@ -54,6 +83,46 @@ class SelfConsumptionProbabilityInterpolator:
probabilities = self.interpolator(points)
return float(np.clip(probabilities.sum(), 0.0, 1.0))
@cache_energy_management
def calculate_expected_direct_consumption(
self, mean_load_power_w: float, pv_power_w: float
) -> float:
"""Calculate expected direct PV-to-load power in watts.
For conditional minute-load probabilities ``p_i`` and load-bin powers
``L_i``, the expected direct consumption is
``sum(p_i * min(L_i, pv_power_w))``.
The tabulated load-bin powers are rescaled to preserve the supplied
forecast mean exactly. This compensates for discretization and the
finite upper table boundary while retaining the distribution shape.
Args:
mean_load_power_w: Mean load power of the forecast interval [W].
pv_power_w: Mean PV power of the forecast interval [W].
Returns:
Expected direct PV-to-load power [W].
"""
mean_load_power_w = max(float(mean_load_power_w), 0.0)
pv_power_w = max(float(pv_power_w), 0.0)
if mean_load_power_w == 0.0 or pv_power_w == 0.0:
return 0.0
load_levels_w, probabilities = self._load_distribution(mean_load_power_w)
modeled_mean_load_w = float(np.dot(probabilities, load_levels_w))
if modeled_mean_load_w <= 0.0:
return 0.0
# Preserve the requested mean load while keeping the conditional shape
# from the probability table.
normalized_load_levels_w = load_levels_w * (mean_load_power_w / modeled_mean_load_w)
expected_direct_power_w = float(
np.dot(probabilities, np.minimum(normalized_load_levels_w, pv_power_w))
)
return float(np.clip(expected_direct_power_w, 0.0, min(mean_load_power_w, pv_power_w)))
# def calculate_self_consumption(self, load_1h_power: float, pv_power: float) -> float:
# """Calculate the PV self-consumption rate using RegularGridInterpolator.