Handling Uncertainty: From Deterministic to Robust

Every lesson so far assumed we know prices perfectly. We don't. Forecasts are wrong. Wind drops unexpectedly. A plant trips. This lesson teaches you three approaches to handling price uncertainty — from simplest to most sophisticated.

The problem with deterministic optimization

Your Lesson 3 optimizer takes a single price vector and finds the optimal schedule for those exact prices. If actual prices differ, your "optimal" schedule may be suboptimal or even loss-making.

# Example: forecast says hour 18 = €90/MWh → optimizer discharges
# Reality: hour 18 = €45/MWh (cloudy day, solar overproduced)
# Result: you sold at €45 instead of waiting for a real peak later

Approach 1: Sensitivity analysis (simplest)

Run your optimizer multiple times with different price scenarios and see how the schedule changes:

"""
sensitivity.py — How sensitive is the schedule to price changes?
"""
import numpy as np

# Base forecast
prices_base = [28, 25, 22, 20, 19, 21, 35, 55,
               62, 58, 52, 48, 45, 44, 47, 53,
               68, 85, 92, 78, 60, 45, 35, 30]

# Generate scenarios: base ± noise
np.random.seed(42)
n_scenarios = 20
scenarios = []
for _ in range(n_scenarios):
    noise = np.random.normal(0, 8, 24)  # ±€8/MWh std dev
    scenario = np.clip(np.array(prices_base) + noise, -50, 200)
    scenarios.append(scenario.tolist())

# Run optimizer for each scenario, collect schedules
# (using solve_battery function from previous lessons)
# results = [solve_battery(s) for s in scenarios]

# Analyze: which hours are ALWAYS charge/discharge regardless of scenario?
# Those are robust decisions. Hours that flip → risky, bid conservatively.

Approach 2: Scenario-based stochastic optimization

Instead of optimizing for one price path, optimize across multiple scenarios simultaneously. The objective becomes expected revenue across all scenarios:

"""
stochastic.py — Two-stage stochastic optimization
"""
import pyomo.environ as pyo
import numpy as np

T = 24
P_max = 10.0
E_max = 20.0
eta = 0.90
eta_sqrt = np.sqrt(eta)
soc_init = 10.0
deg_cost = 5.0

# Price scenarios with probabilities
scenarios = {
    'base':   {'prob': 0.5, 'prices': [28,25,22,20,19,21,35,55,62,58,52,48,45,44,47,53,68,85,92,78,60,45,35,30]},
    'high':   {'prob': 0.25, 'prices': [35,30,28,25,24,28,45,70,80,75,68,62,58,55,60,68,88,110,120,95,75,55,42,35]},
    'low':    {'prob': 0.25, 'prices': [22,20,18,16,15,17,28,42,48,45,40,38,35,34,36,40,52,65,70,58,45,35,28,22]},
}

model = pyo.ConcreteModel()
model.T = pyo.RangeSet(0, T - 1)
model.S = pyo.Set(initialize=scenarios.keys())

# FIRST-STAGE decisions: bid quantities (committed before knowing actual prices)
model.p_charge = pyo.Var(model.T, bounds=(0, P_max))
model.p_discharge = pyo.Var(model.T, bounds=(0, P_max))
model.soc = pyo.Var(model.T, bounds=(0.1*E_max, 0.9*E_max))
model.y = pyo.Var(model.T, within=pyo.Binary)

# Objective: expected revenue across scenarios
def obj_rule(m):
    total = 0
    for s_name, s_data in scenarios.items():
        prob = s_data['prob']
        prices = s_data['prices']
        revenue = sum(prices[t] * (m.p_discharge[t] - m.p_charge[t]) for t in m.T)
        deg = sum(deg_cost * (m.p_charge[t] + m.p_discharge[t]) for t in m.T)
        total += prob * (revenue - deg)
    return total
model.obj = pyo.Objective(rule=obj_rule, sense=pyo.maximize)

# Constraints (same structure as before — applied once since schedule is fixed)
def soc_rule(m, t):
    prev = soc_init if t == 0 else m.soc[t-1]
    return m.soc[t] == prev + eta_sqrt * m.p_charge[t] - m.p_discharge[t] / eta_sqrt
model.soc_bal = pyo.Constraint(model.T, rule=soc_rule)

def nc(m, t): return m.p_charge[t] <= P_max * (1 - m.y[t])
def nd(m, t): return m.p_discharge[t] <= P_max * m.y[t]
model.nc = pyo.Constraint(model.T, rule=nc)
model.nd = pyo.Constraint(model.T, rule=nd)

solver = pyo.SolverFactory('highs')
result = solver.solve(model, tee=False)
print(f"Expected revenue: €{pyo.value(model.obj):.2f}")

Approach 3: Robust optimization (worst-case protection)

Instead of optimizing expected value, guarantee a minimum performance under the worst-case scenario:

"""
robust.py — Maximize worst-case revenue across scenarios
"""
import pyomo.environ as pyo
import numpy as np

# Same setup as above...
model = pyo.ConcreteModel()
model.T = pyo.RangeSet(0, T - 1)

model.p_charge = pyo.Var(model.T, bounds=(0, P_max))
model.p_discharge = pyo.Var(model.T, bounds=(0, P_max))
model.soc = pyo.Var(model.T, bounds=(0.1*E_max, 0.9*E_max))
model.y = pyo.Var(model.T, within=pyo.Binary)

# Auxiliary variable: worst-case revenue
model.worst_case = pyo.Var()

# Objective: maximize the worst case
model.obj = pyo.Objective(expr=model.worst_case, sense=pyo.maximize)

# Worst case must be ≤ revenue in every scenario
def scenario_bound(m, s_name):
    prices = scenarios[s_name]['prices']
    revenue = sum(prices[t] * (m.p_discharge[t] - m.p_charge[t]) for t in m.T)
    deg = sum(deg_cost * (m.p_charge[t] + m.p_discharge[t]) for t in m.T)
    return m.worst_case <= revenue - deg

model.S = pyo.Set(initialize=scenarios.keys())
model.scen_bounds = pyo.Constraint(model.S, rule=scenario_bound)

# ... same SoC and binary constraints as before ...

The robust solution is more conservative: it won't chase high-risk peaks that only appear in one scenario. Instead, it finds schedules that perform decently across all scenarios.

When to use which approach

Check your understanding

What's next

You now have the complete picture: a battery optimizer that can bid into day-ahead, reoptimize intraday, co-optimize with FCR, and handle uncertainty. Lesson 8 brings everything together into a production architecture.

Rechtliches: Impressum