#!/usr/bin/env python3 """ Fáze 0 – ekonomický report: skutečný provoz vs. oracle LP (perfect hindsight). Pro každý pražský den v zadaném rozsahu spočítá: 1. ACTUAL – skutečný cashflow z ems.audit_interval (import × eff. buy − export × eff. sell), 2. ORACLE – dolní mez nákladů: malý ČISTÝ MILP (jen reálné peníze: nákup − prodej − degradace, žádné heuristické penalty) nad SKUTEČNOU PV výrobou, SKUTEČNOU spotřebou a skutečnými cenami dne (perfect foresight), 3. GAP – actual − oracle, férově očištěno o rozdíl koncového SoC (koncová energie oceněna průměrnou denní nákupní cenou). GAP = forecast error + neefektivita dispatche. Oracle nelze v reálu dosáhnout (zná budoucnost), ale trend GAPu je objektivní metrika „neekonomického provozu“ a regresní metr pro Fázi 2/3 (čistý plánovač). Oracle LP je zároveň zárodek čistého jádra solveru. Model oracle (15min sloty, Wh): pv_a_used + pv_b + gi + bd = load + bc + ge (bilance sběrnice) pv_a_used ≤ pv_a_actual (curtailment jen pole A) soc[t] = soc[t-1] + bc·η_c − bd/η_d (SoC dynamika) min_soc ≤ soc ≤ soc_max; výkonové stropy baterie i sítě sell < 0 ∧ block_export_on_negative_sell → ge = 0 (hard pravidlo KV1) binárka: zákaz současného importu a exportu objective: Σ gi·buy − ge·sell + ½(bc+bd)·degradace − soc[T]·avg_buy Zjednodušení (dokumentovaná): spotřeba je fixní (EV/TČ se nepřesouvá), zelený bonus PV B vyloučen z obou stran, konfigurace baterie = aktuální stav DB. Použití: EMS_DB_DSN=postgresql://ems_user:***@10.200.200.1:5432/ems \ python3 scripts/harness/economics_report.py --site-code home-01 --from 2026-05-12 --to 2026-06-09 """ from __future__ import annotations import argparse import asyncio import json import os import sys from dataclasses import dataclass from datetime import datetime, timedelta from zoneinfo import ZoneInfo import asyncpg import pulp PRAGUE = ZoneInfo("Europe/Prague") INTERVAL_H = 0.25 SLOT_WH_TO_W = 4 # Wh za 15 min → W def _build_dsn(args: argparse.Namespace) -> str: if args.dsn: return args.dsn env_dsn = os.environ.get("EMS_DB_DSN") if env_dsn: return env_dsn host = os.environ.get("DB_HOST", "127.0.0.1") port = os.environ.get("DB_PORT", "5432") user = os.environ.get("DB_USER", "ems_user") password = os.environ.get("DB_PASSWORD", "") name = os.environ.get("DB_NAME", "ems") return f"postgresql://{user}:{password}@{host}:{port}/{name}" @dataclass class DaySlot: interval_start: datetime buy: float # Kč/kWh efektivní nákup sell: float # Kč/kWh efektivní prodej pv_a_wh: float # skutečná výroba pole A (curtailable) pv_b_wh: float # skutečná výroba pole B (fixní) load_wh: float # skutečná spotřeba (vč. EV/TČ) grid_import_wh: float grid_export_wh: float soc_pct: float | None @dataclass class BatteryParams: usable_wh: float min_soc_wh: float soc_max_wh: float charge_eff: float discharge_eff: float max_charge_w: float max_discharge_w: float degradation_czk_kwh: float async def _load_battery_and_grid(conn: asyncpg.Connection, site_id: int) -> tuple[BatteryParams, dict]: raw = await conn.fetchval("select ems.fn_planning_site_context($1::int)", site_id) ctx = raw if isinstance(raw, dict) else json.loads(raw) b = ctx["battery"] bat = BatteryParams( usable_wh=float(b["usable_capacity_wh"]), min_soc_wh=float(b["min_soc_wh"]), soc_max_wh=float(b.get("planner_soc_max_wh", b["soc_max_wh"])), charge_eff=float(b["charge_efficiency"]), discharge_eff=float(b["discharge_efficiency"]), max_charge_w=float(b["max_charge_power_w"]), max_discharge_w=float(b["max_discharge_power_w"]), degradation_czk_kwh=float(b["degradation_cost_czk_kwh"]), ) g = ctx["grid"] grid = { "max_import_w": float(g["max_import_power_w"]), "max_export_w": float(g["max_export_power_w"]), "block_export_on_negative_sell": bool(g.get("block_export_on_negative_sell") or False), } return bat, grid async def _load_day( conn: asyncpg.Connection, site_id: int, day_start: datetime ) -> list[DaySlot]: day_end = day_start + timedelta(days=1) rows = await conn.fetch( """ select a.interval_start, p.effective_buy_price_czk_kwh as buy, p.effective_sell_price_czk_kwh as sell, coalesce(a.actual_pv_production_wh, 0) as pv_wh, coalesce(a.pv_b_production_wh, 0) as pv_b_wh, coalesce(a.actual_load_consumption_wh, 0) as load_wh, coalesce(a.actual_grid_import_wh, 0) as gi_wh, coalesce(a.actual_grid_export_wh, 0) as ge_wh, a.actual_battery_soc_pct as soc_pct from ems.audit_interval a join ems.vw_site_effective_price p on p.site_id = a.site_id and p.interval_start = a.interval_start where a.site_id = $1 and a.interval_start >= $2 and a.interval_start < $3 order by a.interval_start """, site_id, day_start, day_end, ) return [ DaySlot( interval_start=r["interval_start"], buy=float(r["buy"]), sell=float(r["sell"]), pv_a_wh=max(0.0, float(r["pv_wh"]) - float(r["pv_b_wh"])), pv_b_wh=float(r["pv_b_wh"]), load_wh=float(r["load_wh"]), grid_import_wh=float(r["gi_wh"]), grid_export_wh=float(r["ge_wh"]), soc_pct=float(r["soc_pct"]) if r["soc_pct"] is not None else None, ) for r in rows ] def _actual_cashflow_czk(slots: list[DaySlot]) -> float: return sum( s.grid_import_wh / 1000.0 * s.buy - s.grid_export_wh / 1000.0 * s.sell for s in slots ) def solve_oracle( slots: list[DaySlot], bat: BatteryParams, grid: dict, soc_start_wh: float, terminal_value_czk_kwh: float, ) -> tuple[float, float]: """Vrátí (cash_czk, soc_end_wh) optimálního dispatche s perfektní znalostí dne.""" n = len(slots) prob = pulp.LpProblem("oracle_day", pulp.LpMinimize) max_chg_wh = bat.max_charge_w * INTERVAL_H max_dis_wh = bat.max_discharge_w * INTERVAL_H max_gi_wh = grid["max_import_w"] * INTERVAL_H max_ge_wh = grid["max_export_w"] * INTERVAL_H gi = [pulp.LpVariable(f"gi_{t}", 0, max_gi_wh) for t in range(n)] ge = [pulp.LpVariable(f"ge_{t}", 0, max_ge_wh) for t in range(n)] bc = [pulp.LpVariable(f"bc_{t}", 0, max_chg_wh) for t in range(n)] bd = [pulp.LpVariable(f"bd_{t}", 0, max_dis_wh) for t in range(n)] pa = [pulp.LpVariable(f"pa_{t}", 0, slots[t].pv_a_wh) for t in range(n)] soc = [pulp.LpVariable(f"soc_{t}", bat.min_soc_wh, bat.soc_max_wh) for t in range(n)] z_imp = [pulp.LpVariable(f"zi_{t}", cat="Binary") for t in range(n)] for t in range(n): s = slots[t] # bilance sběrnice (Wh ve slotu) prob += pa[t] + s.pv_b_wh + gi[t] + bd[t] == s.load_wh + bc[t] + ge[t] # SoC dynamika prev = soc_start_wh if t == 0 else soc[t - 1] prob += soc[t] == prev + bc[t] * bat.charge_eff - bd[t] / bat.discharge_eff # zákaz současného importu a exportu prob += gi[t] <= max_gi_wh * z_imp[t] prob += ge[t] <= max_ge_wh * (1 - z_imp[t]) # tvrdé pravidlo: záporný sell → žádný export (kde konfigurováno) if s.sell < 0 and grid["block_export_on_negative_sell"]: prob += ge[t] == 0 cash = pulp.lpSum( gi[t] / 1000.0 * slots[t].buy - ge[t] / 1000.0 * slots[t].sell for t in range(n) ) degradation = pulp.lpSum( 0.5 * (bc[t] + bd[t]) / 1000.0 * bat.degradation_czk_kwh for t in range(n) ) terminal = soc[n - 1] / 1000.0 * terminal_value_czk_kwh prob += cash + degradation - terminal solver = pulp.HiGHS_CMD(msg=False) if pulp.HiGHS_CMD().available() else pulp.PULP_CBC_CMD(msg=False) prob.solve(solver) if pulp.LpStatus[prob.status] != "Optimal": raise RuntimeError(f"Oracle LP není Optimal: {pulp.LpStatus[prob.status]}") cash_val = sum( gi[t].value() / 1000.0 * slots[t].buy - ge[t].value() / 1000.0 * slots[t].sell for t in range(n) ) return cash_val, float(soc[n - 1].value()) async def run_report(args: argparse.Namespace) -> None: dsn = _build_dsn(args) conn = await asyncpg.connect(dsn) try: site_row = await conn.fetchrow("select id from ems.site where code = $1", args.site_code) if site_row is None: raise SystemExit(f"Site code '{args.site_code}' nenalezen") site_id = int(site_row["id"]) bat, grid = await _load_battery_and_grid(conn, site_id) d_from = datetime.strptime(getattr(args, "from"), "%Y-%m-%d").replace(tzinfo=PRAGUE) d_to = datetime.strptime(args.to, "%Y-%m-%d").replace(tzinfo=PRAGUE) print(f"# Ekonomický report — {args.site_code} (site_id={site_id})") print(f"# Okno: {getattr(args, 'from')} … {args.to} (Prague dny), baterie {bat.usable_wh/1000:.1f} kWh") print() header = ( f"{'den':<11} {'actual':>9} {'oracle':>9} {'gap':>8} {'gap%':>6} " f"{'soc0%':>5} {'socT_a%':>7} {'socT_o%':>7} {'avg_buy':>7}" ) print(header) print("-" * len(header)) tot_actual = tot_oracle = tot_gap = 0.0 days_ok = 0 day = d_from while day <= d_to: slots = await _load_day(conn, site_id, day) if len(slots) < 90 or all(s.soc_pct is None for s in slots): print(f"{day.date()!s:<11} — přeskočeno (slotů: {len(slots)})") day += timedelta(days=1) continue soc0_pct = next(s.soc_pct for s in slots if s.soc_pct is not None) socT_pct = next(s.soc_pct for s in reversed(slots) if s.soc_pct is not None) soc_start_wh = soc0_pct / 100.0 * bat.usable_wh soc_end_actual_wh = socT_pct / 100.0 * bat.usable_wh avg_buy = sum(s.buy for s in slots) / len(slots) actual_cash = _actual_cashflow_czk(slots) oracle_cash, soc_end_oracle_wh = solve_oracle(slots, bat, grid, soc_start_wh, avg_buy) # férové očištění: koncový SoC obou stran oceněn avg_buy dne actual_adj = actual_cash - soc_end_actual_wh / 1000.0 * avg_buy oracle_adj = oracle_cash - soc_end_oracle_wh / 1000.0 * avg_buy gap = actual_adj - oracle_adj gap_pct = (gap / abs(oracle_adj) * 100.0) if abs(oracle_adj) > 1e-6 else float("nan") print( f"{day.date()!s:<11} {actual_adj:>9.1f} {oracle_adj:>9.1f} {gap:>8.1f} {gap_pct:>5.0f}% " f"{soc0_pct:>5.0f} {socT_pct:>7.0f} {soc_end_oracle_wh / bat.usable_wh * 100:>7.0f} {avg_buy:>7.2f}" ) tot_actual += actual_adj tot_oracle += oracle_adj tot_gap += gap days_ok += 1 day += timedelta(days=1) print("-" * len(header)) if days_ok: print( f"{'CELKEM':<11} {tot_actual:>9.1f} {tot_oracle:>9.1f} {tot_gap:>8.1f}" f" ({days_ok} dní; Kč, SoC-adjusted; gap = forecast error + neefektivita dispatche)" ) else: print("Žádný den s kompletním auditem.") finally: await conn.close() def main() -> None: p = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) p.add_argument("--site-code", required=True) p.add_argument("--from", required=True, help="YYYY-MM-DD (Prague)") p.add_argument("--to", required=True, help="YYYY-MM-DD (Prague), včetně") p.add_argument("--dsn", default=None) args = p.parse_args() asyncio.run(run_report(args)) if __name__ == "__main__": sys.exit(main())