"""Exact rectangular-inflow response, all external times in minutes."""
import math,csv
from pathlib import Path

def sample(k,t,amplitude=100,duration=30):
    if not all(math.isfinite(x) for x in (k,t,amplitude,duration)) or k<=0 or t<0 or amplitude<0 or duration<=0:
        raise ValueError('Invalid physical inputs')
    rise=amplitude*(-math.expm1(-min(t,duration)/k))
    q=rise*math.exp(-max(0,t-duration)/k)
    storage=k*60*q
    vin=amplitude*min(t,duration)*60
    return dict(t_minutes=t,k_minutes=k,q_m3s=q,storage_m3=storage,in_m3=vin,out_m3=vin-storage)

def main():
    rows=[]
    for k in (10,30):
        for t in (0,10,30,60,120,360):
            r=sample(k,t);rows.append(r)
            assert r['storage_m3']>=0 and r['out_m3']>=-1e-8
            assert abs(r['in_m3']-r['out_m3']-r['storage_m3'])<1e-8
            print(f"K={k:2} min, t={t:3} min, Q={r['q_m3s']:.6f} m3/s, S={r['storage_m3']:.3f} m3")
        assert math.isclose(sample(k,30+k)['q_m3s']/sample(k,30)['q_m3s'],math.exp(-1),rel_tol=1e-12)
    assert math.isclose(sample(10,30)['q_m3s'],95.0212931632136,rel_tol=1e-12)
    assert math.isclose(sample(30,30)['q_m3s'],63.2120558828558,rel_tol=1e-12)
    with (Path(__file__).parent/'reference.csv').open('w',newline='',encoding='utf-8-sig') as f:
        writer=csv.DictWriter(f,fieldnames=rows[0].keys());writer.writeheader();writer.writerows(rows)
    print('Exact reference values, conservation and exponential recession verified.')

if __name__=='__main__':main()