Files
Enelix-EMS/services/netplan-v4/netplan_v4/economic_replay.py
T

208 lines
14 KiB
Python

"""Causal daily economic family comparison using frozen plans and later physical data.
This is a labelled *simulation*, not savings on an invoice. A cohort is captured
near local midnight, uses only prices/forecasts then available, and all families
share initial storage, physical outcomes, tariff and constraints. No device calls.
Supported v1: one grid-charge-enabled battery; other topologies remain explicit.
"""
from dataclasses import asdict, replace
from datetime import datetime, timedelta
from hashlib import sha256
from math import sqrt
import json
from .domain import ZURICH, utc, month_key, quarter_start
from .selection import ReplayScore
from . import measurement_pipeline as m
def schema(con):
con.executescript('''
CREATE TABLE IF NOT EXISTS planner_economic_cohorts(
plant TEXT NOT NULL, local_day TEXT NOT NULL, issued_at INTEGER NOT NULL,
ends_at INTEGER NOT NULL, policy_id TEXT NOT NULL, value TEXT NOT NULL,
PRIMARY KEY(plant,local_day));
CREATE TABLE IF NOT EXISTS planner_economic_results(
plant TEXT NOT NULL, local_day TEXT NOT NULL, evaluated_at INTEGER NOT NULL,
policy_id TEXT NOT NULL, value TEXT NOT NULL, PRIMARY KEY(plant,local_day));
CREATE TABLE IF NOT EXISTS planner_economic_attempts(
plant TEXT NOT NULL, local_day TEXT NOT NULL, checked_at INTEGER NOT NULL,
reason TEXT NOT NULL, PRIMARY KEY(plant,local_day));
''')
def _policy(data, cfg):
batteries=[]
for b in data['batteries']:
v=asdict(b)
for k in ('soc_percent','measured_at','discharge_blocked'):v.pop(k)
batteries.append(v)
return sha256(m.canonical({'battery':batteries,'limits':asdict(data['limits']),
'mapping':cfg['mappingSha256'],'formula':cfg['formula'],'dataset':cfg['datasetId'],
'method':'frozen_day_ahead_grid_tracking_v1','external':'sdl_request_estimate','peakTariffs':data['peak_prices']}).encode()).hexdigest()
def capture(store, plant, now, assemble, optimizer):
settings=store.settings(plant)
if settings.get('forecastSource')!='corrected_profile': return
local=utc(now).astimezone(ZURICH); day=local.date().isoformat()
# Freeze only at the beginning of a local day. Never reconstruct an old forecast from hindsight.
if local.hour!=0 or local.minute>=15:return
if store.con.execute('SELECT 1 FROM planner_economic_cohorts WHERE plant=? AND local_day=?',(plant,day)).fetchone():return
attempted=store.con.execute('SELECT checked_at FROM planner_economic_attempts WHERE plant=? AND local_day=?',(plant,day)).fetchone()
if attempted and now.timestamp()-attempted[0]<300:return
reason='awaiting_comparable_snapshot'
try:
cfg=m.configuration(store.con,plant,settings['measurementDataset']);plans={}; common=None;policy=None
end=datetime.combine(local.date()+timedelta(days=1),datetime.min.time(),tzinfo=ZURICH).astimezone(utc(now).tzinfo)
for family in store.registry.entries():
data,_,quality=assemble(store,plant,family.key,now)
data['batteries']=[replace(b,roundtrip_efficiency=settings['roundtripEfficiency']) for b in data['batteries']]
if len(data['batteries'])!=1 or not data['batteries'][0].grid_charging:
raise ValueError('unsupported_replay_topology')
if data['steps'][-1].end<end:raise ValueError('published_prices_do_not_cover_day')
p=optimizer(**data,config_revision=settings['revision'],family=family.key,timeout_seconds=1.0)
if not p.get('executable'):raise ValueError('candidate_not_feasible')
points=[x for x in p['points'] if utc(x['time'])<end]
p['points']=points
if utc(points[-1]['validUntil'])!=end:raise ValueError('day_boundary_not_covered')
this=_policy(data,cfg)
if policy is not None and this!=policy:raise ValueError('inconsistent_candidate_context')
policy=this
initial=asdict(data['batteries'][0]);initial['measured_at']=utc(initial['measured_at']).isoformat()
context={'battery':initial,'limits':asdict(data['limits']),'peaks':data['observed_peaks'],
'peakTariffs':data['peak_prices'],'quarterPast':{utc(k).isoformat():asdict(v) for k,v in data['quarter_history'].items()}}
if common is not None and m.canonical(common)!=m.canonical(context):raise ValueError('different_initial_conditions')
common=context;plans[family.key]=p
billing=[[(p['time'],p['validUntil'],p['importPriceChfKwh'],p['exportPriceChfKwh']) for p in plan['points']] for plan in plans.values()]
if any(v!=billing[0] for v in billing):raise ValueError('Different priced intervals between families')
value={'datasetId':cfg['datasetId'],'start':next(iter(plans.values()))['points'][0]['time'],
'end':end.isoformat(),'context':common,'plans':plans,'policyId':policy,
'method':'frozen_day_ahead_grid_tracking_v1','actuation':False}
with store.con:store.con.execute('INSERT INTO planner_economic_cohorts VALUES(?,?,?,?,?,?)',
(plant,day,int(now.timestamp()),int(end.timestamp()),policy,m.canonical(value)))
reason='captured'
except (ValueError,KeyError,TypeError):
# No exception text copied from arbitrary data. Retry bounded to one attempt per five minutes.
reason='awaiting_comparable_snapshot'
with store.con:store.con.execute('INSERT INTO planner_economic_attempts VALUES(?,?,?,?) ON CONFLICT(plant,local_day) DO UPDATE SET checked_at=excluded.checked_at,reason=excluded.reason',
(plant,day,int(now.timestamp()),reason))
def actuals(records,cfg):
required=[s['key'] for s in cfg['sources'] if s['role'] in ('grid','physical_storage','sdl_request')]
if sum(s['role']=='sdl_request' for s in cfg['sources'])!=1:raise ValueError('Missing SDL outcome source')
def residual(values,c):
grid=storage=external=0.0
for s in c['sources']:
role=s['role']
if role not in ('grid','physical_storage','sdl_request'):continue
x=values[s['key']]*s['factorToW']
if role=='grid':grid+=x
elif role=='physical_storage':storage+=x
else:external+=x
# External SDL is an explicitly estimated historical contribution, never called a meter.
return grid-storage+external
return m.reconstruct(records,cfg,projection={'keys':required,'calculate':residual})
def simulate(cohort,family,outcomes):
plan=cohort['plans'][family];ctx=cohort['context'];b=ctx['battery'];limits=ctx['limits']
capacity=b['capacity_kwh'];initial=energy=capacity*b['soc_percent']/100
low=capacity*b['min_soc_percent']/100; high=capacity*b['max_soc_percent']/100
eta=sqrt(b['roundtrip_efficiency']);peaks=dict(ctx['peaks']);quarter_max={};quarters={};money=wear=0.; coverage=[];breaches=0
blocked=b['discharge_blocked'];rearm=capacity*(b['rearm_soc_percent'] or b['min_soc_percent'])/100
for stamp,q in ctx['quarterPast'].items():quarters[stamp]=[q['import_kwh'],q['measured_seconds']]
for p in plan['points']:
start,end=utc(p['time']),utc(p['validUntil']);t=int(start.timestamp());seconds=int((end-start).total_seconds());dt=seconds/3600
w=outcomes.get(t//300*300)
if not w or not w['profileUsable'] or not m.numeric(w['loadW']):raise ValueError('Incomplete actual outcome')
residual=w['loadW'];coverage.append(w['coverage'])
charge=min(b['max_charge_w'],max(0.,(high-energy)/eta/dt*1000))
if blocked and energy>=rearm-1e-9:blocked=False
discharge=0. if blocked else min(b['max_discharge_w'],max(0.,(energy-low)*eta/dt*1000))
target=p['gridTargetW'];cap=limits['import_w'];monthly=limits['manager_month_limits_w'].get(str(start.astimezone(ZURICH).month),limits['manager_month_limits_w'].get(start.astimezone(ZURICH).month))
if monthly is not None:cap=monthly if cap is None else min(cap,monthly)
if cap is not None:target=min(target,cap)
if limits['export_w'] is not None:target=max(target,-limits['export_w'])
battery=min(charge,max(-discharge,target-residual));grid=residual+battery
energy+=battery/1000*dt*(eta if battery>=0 else 1/eta)
if energy<low-1e-6 or energy>high+1e-6:raise ValueError('Replay SOC invariant')
if energy<=low+1e-9:blocked=True
if cap is not None and grid>cap+1.:breaches+=1
if limits['export_w'] is not None and grid < -limits['export_w']-1.:breaches+=1
money+=(max(grid,0)*p['importPriceChfKwh']-max(-grid,0)*p['exportPriceChfKwh'])/1000*dt
wear+=abs(battery)/1000*dt*b['throughput_chf_kwh']
q=quarter_start(start).isoformat();v=quarters.setdefault(q,[0.,0]);v[0]+=max(grid,0)/1000*dt;v[1]+=seconds
if end==quarter_start(start)+timedelta(minutes=15):
if v[1]!=900:raise ValueError('Incomplete simulated billing quarter')
month=month_key(start);peaks[month]=max(peaks[month],v[0]/.25);quarter_max[month]=max(quarter_max.get(month,0.),v[0]/.25)
additional=sum(max(0,v-ctx['peaks'][month])*ctx['peakTariffs'][month] for month,v in peaks.items())
# Identical terminal valuation for all families, known at cohort creation; separate from cash.
last=plan['points'][-1];buy=max(0.,last['importPriceChfKwh']);sell=max(0.,min(buy,last['exportPriceChfKwh']))
terminal=max(initial-energy,0)/eta*buy-max(energy-initial,0)*eta*sell
return {'cashCostChf':money+additional,'energyCostChf':money,'peakCostChf':additional,'throughputCostChf':wear,
'terminalAdjustmentChf':terminal,'costChf':money+additional+wear+terminal,
'initialEnergyKwh':initial,'finalEnergyKwh':energy,'coverage':min(coverage),
'quarterMaximaKw':quarter_max,'initialPeaksKw':ctx['peaks'],'peakTariffs':ctx['peakTariffs'],
'constraintBreaches':breaches,'terminalNormalized':True,
'terminalMethod':'common_known_end_price_inventory_valuation_not_physical_restoration'}
def advance(store,plant,now):
timestamp=int(now.timestamp())
rows=store.con.execute('SELECT c.* FROM planner_economic_cohorts c LEFT JOIN planner_economic_results r USING(plant,local_day) WHERE c.plant=? AND c.ends_at<=? AND c.ends_at>=? AND r.local_day IS NULL AND NOT EXISTS(SELECT 1 FROM planner_economic_attempts a WHERE a.plant=c.plant AND a.local_day=c.local_day AND a.checked_at>?) ORDER BY c.issued_at LIMIT 1',(plant,timestamp-120,timestamp-90*86400,timestamp-300)).fetchall()
for row in rows:
with store.con:store.con.execute('INSERT INTO planner_economic_attempts VALUES(?,?,?,?) ON CONFLICT(plant,local_day) DO UPDATE SET checked_at=excluded.checked_at,reason=excluded.reason',(plant,row['local_day'],timestamp,'evaluating_actuals'))
try:
cohort=json.loads(row['value']);cfg=m.configuration(store.con,plant,cohort['datasetId']);source=cfg.get('sourceDatasetId',cfg['datasetId']);start=m.epoch(cohort['start']);end=m.epoch(cohort['end'])
observed=store.con.execute('SELECT value FROM planner_observations WHERE plant=? AND dataset=? AND captured_at>=? AND captured_at<=? AND received_at<=? ORDER BY captured_at',(plant,source,start-600,end+600,timestamp))
windows=actuals([json.loads(r[0]) for r in observed],cfg);outcomes={w['start']:w for w in windows}
results={f:simulate(cohort,f,outcomes) for f in cohort['plans']}
value={'status':'evaluated','results':results,'start':cohort['start'],'end':cohort['end'],
'evaluationBasis':'physical_balance_with_sdl_request_estimate','method':cohort['method'],
'actualCashSavings':False,'controlEnabled':False}
with store.con:store.con.execute('INSERT INTO planner_economic_results VALUES(?,?,?,?,?)',(plant,row['local_day'],timestamp,row['policy_id'],m.canonical(value)))
except (ValueError,KeyError,TypeError):
# Retain unscored cohort. Lack of actual data must never become a zero cost.
pass
def scores(store,plant,now):
settings=store.settings(plant)
if settings.get('forecastSource')!='corrected_profile':return []
since=int((utc(now)-timedelta(days=settings['autoLookbackDays'])).timestamp())
recent=store.con.execute('SELECT policy_id,value FROM planner_economic_cohorts WHERE plant=? ORDER BY issued_at DESC LIMIT 1',(plant,)).fetchone()
if not recent or json.loads(recent[1])['datasetId']!=settings.get('measurementDataset'):return []
rows=store.con.execute('SELECT c.issued_at,c.ends_at,r.evaluated_at,r.value FROM planner_economic_results r JOIN planner_economic_cohorts c USING(plant,local_day) WHERE r.plant=? AND r.policy_id=? AND c.issued_at>=? AND r.evaluated_at<=? ORDER BY c.issued_at',(plant,recent[0],since,int(now.timestamp()))).fetchall()
if not rows:return []
keys=[f.key for f in store.registry.entries()];valid=[]
for r in rows:
v=json.loads(r['value']);out=v['results']
if set(out)!=set(keys) or any(out[k]['constraintBreaches'] for k in keys):continue
valid.append((r,v))
if not valid:return []
result=[]
for k in keys:
# Monthly demand cost is paid ONCE for the maximum, not once per replay day.
total=sum(v['results'][k]['energyCostChf']+v['results'][k]['throughputCostChf']+v['results'][k]['terminalAdjustmentChf'] for _,v in valid)
bases={}; maxima={}; rates={}
for _,v in valid:
d=v['results'][k]
for month,peak in d['quarterMaximaKw'].items():
bases.setdefault(month,d['initialPeaksKw'][month])
maxima[month]=max(maxima.get(month,0.),peak);rates[month]=d['peakTariffs'][month]
total+=sum(max(0.,maxima[month]-bases[month])*rates[month] for month in maxima)
cover=min(v['results'][k]['coverage'] for _,v in valid)
first=utc(valid[0][1]['start']);end=utc(valid[-1][1]['end']);available=datetime.fromtimestamp(max(r['evaluated_at'] for r,_ in valid),utc(now).tzinfo)
result.append(ReplayScore(k,recent[0],first,end,first,available,total,cover,len(valid)))
return result
def status(store,plant):
rows=store.con.execute('SELECT local_day,value FROM planner_economic_results WHERE plant=? ORDER BY local_day DESC LIMIT 14',(plant,))
evaluated=[{'day':r[0],**json.loads(r[1])} for r in rows]
count=store.con.execute('SELECT COUNT(*) FROM planner_economic_cohorts WHERE plant=?',(plant,)).fetchone()[0]
return {'connected':True,'method':'frozen_day_ahead_grid_tracking_v1','cohorts':count,'completedComparisons':evaluated,
'evaluationBasis':'physical_balance_with_sdl_request_estimate','isBillingEvidence':False,
'limitations':['one_grid_charging_battery','frozen_daily_plan_not_receding_horizon_field_replay']}