diff --git a/openapi.yaml b/openapi.yaml index 7447c542e..3e416c342 100644 --- a/openapi.yaml +++ b/openapi.yaml @@ -269,6 +269,39 @@ components: maximum: 2 default: 0 description: Charging and discharging priority 0..2 compared to other batteries. 2 = highest priority. + charge_profile: + type: string + description: | + Named charge profile for SoC-dependent charge power limiting + (CC-CV taper). When set, the MILP solver uses iterative constraint + tightening to enforce realistic charge power limits at high SoC. + Available profiles: lifepo4_conservative, lifepo4_moderate. + Individual parameters can be overridden with charge_knee, + charge_k, charge_c_rate_float. + example: lifepo4_conservative + charge_knee: + type: number + minimum: 0 + maximum: 100 + description: | + SoC percentage where charge power tapering begins. + With charge_profile: overrides the profile default. + Without charge_profile: requires charge_k to build a custom profile. + example: 65 + charge_k: + type: number + minimum: 0 + description: | + Exponential decay constant for the CC-CV taper curve. + Higher values mean steeper power reduction above the knee. + example: 0.052 + charge_c_rate_float: + type: number + minimum: 0 + description: | + Minimum C-rate at the tail end of charging (float/balancing phase). + Default: 0.01 (C/100). + example: 0.01 TimeSeries: type: object diff --git a/src/optimizer/app.py b/src/optimizer/app.py index 0efdd660a..16180b513 100644 --- a/src/optimizer/app.py +++ b/src/optimizer/app.py @@ -5,6 +5,7 @@ from flask_restx import Api, Resource, fields from werkzeug.exceptions import BadRequest +from .charging_profiles import get_profile from .optimizer import BatteryConfig, GridConfig, OptimizationStrategy, Optimizer, TimeSeriesData app = Flask(__name__) @@ -80,7 +81,11 @@ def handle_validation_error(error): 'c_max': fields.Float(required=True, description='Maximum charge power (W)'), 'd_max': fields.Float(required=True, description='Maximum discharge power (W)'), 'p_a': fields.Float(required=True, description='Monetary value per Wh at end of the optimization horizon'), - 'c_priority': fields.Integer(required=False, description='Charging and discharging priority compared to other batteries. 2 = highest priority.') + 'c_priority': fields.Integer(required=False, description='Charging and discharging priority compared to other batteries. 2 = highest priority.'), + 'charge_profile': fields.String(required=False, description='Named CC-CV taper profile (e.g. lifepo4_conservative).'), + 'charge_knee': fields.Float(required=False, description='SoC% where taper begins. Overrides profile default.'), + 'charge_k': fields.Float(required=False, description='Exponential decay constant for the CC-CV taper curve.'), + 'charge_c_rate_float': fields.Float(required=False, description='Minimum C-rate at float / tail end of charging (default 0.01).'), }) time_series_model = api.model('TimeSeries', { @@ -157,20 +162,50 @@ def post(self): # Parse battery configurations batteries = [] for bat_data in data['batteries']: + # Build optional CC-CV taper profile + profile = None + profile_name = bat_data.get('charge_profile') + charge_knee = bat_data.get('charge_knee') + charge_k = bat_data.get('charge_k') + charge_c_rate_float = bat_data.get('charge_c_rate_float') + + s_capacity = bat_data.get('s_capacity', bat_data['s_max']) + c_max = bat_data['c_max'] + + if profile_name: + # Named profile with optional overrides + overrides = {} + if charge_knee is not None: + overrides['knee'] = charge_knee + if charge_k is not None: + overrides['k'] = charge_k + if charge_c_rate_float is not None: + overrides['c_rate_float'] = charge_c_rate_float + profile = get_profile(profile_name, **overrides) + elif charge_knee is not None and charge_k is not None: + # Custom profile: derive c_rate_max from c_max / s_capacity + profile = { + 'c_rate_max': c_max / s_capacity if s_capacity > 0 else 0.25, + 'knee': charge_knee, + 'k': charge_k, + 'c_rate_float': charge_c_rate_float if charge_c_rate_float is not None else 0.01, + } + batteries.append(BatteryConfig( charge_from_grid=bat_data.get('charge_from_grid', False), discharge_to_grid=bat_data.get('discharge_to_grid', False), - s_capacity=bat_data.get('s_capacity', bat_data['s_max']), + s_capacity=s_capacity, s_min=bat_data['s_min'], s_max=bat_data['s_max'], s_initial=bat_data['s_initial'], p_demand=bat_data.get('p_demand'), s_goal=bat_data.get('s_goal'), c_min=bat_data['c_min'], - c_max=bat_data['c_max'], + c_max=c_max, d_max=bat_data['d_max'], p_a=bat_data['p_a'], c_priority=bat_data.get('c_priority', 0), + charge_profile=profile, )) # Parse time series data diff --git a/src/optimizer/charging_profiles.py b/src/optimizer/charging_profiles.py new file mode 100644 index 000000000..960ea4c6b --- /dev/null +++ b/src/optimizer/charging_profiles.py @@ -0,0 +1,132 @@ +""" +Battery charging profiles by cell chemistry. + +Each profile defines how the maximum charge power decreases as a function +of State of Charge (SoC). This models the real-world CC-CV charging behavior +where the Battery Management System (BMS) reduces the charge current at +higher SoC to protect cells and minimize resistive losses. + +Usage: + from charging_profiles import get_profile, max_charge_power + + profile = get_profile("lifepo4_conservative") + watts = max_charge_power(profile, capacity_wh=31000, soc_percent=75) + +To find values for your own system: + - Set c_rate_max to the highest C-rate your battery sustains + without the BMS reducing the charge current + - Set knee to the SoC percentage where the BMS starts + reducing the charge current below c_rate_max + - Leave k and c_rate_float at their defaults unless you have + your own measurements (fit with Wolfram Alpha: + "exponential fit {SoC1,I1}, {SoC2,I2}, ...") + +General advice: always stay slightly below your measured values. +If your BMS starts tapering at 66%, set knee to 65%. If your +inverter sustains 160A, set c_rate_max to match C/4 or lower. +This avoids triggering BMS intervention, reduces resistive losses, +and extends cell life. +""" + +import math + +# --------------------------------------------------------------------------- +# Chemistry profiles +# --------------------------------------------------------------------------- + +CHEMISTRY_PROFILES = { + "lifepo4_conservative": { + # Charging profile for LiFePO4 cells optimized for longevity and + # minimal resistive losses. Derived from real-world measurements + # on a mixed environment of aged LiYFePO4 cells (200Ah + 90Ah) + # and newer LiFePO4 cells (314Ah) with Victron inverters. + # + # The charging curve has two phases: + # 1. Below 'knee': constant power at 'c_rate_max' rate + # 2. Above 'knee': power drops exponentially toward 'c_rate_float' + # + # Formula above knee: + # P = max(c_rate_float, c_rate_max * e^(-k * (SoC - knee))) + + "c_rate_max": 0.25, + # Maximum charge rate below knee, as fraction of capacity. + # C/4 means a 31 kWh battery charges at up to 7750 W. + + "knee": 65, + # SoC percentage where charge current begins to taper. + + "k": 0.052, + # Exponential decay constant (R^2 = 0.996 from BMS measurements). + + "c_rate_float": 0.01, + # Minimum charge rate (maintenance/balancing phase). + # C/100 means a 31 kWh battery trickle-charges at 310 W. + }, + + "lifepo4_moderate": { + # Moderate profile for newer LiFePO4 cells in good condition. + + "c_rate_max": 0.33, # C/3 + "knee": 75, + "k": 0.045, + "c_rate_float": 0.02, # C/50 + }, +} + + +# --------------------------------------------------------------------------- +# Public API +# --------------------------------------------------------------------------- + +def list_profiles(): + """Return list of available profile names.""" + return list(CHEMISTRY_PROFILES.keys()) + + +def get_profile(name, **overrides): + """ + Get a charging profile by name with optional overrides. + + Args: + name: Profile name (e.g. "lifepo4_conservative") + **overrides: Override individual values (e.g. knee=70, c_rate_max=0.2) + + Returns: + dict with keys: c_rate_max, knee, k, c_rate_float + """ + if name not in CHEMISTRY_PROFILES: + available = ", ".join(CHEMISTRY_PROFILES.keys()) + raise ValueError(f"Unknown profile '{name}'. Available: {available}") + + profile = dict(CHEMISTRY_PROFILES[name]) + for key, value in overrides.items(): + if key not in profile: + raise ValueError(f"Unknown parameter '{key}'. Valid: {', '.join(profile.keys())}") + profile[key] = value + + return profile + + +def max_charge_power(profile, capacity_wh, soc_percent): + """ + Calculate maximum charge power in watts for a given SoC. + + Args: + profile: Charging profile dict from get_profile() + capacity_wh: Battery capacity in Wh + soc_percent: Current state of charge (0-100) + + Returns: + Maximum charge power in watts + """ + c_rate_max = profile["c_rate_max"] + knee = profile["knee"] + k = profile["k"] + c_rate_float = profile["c_rate_float"] + + if soc_percent <= knee: + c_rate = c_rate_max + else: + c_rate = max(c_rate_float, c_rate_max * math.exp(-k * (soc_percent - knee))) + + return capacity_wh * c_rate diff --git a/src/optimizer/optimizer.py b/src/optimizer/optimizer.py index b878ec8e6..f8007baf0 100644 --- a/src/optimizer/optimizer.py +++ b/src/optimizer/optimizer.py @@ -1,3 +1,4 @@ +import math from dataclasses import dataclass from tempfile import TemporaryDirectory from typing import Dict, List, Optional @@ -5,6 +6,7 @@ import numpy as np import pulp +from .charging_profiles import max_charge_power from .settings import OptimizerSettings @@ -36,6 +38,7 @@ class BatteryConfig: p_demand: Optional[List[float]] = None # Minimum charge demand (Wh) s_goal: Optional[List[float]] = None # Goal state of charge (Wh) c_priority: int = 0 + charge_profile: Optional[dict] = None # CC-CV taper profile (charging_profiles format) @dataclass @@ -115,6 +118,7 @@ def create_model(self): self._setup_target_function() self._add_energy_balance_constraints() self._add_battery_constraints() + self._add_charge_profile_constraints() def _setup_variables(self): """ @@ -491,10 +495,80 @@ def _add_battery_constraints(self): # Charge constraint self.problem += self.variables['c'][i][t] <= self.M * (1 - self.variables['z_cd'][i][t]) + def _add_charge_profile_constraints(self): + """Add piecewise-linear upper bounds on charge power from CC-CV taper profiles. + + For each battery with a charge_profile, the taper curve P_max(SoC) is + a convex decreasing function above the knee. We approximate it with + tangent lines at 4 breakpoints. Each tangent is a valid linear upper + bound (convexity guarantees the tangent never exceeds the curve), so + the MILP enforces realistic charge power limits in a single solve. + + The constraint for slot t uses s[i][t-1] (SoC at slot start). For + t=0, s_initial is a constant so the bound is applied directly to + the variable's upper bound. + """ + for i, bat in enumerate(self.batteries): + if bat.charge_profile is None: + continue + + profile = bat.charge_profile + cap = bat.s_capacity + if cap <= 0: + continue + + knee = profile['knee'] + k = profile['k'] + c_rate_max = profile['c_rate_max'] + + # 4 breakpoints from knee to 100% SoC + n_points = 4 + breakpoints = [] + for j in range(n_points): + soc_pct = knee + (100.0 - knee) * j / (n_points - 1) + s_wh = soc_pct / 100.0 * cap + p_w = max_charge_power(profile, cap, soc_pct) + # derivative of P_max w.r.t. SoC in Wh: + # dP/ds = dP/d(soc%) * d(soc%)/ds + # = (-k * c_rate_max * cap * e^(-k*(soc%-knee))) * (100/cap) + # = -k * c_rate_max * 100 * e^(-k*(soc%-knee)) + dp_ds = -k * c_rate_max * 100.0 * math.exp(-k * (soc_pct - knee)) + breakpoints.append((s_wh, p_w, dp_ds)) + + for s_b, p_b, dp_ds in breakpoints: + for t in self.time_steps: + dt_h = self.time_series.dt[t] / 3600.0 + + if t == 0: + # s_initial is a constant: compute a fixed upper bound + p_limit = p_b + dp_ds * (bat.s_initial - s_b) + e_limit = p_limit * dt_h + e_var_ub = bat.c_max * dt_h + # only add if tighter than the existing variable bound + if e_limit < e_var_ub - 1.0: + self.problem += ( + self.variables['c'][i][t] <= e_limit, + f"taper_{i}_{t}_bp{breakpoints.index((s_b, p_b, dp_ds))}", + ) + else: + # linear constraint: c[i][t] <= (p_b - dp_ds*s_b)*dt_h + dp_ds*dt_h * s[i][t-1] + s_prev = self.variables['s'][i][t - 1] + intercept = (p_b - dp_ds * s_b) * dt_h + slope = dp_ds * dt_h + self.problem += ( + self.variables['c'][i][t] <= intercept + slope * s_prev, + f"taper_{i}_{t}_bp{breakpoints.index((s_b, p_b, dp_ds))}", + ) + def solve(self) -> Dict: """ Creates the MILP model if none exists and solves the optimization problem. - Returns a dictionary with the optimization results + + When batteries have a charge_profile (CC-CV taper), piecewise-linear + upper bounds on charge power are part of the model, so a single solve + enforces realistic SoC-dependent charge limits. + + Returns a dictionary with the optimization results. """ if self.problem is None: