"""A conventional cold-compression / warm-expansion closed nitrogen loop.

Cold nitrogen is an explicitly consumed resource. This is a steady component
screen with assumed machine maps; it is not a new source of energy or a
transient prediction of a self-starting pressure chamber.
"""
from pathlib import Path
import os,json,math
os.environ['OPENBLAS_NUM_THREADS']='1'
import numpy as np
from CoolProp.CoolProp import PropsSI
R=Path(__file__).resolve().parents[1]
P0=101325.;TA=293.15;DUTY=20.;M=2.;MAIR=30.;CAIR=1005.;RAIR=287.05
PL=55000.;PH=110000.;DP=2000.;EG=.9;EM=.9;FANETA=.65;CONTROLS=100.
hl=PropsSI('H','P',P0,'Q',0,'Nitrogen');sl=PropsSI('S','P',P0,'Q',0,'Nitrogen')
hv=PropsSI('H','P',P0,'Q',1,'Nitrogen');sv=PropsSI('S','P',P0,'Q',1,'Nitrogen')
tb=PropsSI('T','P',P0,'Q',0,'Nitrogen');ha=PropsSI('H','P',P0,'T',TA,'Nitrogen');sa=PropsSI('S','P',P0,'T',TA,'Nitrogen')
def availability(h,s):return h-ha-TA*(s-sa)

def cycle(cold=90.,hot=290.,source=300.,et=.8,ec=.8,air_dp=200.,coupling='electric'):
    # Compressor inlet, compressor outlet, warm turbine inlet, turbine outlet.
    h1=PropsSI('H','P',PL,'T',cold,'Nitrogen');s1=PropsSI('S','P',PL,'T',cold,'Nitrogen')
    h2s=PropsSI('H','P',PH+DP,'S',s1,'Nitrogen');h2=h1+(h2s-h1)/ec
    s2=PropsSI('S','P',PH+DP,'H',h2,'Nitrogen')
    h3=PropsSI('H','P',PH,'T',hot,'Nitrogen');s3=PropsSI('S','P',PH,'T',hot,'Nitrogen')
    h4s=PropsSI('H','P',PL+DP,'S',s3,'Nitrogen');h4=h3-et*(h3-h4s)
    s4=PropsSI('S','P',PL+DP,'H',h4,'Nitrogen')
    assert s2>=s1 and s4>=s3
    wc=M*(h2-h1);wt=M*(h3-h4);qhot=M*(h3-h2);qcold=M*(h4-h1)
    assert qhot>0 and qcold>0
    fan_shaft=MAIR*CAIR*source*((1+air_dp/P0)**(RAIR/CAIR)-1)/FANETA
    fan_dc=fan_shaft/EM;hotin=source+fan_shaft/(MAIR*CAIR);hotout=hotin-qhot/(MAIR*CAIR)
    if coupling=='shaft':
        drive=wc/.96;surplus=wt-drive
        generated=max(0.,surplus)*EG;pump_dc=max(0.,-surplus)/EM
        drive_loss=drive-wc+max(0.,surplus)*(1-EG)+max(0.,-surplus)*(1/EM-1)
    else:
        pump_dc=wc/EM;generated=wt*EG
        drive_loss=wt*(1-EG)+pump_dc-wc
    control=CONTROLS*DUTY;net=generated-pump_dc-fan_dc-control
    boil=qcold/(hv-hl);remaining=6-boil
    reject=drive_loss+(fan_dc-fan_shaft)+control
    first=M*(qhot/M-qcold/M)-(wt-wc)
    complete=qcold+MAIR*CAIR*(hotout-source)-(pump_dc+fan_dc+control-generated-reject)
    air_ex=lambda T:MAIR*CAIR*((T-TA)-TA*math.log(T/TA))
    cold_used=boil*(availability(hl,sl)-availability(hv,sv))
    resources=cold_used+air_ex(source)-air_ex(hotout)
    destroyed=resources-net
    entropy=boil*(sv-sl)+MAIR*CAIR*math.log(hotout/source)+reject/TA
    identity=destroyed-TA*entropy
    assert abs(first)<1e-6 and abs(complete)<1e-6 and abs(identity)<1e-6 and destroyed>=0 and entropy>=0
    hs=np.linspace(h2,h3,121);tc=np.array([PropsSI('T','P',PH,'H',h,'Nitrogen') for h in hs]);th=hotout+(hs-h2)/(h3-h2)*(hotin-hotout);delta=th-tc
    hc=np.linspace(h1,h4,121);tcr=np.array([PropsSI('T','P',PL,'H',h,'Nitrogen') for h in hc]);dc=tcr-tb
    def ua(h,d):
        a,b=d[:-1],d[1:];lm=np.divide(b-a,np.log(b/a),out=(a+b)/2,where=abs(b-a)>1e-9)
        return float(sum(M*np.diff(h)/DUTY/lm))
    feasible=min(delta)>0 and min(dc)>0 and remaining>.35
    return {'coupling':coupling,'shaft_transmission_efficiency_target':.96 if coupling=='shaft' else None,'cold_compressor_inlet_K':cold,'warm_turbine_inlet_K':hot,'heat_source_initial_K':source,'heat_source_final_K':hotout,
        'turbine_isentropic_efficiency_target':et,'compressor_isentropic_efficiency_target':ec,'hot_air_pressure_loss_Pa_target':air_dp,
        'working_nitrogen_throughput_kg':M,'duty_s':DUTY,'working_loop_pressure_low_Pa':PL,'working_loop_pressure_high_Pa':PH,'each_core_N2_pressure_loss_Pa_target':DP,
        'compressor_shaft_J':wc,'turbine_shaft_J':wt,'compressor_DC_J':pump_dc,'generator_DC_J':generated,'warm_air_blower_DC_J':fan_dc,'controls_DC_J':control,'net_steady_DC_J':net,
        'heat_supplied_J':qhot,'cold_bath_heat_received_J':qcold,'coolant_nitrogen_boiloff_kg':boil,'coolant_remaining_before_unknown_startup_kg':remaining,
        'nitrogen_bath_temperature_K':tb,'cold_resource_availability_used_J':cold_used,'cold_vapour_retained_availability_J':boil*availability(hv,sv),
        'total_resource_availability_used_J':resources,'exergy_destroyed_J':destroyed,'total_entropy_generation_J_K':entropy,
        'first_law_loop_residual_J':first,'complete_first_law_residual_J':complete,'exergy_entropy_identity_residual_J':identity,
        'warm_core_minimum_pinch_K':float(min(delta)),'cold_core_minimum_pinch_K':float(min(dc)),
        'warm_core_required_UA_W_K':ua(hs,delta) if min(delta)>0 else None,'cold_core_required_UA_W_K':ua(hc,dc) if min(dc)>0 else None,
        'positive_net_DC':net>0,'thermal_inventory_gate':feasible,
        'nodes':[{'name':n,'pressure_Pa':p,'enthalpy_J_kg':h,'temperature_K':PropsSI('T','P',p,'H',h,'Nitrogen')} for n,p,h in [('Cold compressor inlet',PL,h1),('Compressor outlet',PH+DP,h2),('Warm turbine inlet',PH,h3),('Turbine outlet',PL+DP,h4)]]}

