"""Два переходных опыта; скорость фиксирована, экстраполяция карты запрещена."""
import sys,json,math
ROOT=globals().get('TVV_ROOT','')
if ROOT not in sys.path:sys.path.insert(0,ROOT)
from run_record import start,finish,sha,save_json
from numerics import FluxMap,newton_bounded
from energy_map import from_record
from flux_dynamics import implicit_flux_step,open_field_step
from map_processing import magnetic_energy
from regimes import electrical
source=globals().get('TVV_MAGNETIC_MAP','')
if not source:raise ValueError('Задайте TVV_MAGNETIC_MAP: проверенный magnetic-map.json')
with open(source) as f:data=json.loads(f.read())
if not data['validated']:raise ValueError('Карта не прошла независимую проверку')
folder,record=start(ROOT,globals().get('TVV_OUTPUT_ROOT',ROOT+'/runs'),'dynamics',['run_dynamics.py','flux_dynamics.py','numerics.py','map_processing.py','energy_map.py','regimes.py','run_record.py','machine.json','dynamic-plan.json'],'two fixed-speed nonlinear electrical transients')
record['map_source']={'path':source,'sha256':sha(source)};save_json(folder+'/magnetic-map.json',data)
with open(folder+'/snapshot/machine.json') as f:machine=json.loads(f.read())
with open(folder+'/snapshot/dynamic-plan.json') as f:plan=json.loads(f.read())
magnetic_map=from_record(data);rs=machine['resistance']['stator_phase_ohm'];rf=machine['resistance']['field_ohm'];omega=2*math.pi*50
anchor=data['anchor_current_A'];anchor_energy=data['anchor_energy_J'];voltage0=rf*1500.
def steady_load(resistance):
    def residual(current):
        flux,jac=magnetic_map.evaluate(current);v=electrical(current,flux,rs)
        values=[v['ud_V']+resistance*current[0],v['uq_V']+resistance*current[1],voltage0-rf*current[2]]
        matrix=[[(rs+resistance if j==0 else 0)-omega*jac[1][j] for j in range(3)],
                [(rs+resistance if j==1 else 0)+omega*jac[0][j] for j in range(3)], [0.,0.,-rf]]
        return values,matrix
    result=newton_bounded(residual,[-20.,-400.,1500.],[(a[0],a[-1]) for a in magnetic_map.axes],[20000.,20000.,200.],1e-10)
    if not result['converged']:raise ValueError('Начальная нагрузка не находится в карте')
    return result['x']
def simulate(kind,factor):
    current=[0.,0.,1500.] if kind=='field-step' else steady_load(50.)
    flux,_=magnetic_map.evaluate(current);time=0.;duration=plan['field_duration_s'] if kind=='field-step' else plan['load_duration_s'];rows=[];input_work=0.;event=plan['event_s']
    energy0=magnetic_energy(magnetic_map,current,anchor,anchor_energy)
    step=factor*(plan['field_step_s'] if kind=='field-step' else plan['load_initial_step_s'])
    while True:
        active=time>=event-1e-12;voltage=voltage0*(1.02 if active and kind=='field-step' else 1.)
        resistance=(40. if active else 50.) if kind=='load-step' else None
        torque=1.5*(flux[0]*current[1]-flux[1]*current[0]);cu=1.5*rs*(current[0]**2+current[1]**2)+rf*current[2]**2
        load=1.5*resistance*(current[0]**2+current[1]**2) if resistance is not None else 0.
        power=voltage*current[2]-cu-load-omega*torque
        energy=magnetic_energy(magnetic_map,current,anchor,anchor_energy)
        rows.append({'time_s':time,'currents_A':current[:],'flux_Wb':flux[:],'energy_J':energy,'input_work_J':input_work,
                     'energy_balance_residual_J':energy-energy0-input_work,'field_voltage_V':voltage,'load_ohm':resistance,'load_power_W':load,'torque_Nm':torque})
        if time>=duration-1e-12:break
        limit=factor*(plan['field_step_s'] if kind=='field-step' else plan['load_slow_step_s'] if time>=event+plan['load_fast_window_s'] else plan['load_fast_step_s'])
        dt=min(step,limit,duration-time)
        if time<event and time+dt>=event:dt=event-time
        target_time=time+dt
        # Use the old side at the exact event endpoint; the next step starts with
        # the new voltage/load, so no impulse of work is inserted by averaging.
        vnew=voltage0*(1.02 if kind=='field-step' and time>=event-1e-12 else 1.)
        rnew=40. if time>=event-1e-12 else 50.
        if kind=='field-step':new_flux,new_current=open_field_step(magnetic_map,flux[2],current[2],dt,rf,vnew)
        else:new_flux,new_current,_=implicit_flux_step(magnetic_map,flux,current,dt,rs,rf,omega,rnew,vnew)
        new_torque=1.5*(new_flux[0]*new_current[1]-new_flux[1]*new_current[0])
        new_power=vnew*new_current[2]-rf*new_current[2]**2-omega*new_torque
        if kind=='load-step':new_power-=1.5*(rs+rnew)*(new_current[0]**2+new_current[1]**2)
        input_work+=dt*(power+new_power)/2
        time=target_time;current=new_current;flux=new_flux
        step=factor*plan['load_initial_step_s'] if abs(time-event)<1e-12 and kind=='load-step' else min(limit,step*plan['step_growth'])
    return rows
try:
    results={}
    for kind in ('field-step','load-step'):
        trajectories={}
        for factor in (1.,.5):
            rows=simulate(kind,factor);label=kind+('-half-step' if factor==.5 else '')
            trajectories[factor]=rows
            save_json(folder+'/'+label+'.json',rows);results[label]={'steps':len(rows)-1,'final':rows[-1],
                'energy_balance_scaled_error':max(abs(r['energy_balance_residual_J']) for r in rows)/max(plan['energy_absolute_scale_J'],max(abs(r['energy_J']-rows[0]['energy_J']) for r in rows))}
        coarse=trajectories[1.];fine=trajectories[.5];index=0;changes=[0.,0.,0.]
        scales=[max(1.,max(abs(r['currents_A'][j]) for r in fine)) for j in range(3)]
        for row in coarse:
            t=row['time_s']
            while index<len(fine)-2 and fine[index+1]['time_s']<t:index+=1
            a,b=fine[index],fine[index+1];weight=(t-a['time_s'])/(b['time_s']-a['time_s'])
            for j in range(3):
                reference=a['currents_A'][j]+weight*(b['currents_A'][j]-a['currents_A'][j])
                changes[j]=max(changes[j],abs(row['currents_A'][j]-reference)/scales[j])
        results[kind]['step_halving_current_errors']=changes
        if max(changes)>plan['relative_tolerance'] or results[kind]['energy_balance_scaled_error']>plan['relative_tolerance']:raise ValueError('Динамика не прошла шаг/энергию')
    save_json(folder+'/results.json',results);finish(folder,record,'passed_time_step_and_energy_pending_independent_trajectory_FEMM')
except Exception as error:record['error']=str(error);save_json(folder+'/results.json',locals().get('results',{}));finish(folder,record,'failed');raise
print('TVV_DYNAMICS_RESULT',folder)
