"""Reproduce the hypothetical MHD bypass ledger; not an engine model."""
import csv
from pathlib import Path
ROOT=Path(__file__).resolve().parent
FLOW=10.0; AVAILABLE_KJ_KG=100.0; AUX_KW=20.0
with (ROOT/"mhd-bypass-ledger.csv").open("w",newline="",encoding="utf-8") as f:
 w=csv.writer(f); w.writerow(["mass_flow_kg_s","assumed_available_energy_kJ_kg","extraction_fraction","generator_efficiency","extracted_stream_power_kW","generated_electric_power_kW","assumed_other_auxiliary_power_kW","assumed_field_excitation_power_kW","unmet_auxiliary_and_field_power_kW","accelerator_efficiency","accelerator_flow_work_kW","net_stream_energy_change_kW"])
 for frac in (0,.05,.10,.20,.40,.60):
  for eta_g in (.40,.60,.80):
   for eta_a in (.40,.60,.80):
    for field in (0,15,30,60,100):
     extracted=FLOW*AVAILABLE_KJ_KG*frac; generated=extracted*eta_g; to_accel=max(0.0,generated-AUX_KW-field); shortfall=max(0.0,AUX_KW+field-generated); work=to_accel*eta_a
     w.writerow([f"{FLOW:.1f}",f"{AVAILABLE_KJ_KG:.1f}",f"{frac:.2f}",f"{eta_g:.2f}",f"{extracted:.3f}",f"{generated:.3f}",f"{AUX_KW:.3f}",f"{field:.3f}",f"{shortfall:.3f}",f"{eta_a:.2f}",f"{work:.3f}",f"{-extracted+work:.3f}"])
extracted=200.0; generated=120.0; to_accel=max(0,generated-20.0-15.0); work=to_accel*.6
assert (to_accel,work,-extracted+work)==(85.0,51.0,-149.0)
assert max(0,120.0-200.0-15.0)==0 and max(0,20.0+15.0-0)==35.0
with (ROOT/"mhd-bypass-ledger.csv").open(newline="",encoding="utf-8") as f:
 rows=list(csv.DictReader(f))
assert len(rows)==270
for row in rows:
 flow=float(row["mass_flow_kg_s"]); q=float(row["assumed_available_energy_kJ_kg"])
 frac=float(row["extraction_fraction"]); eta_g=float(row["generator_efficiency"])
 aux=float(row["assumed_other_auxiliary_power_kW"]); field=float(row["assumed_field_excitation_power_kW"])
 eta_a=float(row["accelerator_efficiency"])
 removed=flow*q*frac; generated=removed*eta_g
 available=max(0.0,generated-aux-field); shortfall=max(0.0,aux+field-generated)
 returned=available*eta_a; delta=-removed+returned
 close=lambda actual,expected: abs(actual-expected)<1e-8
 assert close(float(row["extracted_stream_power_kW"]),removed)
 assert close(float(row["generated_electric_power_kW"]),generated)
 assert close(float(row["unmet_auxiliary_and_field_power_kW"]),shortfall)
 assert close(float(row["accelerator_flow_work_kW"]),returned)
 assert close(float(row["net_stream_energy_change_kW"]),delta)
 assert removed>=0 and generated>=0 and returned>=0 and returned<=generated+1e-8 and delta<=1e-8
boundary_cases=0
for flow in (1,10,50):
 for q in (20,100,300):
  for eg in (0,.4,.8,1):
   for ea in (0,.4,.8,1):
    for field in (0,15,100,1000):
     for frac in (0,.2,.6,1):
      removed=flow*q*frac; generated=removed*eg
      returned=max(0,generated-20-field)*ea
      boundary_cases+=1
      assert removed>=0 and 0<=returned<=generated+1e-8 and -removed+returned<=1e-8
     changes=[-flow*q*frac+max(0,flow*q*frac*eg-20-field)*ea for frac in (0,.2,.6,1)]
     assert all(b<=a+1e-8 for a,b in zip(changes,changes[1:]))
assert boundary_cases==2304
print(f"MHD ledger reproduced: {len(rows)} rows; {len(rows)} row-wise accounting checks and {boundary_cases} boundary/property cases passed.")