Review index · Download raw source
reference/cloud/run_grid.py
1"""Run pinned independent finite-cylinder grid; resumable per-case outputs."""2from pathlib import Path3import sys,json,hashlib,time,datetime,argparse4import numpy as np5from cylinder_mc import Geometry,load_phases,simulate,phase_qa,intervals,segment_inside,sha6def main():7 ap=argparse.ArgumentParser()8 ap.add_argument("--lut",required=True);ap.add_argument("--out",required=True)9 ap.add_argument("--packets",type=int,default=50000);ap.add_argument("--batches",type=int,default=6)10 args=ap.parse_args();out=Path(args.out);out.mkdir(exist_ok=True)11 for name in ["batch_info","batch_kernels"]: (out/name).mkdir(exist_ok=True)12 phases,meta=load_phases(args.lut);geo=Geometry()13 q={"phase_input":meta,"phases":{k:phase_qa(v)for k,v in phases.items()}}14 q["geometry_tests"]=[]15 tests=[([0,0,15000],[0,0,1],0,2000,True),16 ([0,0,15000],[0,0,-1],0,1000,True),17 ([0,0,15000],[1,0,0],0,30000,True),18 ([0,0,18000],[0,0,-1],1000,4000,True),19 ([0,0,18000],[0,0,1],None,None,False),20 ([0,0,17000],[0,0,1],None,None,False),21 ([0,0,17000],[0,0,-1],0,3000,True),22 ([40000,0,15000],[-1,0,0],10000,70000,True)]23 for pos,vec,exlo,exhi,hit in tests:24 lo,hi,h=intervals(np.array([pos],float),np.array([vec],float),geo)25 passed=bool(h[0])==hit and(not hit or(abs(lo[0]-exlo)<1e-8 and abs(hi[0]-exhi)<1e-8))26 q["geometry_tests"].append({"position_m":pos,"direction":vec,"passed":bool(passed)})27 q["all_geometry_tests_pass"]=all(x["passed"]for x in q["geometry_tests"])28 if not q["all_geometry_tests_pass"] or not all(x["pass"]for x in q["phases"].values()):raise RuntimeError("Preflight QA failed")29 (out/"phase_and_geometry_QA.json").write_text(json.dumps(q,indent=2),encoding="utf-8")30 manifest={"status":"RUNNING","started_utc":datetime.datetime.now(datetime.timezone.utc).isoformat(),31 "packets_per_batch":args.packets,"independent_batches":args.batches,32 "total_tau":[0,.3,1,3,10,30,45],"source_altitude_km":[15,16,17,18],33 "geometry":geo.__dict__,"dt_s":2e-6,"response_window_s":.012,"max_scatter":20000,34 "phase_input":meta,"engine_sha256":sha(Path(__file__).parent/"cylinder_mc.py"),35 "runner_sha256":sha(__file__),"case_ids":[]}36 manifest_path=out/"run_manifest.json"37 manifest_path.write_text(json.dumps(manifest,indent=2),encoding="utf-8")38 start=time.perf_counter()39 for phase_i,(phase_name,phase)in enumerate(phases.items()):40 for alt_i,source in enumerate([15000.,16000.,17000.,18000.]):41 for ti,tau in enumerate([0.,.3,1.,3.,10.,30.,45.]):42 case=f"{phase_name}_z{int(source/1000)}_tau{str(tau).replace('.','p')}"43 info_path=out/"batch_info"/(case+".json");kernel_path=out/"batch_kernels"/(case+".npz")44 if info_path.exists()and kernel_path.exists():45 print("REUSE",case,flush=True)46 else:47 infos=[];hists=[]48 for bi in range(args.batches if tau else 1):49 seed=2026102200+phase_i*100000+alt_i*10000+ti*100+bi50 n=args.packets if tau else 151 info,h=simulate(phase,source,tau,n,seed,geo,importance=.5)52 info["phase_id"]=phase_name53 infos.append(info);hists.append(h)54 info_path.write_text(json.dumps({"case_id":case,"geometry":geo.__dict__,"batches":infos},indent=2),encoding="utf-8")55 np.savez_compressed(kernel_path,batch_mass_relative=np.array(hists),dt_s=np.array(2e-6),window_s=np.array(.012))56 a=np.array([x["total_relative_to_vacuum"]for x in infos])57 se=float(a.std(ddof=1)/np.sqrt(len(a)))if len(a)>1 else 0.58 print(json.dumps({"case":case,"mean_gain":float(a.mean()),"batch_SE":se,"seconds":sum(x["elapsed_s"]for x in infos)}),flush=True)59 manifest["case_ids"].append(case)60 manifest["elapsed_s"]=time.perf_counter()-start61 manifest_path.write_text(json.dumps(manifest,indent=2),encoding="utf-8")62 manifest["status"]="COMPLETE";manifest["ended_utc"]=datetime.datetime.now(datetime.timezone.utc).isoformat()63 manifest["elapsed_s"]=time.perf_counter()-start64 manifest["nonvacuum_total_packets"]=2*4*6*args.batches*args.packets65 manifest_path.write_text(json.dumps(manifest,indent=2),encoding="utf-8")66 print("COMPLETE",manifest["elapsed_s"],flush=True)67if __name__=="__main__":main()68