def main():
    rows=[cycle(cold=c,et=t,ec=k,coupling=drive) for c in [90.,100.,110.,130.] for t in [.65,.7,.75,.8,.85] for k in [.65,.7,.75,.8,.85] for drive in ['electric','shaft']]
    selected=cycle(coupling='shaft');electrical=cycle();out={'physical_measurements':0,'classification':'STEADY_COLD_PRESSURE_LOOP_OPPORTUNITY_NOT_HARDWARE_VALIDATION','selected':selected,'electrical_conversion_control':electrical,'same_thermal_cycle_net_gain_J':selected['net_steady_DC_J']-electrical['net_steady_DC_J'],'cases':rows,
       'inputs':{'prepared_coolant_N2_kg':6.,'finite_heat_source_air_kg':MAIR,'generator_efficiency_target':EG,'motor_efficiency_target':EM,'air_blower_efficiency_target':FANETA,'startup_inventory_and_wet_mass_closed':False,'full_6kg_liquid_preparation_reversible_floor_J':6*availability(hl,sl)},
       'limits':['A closed pure-nitrogen loop uses cold compression and warm expansion. This is conventional thermal work, not energy from an empty region.','Each heat exchanger has an explicit 2 kPa nitrogen pressure-loss target. Installed machine efficiencies, pressure loss and UA are unmeasured.','The same-shaft target uses 96% transmission efficiency; surplus shaft work alone enters the generator. Shaft-speed matching, mechanical losses and startup motor demand still need hardware evidence.','The bath is an already prepared liquid-nitrogen resource. Boiled cold vapour and unspent liquid remain in the availability ledger.','Startup, working-gas inventory, wall transients, vessel mass, hot/cold seals, compressor maps, controls and real pressure dynamics are not solved.','No electrical energy or vehicle-performance credit is assigned to the driving model.']}
    (R/'output/cold-pressure-loop.json').write_text(json.dumps(out,indent=2)+'\n')
    print(json.dumps({'selected':selected,'positive_cases':sum(z['positive_net_DC'] and z['thermal_inventory_gate'] for z in rows),'case_count':len(rows)},indent=2))
if __name__=='__main__':main()