"""Validate TorsoCAE DGFV/VOF against Lobovsky et al. (2014)."""

from __future__ import annotations

import argparse
import json
import runpy
import sys
import tempfile
import zipfile
from pathlib import Path

import numpy as np


CASE_DIR = Path(__file__).resolve().parent
REPO = CASE_DIR.parents[3]
JOURNAL = REPO / "journals" / "cfd_vof_dam_break1.py"
sys.path.insert(0, str(REPO / "server"))

TANK_LENGTH = 1.610
RESERVOIR_LENGTH = 0.600
INITIAL_WATER_HEIGHT = 0.300
GRAVITY = 9.81
REFERENCE_FRONT_SPEED = 1.56
COARSE_GRID_FRONT_SPEED_ERROR_MAX = 0.15
REFERENCE_CITATION = (
    "Lobovsky, L. et al. (2014), Journal of Fluids and Structures 48, "
    "407-434, doi:10.1016/j.jfluidstructs.2014.03.009"
)


def _project_npz(project: Path, destination: Path) -> Path:
    with zipfile.ZipFile(project) as archive:
        names = [name for name in archive.namelist() if name.startswith("results/fem/") and name.endswith(".npz")]
        if len(names) != 1:
            raise RuntimeError(f"Expected one FEM result in {project}, found {len(names)}")
        destination.write_bytes(archive.read(names[0]))
    return destination


def _history(npz_path: Path) -> dict[str, np.ndarray]:
    with np.load(npz_path, allow_pickle=False) as data:
        required = (
            "vof_history_time",
            "vof_history_phase_volume",
            "vof_history_phase_volume_error",
            "vof_history_phase_upper",
        )
        missing = [key for key in required if key not in data]
        if missing:
            raise RuntimeError(f"VOF result lacks validation history: {missing}")
        return {
            "time": np.asarray(data["vof_history_time"], dtype=np.float64),
            "phase_volume": np.asarray(data["vof_history_phase_volume"], dtype=np.float64),
            "phase_volume_error": np.asarray(
                data["vof_history_phase_volume_error"], dtype=np.float64
            ),
            "phase_upper": np.asarray(data["vof_history_phase_upper"], dtype=np.float64),
        }


def evaluate(npz_path: Path) -> dict:
    history = _history(npz_path)
    time = history["time"]
    front_x = history["phase_upper"][:, 0]
    time_star = time * np.sqrt(GRAVITY / INITIAL_WATER_HEIGHT)
    front_star = (front_x - RESERVOIR_LENGTH) / INITIAL_WATER_HEIGHT

    impact = np.flatnonzero(front_x >= TANK_LENGTH - 1.0e-12)
    if not len(impact):
        raise RuntimeError("The computed phase front did not reach the downstream wall")
    impact_index = int(impact[0])
    fit = (time_star > 1.0) & (np.arange(len(time_star)) < impact_index)
    if np.count_nonzero(fit) < 4:
        raise RuntimeError("Insufficient pre-impact samples after nondimensional time t*=1")
    slope, intercept = np.polyfit(time_star[fit], front_star[fit], 1)
    speed_error = abs(float(slope) - REFERENCE_FRONT_SPEED) / REFERENCE_FRONT_SPEED
    max_volume_error = float(np.max(np.abs(history["phase_volume_error"])))
    passed = bool(
        speed_error <= COARSE_GRID_FRONT_SPEED_ERROR_MAX
        and max_volume_error <= 1.0e-8
    )
    return {
        "case": "Lobovsky dry-bed dam break",
        "citation": REFERENCE_CITATION,
        "geometry": {
            "tank_length_m": TANK_LENGTH,
            "tank_height_m": 0.600,
            "reservoir_length_m": RESERVOIR_LENGTH,
            "initial_water_height_m": INITIAL_WATER_HEIGHT,
            "numerical_depth_m": 0.030,
        },
        "method": {
            "formulation": "Torso staggered DGFV/VOF, p=0",
            "front_definition": "maximum cell-vertex x with cell volume fraction >= 0.5",
            "fit_interval": "t* > 1 until first downstream-wall contact",
        },
        "reference": {"front_speed_over_sqrt_gH": REFERENCE_FRONT_SPEED},
        "computed": {
            "front_speed_over_sqrt_gH": float(slope),
            "front_fit_intercept": float(intercept),
            "front_speed_relative_error_pct": 100.0 * speed_error,
            "impact_time_s": float(time[impact_index]),
            "impact_time_star": float(time_star[impact_index]),
            "max_phase_volume_relative_error": max_volume_error,
            "history_samples": int(len(time)),
        },
        "history": {
            "time_s": time.tolist(),
            "time_star": time_star.tolist(),
            "front_x_m": front_x.tolist(),
            "front_x_over_H": front_star.tolist(),
            "phase_volume_m3": history["phase_volume"].tolist(),
            "phase_volume_relative_error": history["phase_volume_error"].tolist(),
        },
        "acceptance": {
            "front_speed_relative_error_max_pct": 100.0 * COARSE_GRID_FRONT_SPEED_ERROR_MAX,
            "phase_volume_relative_error_max": 1.0e-8,
            "downstream_wall_contact_required": True,
        },
        "pass": passed,
    }


def main() -> None:
    parser = argparse.ArgumentParser()
    parser.add_argument(
        "--project",
        type=Path,
        help="Evaluate an existing .tcae project instead of rerunning the journal",
    )
    args = parser.parse_args()
    with tempfile.TemporaryDirectory(prefix="torsocae-vof-validation-") as directory:
        if args.project:
            npz_path = _project_npz(args.project.resolve(), Path(directory) / "result.npz")
        else:
            namespace = runpy.run_path(str(JOURNAL), run_name="__main__")
            npz_path = Path(namespace["result"]["npz_path"])
        summary = evaluate(npz_path)
        (CASE_DIR / "results.json").write_text(json.dumps(summary, indent=2) + "\n")
        print(json.dumps(summary["computed"], indent=2))
        if not summary["pass"]:
            raise SystemExit("Dam-break validation did not meet its published comparison gates")


if __name__ == "__main__":
    main()