#!/usr/bin/env python3 """Calibrate a read-only model export using file-isolated native EPANET runs.""" from __future__ import annotations import argparse import asyncio import os import sys import traceback from datetime import datetime, timedelta, timezone from pathlib import Path from uuid import uuid4 def main() -> int: parser = argparse.ArgumentParser(description=__doc__) parser.add_argument("--design", type=Path, default=Path("resources/benchmarks/tjwater2-v2-burst-simulation.json")) parser.add_argument("--output-root", type=Path, default=Path("data/benchmark-calibration")) parser.add_argument("--duration-seconds", type=int, help="Override calibration duration, e.g. 18000 for five hours") parser.add_argument("--timeout-seconds", type=int, default=240) parser.add_argument("--all-severities", action="store_true", help="Retain weak-signal sizes even if the severe level fails") args = parser.parse_args() sys.path.insert(0, str(Path(__file__).resolve().parents[1])) os.umask(0o077) report = None report_path = None try: from psycopg.conninfo import make_conninfo from app.infra.db.project_routing import ActiveProjectRouting, activate_project_routing from app.native.wndb.core.connection import close_all_project_pools, project_transaction from app.native.wndb.inp.exporter import dump_inp from scripts.check_burst_simulation_design import readonly, resolve_project, validate_database, validate_tool_inputs from scripts.benchmark.simulation_design import load_simulation_design from scripts.benchmark.calibration import assess, configure_inp, digest, pressure_preview, run_file, save_json, select_calibrated_area design = load_simulation_design(args.design) validate_tool_inputs(design) duration = args.duration_seconds if args.duration_seconds is not None else design.duration_seconds step = design.model_expectations.report_timestep_seconds if not 3600 <= duration <= 18000 or duration % step or not 1 <= args.timeout_seconds <= 1800: raise ValueError("duration must be whole report steps between one and five hours; timeout 1..1800") biz, _ = asyncio.run(resolve_project(design)) with readonly(biz) as connection: validate_database(connection, design, template=False) folder = (args.output_root / (datetime.now(timezone.utc).strftime("%Y%m%dT%H%M%SZ") + "-" + uuid4().hex[:8])).resolve() folder.mkdir(parents=True, mode=0o700) snapshot = folder / "source.inp" readonly_dsn = make_conninfo(biz, options="-c default_transaction_read_only=on -c statement_timeout=20000") routing = ActiveProjectRouting(project_code=design.project_code, business_dsn=readonly_dsn) try: with activate_project_routing(routing), project_transaction(design.project_code) as connection: connection.execute("SET TRANSACTION ISOLATION LEVEL REPEATABLE READ, READ ONLY") # Stateless models permit event-aligned runs and independent warmup history. counts = connection.execute("SELECT (SELECT count(*) FROM network.tanks) tanks, " "(SELECT count(*) FROM network.pumps) pumps, " "(SELECT count(*) FROM network.controls) controls, " "(SELECT count(*) FROM network.rules) rules").fetchone() if any(counts.values()): raise ValueError("stateful model requires an explicit warm-start implementation") pda_options = {row["key"]: row["value"] for row in connection.execute( "SELECT key,value FROM network.simulation_settings WHERE engine_version='legacy' " "AND key=ANY(%s)", (["HEADERROR", "FLOWCHANGE", "MINIMUM PRESSURE", "PRESSURE EXPONENT"],), )} dump_inp(design.project_code, str(snapshot), "2") finally: close_all_project_pools() source = snapshot.read_text(encoding="utf-8") node_ids = sorted({node for pipe in design.pipes for node in pipe.output_node_ids}) link_ids = sorted({link for pipe in design.pipes for link in pipe.output_link_ids}) report = {"version": "burst-calibration@1", "status": "running", "created_at": datetime.now(timezone.utc), "source_effect": "read-only model export; local file simulations only", "model_sha256": digest(snapshot), "design_sha256": digest(args.design), "design": design.model_dump(mode="json"), "duration_seconds": duration, "pda_options": pda_options, "assessment_scope": "strong-signal calibration, not Agent accuracy or absence of a physical burst", "all_severities": args.all_severities, "engine_sha256": digest(Path("app/infra/epanet/linux/runepanet")), "cases": [], "selection": {}, "controls": {}, "limitations": ["No analysis scenario or SCADA rows have been persisted", "No worker, MCP, Agent, or 15-day detection has been executed", "The four-sample pressure replay rule is a calibration preview"]} report_path = folder / "calibration.json" save_json(report_path, report) def simulate(label, start, seconds, emitter=None): path = folder / f"{label}.inp" path.write_text(configure_inp(source, start, seconds, emitter=emitter, pda_options=pda_options), encoding="utf-8") result = run_file(path, node_ids, link_ids, timeout=args.timeout_seconds) save_json(folder / f"{label}.json", result) print(f"{label}: {result['status']}, solver warnings={len(result['warnings'])}", flush=True) return result controls = {} histories = {} for index, at in enumerate(design.event_times): controls[at] = simulate(f"normal-{index}", at, duration) warmup = simulate(f"warmup-{index}", at - timedelta(seconds=4 * step), 3 * step) if warmup.get("status") != "completed" or warmup.get("warnings") or warmup.get("times") != { "report_start": 0, "report_step": step, "sim_duration": 3 * step, "num_periods": 4 }: raise ValueError("normal warmup did not produce four converged samples") histories[at] = warmup["nodes"] report["controls"][at.isoformat()] = {"normal_file": f"normal-{index}.json", "warmup_file": f"warmup-{index}.json"} save_json(report_path, report) for pipe in design.pipes: cases = [case for case in design.cases if case.pipe_id == pipe.pipe_id] passed_levels = [] for severity in design.calibration_order: level_passed = True for case in cases: area = pipe.candidate_orifice_areas_cm2[severity] burst = simulate(f"{case.case_id}-{severity}", case.occurred_at, duration, (pipe.emitter_node_id, area)) normal = controls[case.occurred_at] result = assess(pipe, normal, burst, design.selection_rule, duration=duration, step=step) if normal.get("status") == burst.get("status") == "completed" and not any( reason.endswith("incomplete_or_failed") for reason in result["blockers"] ): for key, current in (("normal_trigger", normal), ("event_trigger", burst)): result[key] = pressure_preview(pipe, case.occurred_at, duration, step, histories[case.occurred_at], current["nodes"], design.selection_rule) if result["normal_trigger"]["raised"]: result["blockers"].append("normal_control_false_positive") if not result["event_trigger"]["raised"]: result["blockers"].append("event_rule_did_not_trigger") coefficient = 0.65 * area / 10000 * (19.6 ** 0.5) * 1000 # PDA demand can change with pressure: demand delta alone is not emitter flow. result["emitter_flow_m3_s"] = [coefficient * max(row["pressure"], 0) ** 0.5 / 1000 for row in burst["nodes"][pipe.emitter_node_id]] result.update(case_id=case.case_id, severity=severity, orifice_area_cm2=area, passed=not result["blockers"]) report["cases"].append(result) level_passed &= result["passed"] save_json(report_path, report) print(f"{case.case_id}-{severity}: physical and rule checks {'PASS' if result['passed'] else 'FAIL'}", flush=True) if level_passed: passed_levels.append(severity) if severity == "severe" and not level_passed and not args.all_severities: break case_ids = {case.case_id for case in cases} report["selection"][pipe.pipe_id] = select_calibrated_area( pipe, passed_levels, [case for case in report["cases"] if case["case_id"] in case_ids], args.all_severities, ) save_json(report_path, report) report["status"] = "completed" save_json(report_path, report) print(f"calibration complete; report: {report_path}", flush=True) return 0 except Exception as exc: frames = traceback.extract_tb(exc.__traceback__) location = frames[-1] if frames else None suffix = f" at {location.filename}:{location.lineno}" if location else "" if report is not None and report_path is not None: report.update(status="failed", failure_type=type(exc).__name__) save_json(report_path, report) print(f"calibration failed ({type(exc).__name__}){suffix}", file=sys.stderr) return 1 if __name__ == "__main__": raise SystemExit(main())