"""Reference 24 V recorder hold-up model. No component values are measured."""
from dataclasses import dataclass, asdict
import math
import numpy as np
from scipy.integrate import solve_ivp

@dataclass(frozen=True)
class Case:
    name: str = 'Nominal'
    capacitance_f: float = 10.0
    esr_ohm: float = .12
    initial_v: float = 24.0
    efficiency: float = .90
    idle_w: float = .03
    cutoff_v: float = 10.5
    current_limit_a: float = 2.0
    peak_w: float = 12.0
    steady_w: float = 6.0
    peak_s: float = 2.0

def operating_point(v, power, resistance):
    """Stable low-current root of P = (V - IR) I."""
    if v <= 0 or power < 0 or resistance < 0 or v*v < 4*resistance*power:
        raise ValueError('No stable constant-power operating point')
    current=2*power/(v+math.sqrt(v*v-4*resistance*power))
    return current, v-resistance*current

def validate(c):
    vals=asdict(c); vals.pop('name')
    if not all(math.isfinite(x) for x in vals.values()):raise ValueError('Finite inputs required')
    if min(c.capacitance_f,c.initial_v,c.cutoff_v,c.current_limit_a,c.peak_w,c.steady_w,c.peak_s)<=0:
        raise ValueError('Positive energy, voltage, load and timing inputs required')
    if not 0<c.efficiency<=1 or c.esr_ohm<0 or c.idle_w<0:raise ValueError('Invalid loss model')
    if c.initial_v<=c.cutoff_v:raise ValueError('Initial voltage must exceed cutoff')
    i,v=operating_point(c.initial_v,c.peak_w/c.efficiency+c.idle_w,c.esr_ohm)
    if v<=c.cutoff_v or i>=c.current_limit_a:raise ValueError('Reference load cannot start')

def simulate(c, max_step=.5):
    """Integrate cap voltage and independently accumulated load/converter/ESR energy.
    The converter supplies a regulated 24 V until input UVLO or current limit.
    The peak/steady discontinuity is an explicit integration boundary.
    """
    validate(c)
    if max_step<=0:raise ValueError('Positive integration step required')
    state=[c.initial_v,0.,0.,0.];rows=[];reason=None
    for start,end,pout in [(0.,c.peak_s,c.peak_w),(c.peak_s,5000.,c.steady_w)]:
        pin=pout/c.efficiency+c.idle_w
        def rhs(t,y):
            current,vin=operating_point(y[0],pin,c.esr_ohm)
            return [-current/c.capacitance_f,pout,pin-pout,current*current*c.esr_ohm]
        def uvlo(t,y):return operating_point(y[0],pin,c.esr_ohm)[1]-c.cutoff_v
        def limit(t,y):return c.current_limit_a-operating_point(y[0],pin,c.esr_ohm)[0]
        uvlo.terminal=limit.terminal=True;uvlo.direction=limit.direction=-1
        sol=solve_ivp(rhs,(start,end),state,rtol=1e-9,atol=1e-10,max_step=max_step,events=[uvlo,limit])
        if not sol.success:raise RuntimeError(sol.message)
        for t,y in zip(sol.t,sol.y.T):
            i,v=operating_point(y[0],pin,c.esr_ohm)
            rows.append(dict(time_s=float(t),capacitor_v=float(y[0]),converter_input_v=v,input_current_a=i,
                             load_power_w=pout,load_energy_j=float(y[1]),converter_loss_j=float(y[2]),esr_loss_j=float(y[3])))
        state=sol.y[:,-1]
        if any(len(x) for x in sol.t_events):
            reason='input UVLO' if len(sol.t_events[0]) else 'input current limit';break
    if reason is None:raise RuntimeError('No cutoff within the 5,000 s model horizon')
    energy=.5*c.capacitance_f*(c.initial_v*c.initial_v-state[0]*state[0])
    residual=energy-sum(state[1:])
    return dict(case=asdict(c),duration_s=rows[-1]['time_s'],cutoff_reason=reason,
                energy_removed_j=energy,load_energy_j=float(state[1]),converter_loss_j=float(state[2]),
                esr_loss_j=float(state[3]),energy_residual_j=float(residual),
                maximum_input_current_a=max(r['input_current_a'] for r in rows),trace=rows)

def ideal_duration(c):
    """Independent zero-ESR energy solution; upper bound if ESR > 0."""
    energy=.5*c.capacitance_f*(c.initial_v**2-c.cutoff_v**2)
    peak=c.peak_w/c.efficiency+c.idle_w
    if energy<=peak*c.peak_s:return energy/peak
    return c.peak_s+(energy-peak*c.peak_s)/(c.steady_w/c.efficiency+c.idle_w)

def transfer(capacitance_f=.0047*.8,resistance=.05,gap_s=.002,power_w=12.,initial_v=24.):
    """Output capacitor supplies constant recorder load during an assumed transfer gap."""
    if capacitance_f<=0 or gap_s<0:raise ValueError('Invalid transfer capacitor/gap')
    def rhs(t,y):return [-operating_point(y[0],power_w,resistance)[0]/capacitance_f]
    if gap_s==0:return operating_point(initial_v,power_w,resistance)[1]
    sol=solve_ivp(rhs,(0,gap_s),[initial_v],rtol=1e-11,atol=1e-11,max_step=max(gap_s/50,1e-8))
    return operating_point(float(sol.y[0,-1]),power_w,resistance)[1]