Upload folder using huggingface_hub

.gitattributes

@@ -1,35 +1,3 @@
-*.7z filter=lfs diff=lfs merge=lfs -text
-*.arrow filter=lfs diff=lfs merge=lfs -text
-*.bin filter=lfs diff=lfs merge=lfs -text
-*.bz2 filter=lfs diff=lfs merge=lfs -text
-*.ckpt filter=lfs diff=lfs merge=lfs -text
-*.ftz filter=lfs diff=lfs merge=lfs -text
-*.gz filter=lfs diff=lfs merge=lfs -text
-*.h5 filter=lfs diff=lfs merge=lfs -text
-*.joblib filter=lfs diff=lfs merge=lfs -text
-*.lfs.* filter=lfs diff=lfs merge=lfs -text
-*.mlmodel filter=lfs diff=lfs merge=lfs -text
-*.model filter=lfs diff=lfs merge=lfs -text
-*.msgpack filter=lfs diff=lfs merge=lfs -text
+*.pt filter=lfs diff=lfs merge=lfs -text
 *.npy filter=lfs diff=lfs merge=lfs -text
 *.npz filter=lfs diff=lfs merge=lfs -text
-*.onnx filter=lfs diff=lfs merge=lfs -text
-*.ot filter=lfs diff=lfs merge=lfs -text
-*.parquet filter=lfs diff=lfs merge=lfs -text
-*.pb filter=lfs diff=lfs merge=lfs -text
-*.pickle filter=lfs diff=lfs merge=lfs -text
-*.pkl filter=lfs diff=lfs merge=lfs -text
-*.pt filter=lfs diff=lfs merge=lfs -text
-*.pth filter=lfs diff=lfs merge=lfs -text
-*.rar filter=lfs diff=lfs merge=lfs -text
-*.safetensors filter=lfs diff=lfs merge=lfs -text
-saved_model/**/* filter=lfs diff=lfs merge=lfs -text
-*.tar.* filter=lfs diff=lfs merge=lfs -text
-*.tar filter=lfs diff=lfs merge=lfs -text
-*.tflite filter=lfs diff=lfs merge=lfs -text
-*.tgz filter=lfs diff=lfs merge=lfs -text
-*.wasm filter=lfs diff=lfs merge=lfs -text
-*.xz filter=lfs diff=lfs merge=lfs -text
-*.zip filter=lfs diff=lfs merge=lfs -text
-*.zst filter=lfs diff=lfs merge=lfs -text
-*tfevents* filter=lfs diff=lfs merge=lfs -text

README.md

@@ -0,0 +1,134 @@
+---
+license: cc-by-4.0
+language:
+- en
+tags:
+- OneScience
+- Earth Science
+- Physics-Informed Neural Network
+frameworks: PyTorch
+---
+
+<p align="center">
+  <strong><span style="font-size: 30px;">PINN-TC</span></strong>
+</p>
+
+# Model Introduction
+
+PINN-TC reconstructs three-dimensional tropical-cyclone wind and pressure-related fields from sparse flight, dropsonde, or simulated observations, using physical equations to constrain unobserved regions for vortex initialization and data-assimilation research.
+
+Paper: Realistic tropical cyclone wind and pressure fields can be reconstructed from sparse data using deep learning  
+https://doi.org/10.1038/s43247-023-01144-2
+
+# Model Description
+
+The method was proposed by research teams at Princeton University, the California Institute of Technology, Stanford University, and NOAA's Geophysical Fluid Dynamics Laboratory. The paper trains and validates the method with Hurricane Ida T-SHiELD simulations and NOAA/AOML flight-level, dropsonde, and Tail Doppler Radar observations. It is intended for sparse tropical-cyclone wind and geopotential reconstruction, physics-consistent constraints, and vortex-initialization research.
+
+# Use Cases
+
+| Use Case | Description |
+| :---: | :--- |
+| Tropical-cyclone reconstruction | Reconstruct wind, geopotential, and pressure-vertical-velocity fields from sparse coordinate observations. |
+| Physics-informed constraints | Constrain unobserved regions with horizontal momentum and pressure-continuity equations. |
+| ModelScope/OneCode execution | Validate sparse sampling, training, full-grid inference, scientific metrics, and visualization in ModelScope or OneCode environments. |
+| Multi-GPU training | Launch distributed data-parallel training through `torchrun`. |
+
+# Usage Instructions
+
+## 1.OneCode
+
+Experience intelligent, one-click AI4S programming through the OneCode online environment:
+
+[Try intelligent, one-click AI4S programming](https://web-2069360198568017922-iaaj.ksai.scnet.cn:58043/home)
+
+## 2. Download and Installation
+
+```bash
+hf download OneScience-Group/PINN-TC --local-dir ./PINN-TC
+cd PINN-TC
+```
+
+### Environment Dependencies
+
+**Hardware Requirements**
+
+- A GPU or DCU is recommended.
+- A CPU can be used to validate the workflow with the default small-sample configuration.
+- DCU users must install DTK in advance. DTK 25.04.2 or later, or the OneScience-recommended version matching the current cluster, is recommended.
+
+**DCU Environment**
+
+```bash
+# Activate DTK and Conda first
+conda create -n onescience311 python=3.11 -y
+conda activate onescience311
+pip install onescience[earth-dcu] -i http://mirrors.onescience.ai:3141/pypi/simple/ --trusted-host mirrors.onescience.ai
+pip install numpy pyyaml matplotlib
+```
+
+**GPU Environment**
+
+```bash
+# Activate Conda first
+conda create -n onescience311 python=3.11 -y libstdcxx-ng=12 libgcc-ng=12 gcc_linux-64=12 gxx_linux-64=12
+conda activate onescience311
+pip install onescience[earth-gpu] -i http://mirrors.onescience.ai:3141/pypi/simple/ --trusted-host mirrors.onescience.ai
+pip install numpy pyyaml matplotlib
+```
+
+### Training Data
+
+Training data consist of tropical-cyclone simulation fields and sparse observations containing storm-centered coordinates, horizontal winds, and geopotential. Model inputs are `[y,x,t,p]` coordinates, and targets are winds and geopotential at the same locations, while vertical velocity is constrained by the physical equations. The data cover a `±400 km` horizontal domain, `150-900 hPa`, and observation times at `-3/0/+3 h`. This repository uses a small number of synthetic samples to validate training, inference, and evaluation; they do not represent the official T-SHiELD or NOAA distribution, training scale, or formal paper performance.
+
+```bash
+python scripts/fake_data.py --force
+```
+
+### Training
+
+For single-device training, use:
+
+```bash
+python scripts/train.py
+```
+
+For multi-GPU training, use:
+
+```bash
+torchrun --nproc_per_node=8 --nnodes=1 --rdzv_id=1000 --rdzv_backend=c10d --max_restarts=0 --master_addr="localhost" --master_port=29500 scripts/train.py
+```
+
+Training jointly optimizes observational data loss and PDE residuals from the beta-plane horizontal momentum equations and the pressure-coordinate continuity equation. The checkpoint at `result/checkpoints/pinn-tc.pt` contains model state, model configuration, data format version, Adam/L-BFGS states, training history, and completed-step state.
+
+### Trained Weights
+
+This repository does not bundle official paper weights, and no confirmed official checkpoint is available in the public resources. Running the training script produces a local engineering checkpoint, but no compatibility with external weights or equivalence to paper training results is claimed.
+
+### Inference
+
+```bash
+python scripts/inference.py
+```
+
+Inference loads the trained checkpoint and generates complete three-dimensional wind, geopotential, and pressure-vertical-velocity fields from sparse tropical-cyclone observation coordinates and physical constraints. Complete numerical results are saved to `result/output/predictions.npz`.
+
+### Evaluation and Visualization
+
+```bash
+python scripts/result.py
+```
+
+Evaluation computes wind-speed RMSE, Pearson correlation, PDE residuals, and radial RMSE as aggregates over the three observation times; it saves neither per-class metrics nor separate per-time metrics. It also generates a four-panel comparison of sparse observations, target wind speed, predicted wind speed, and error. Synthetic-data results validate the engineering workflow only and do not represent formal paper performance.
+
+# Official OneScience Information
+
+| Platform | OneScience Main Repository | Skills Repository |
+| --- | --- | --- |
+| Gitee | https://gitee.com/onescience-ai/onescience | https://gitee.com/onescience-ai/oneskills |
+| GitHub | https://github.com/onescience-ai/OneScience | https://github.com/onescience-ai/oneskills |
+
+# Citation and License
+
+This repository is an independent engineering reproduction of the public PINN-TC specifications.
+
+Use of this repository's code, official model weights, and data remains subject to the licenses and terms of their respective projects.

conf/config.yaml

@@ -0,0 +1,60 @@
+seed: 42
+domain:
+  horizontal_extent_m: 400000
+  time_extent_s: 21600
+  pressure_min_pa: 15000
+  pressure_max_pa: 90000
+data:
+  root: data
+  format_version: pinn_tc_v1
+  input_order: [y, x, t, p]
+  supervised_order: [u, v, geopotential]
+  output_order: [u, v, geopotential, omega]
+  horizontal_extent_m: 400000
+  time_extent_s: 21600
+  pressure_min_pa: 15000
+  pressure_max_pa: 90000
+  grid_points: 251
+  pressure_levels: 44
+  time_steps: 3
+  observation_times_s: [-10800, 0, 10800]
+  observation_mode: Switch
+  observation_line_points: 17
+  observation_pressure_stride: 11
+  boundary_points_per_edge: 5
+  collocation_points: 512
+model:
+  width: 16
+  depth: 2
+physics:
+  f0: 0.00005
+  beta: 0.00000000002
+train:
+  optimizer_sequence: [Adam, L-BFGS]
+  adam_steps: 2
+  adam_lr: 0.001
+  lbfgs_steps: 1
+  lbfgs_lr: 0.5
+  lbfgs_history_size: 10
+  observation_batch: 128
+  collocation_batch: 128
+  wind_speed_weight: 1.0
+  gamma: 0.1
+  momentum_scale: 0.01
+  continuity_scale: 0.0001
+runtime:
+  device: auto
+  ddp: torchrun-compatible
+  inference_chunk_size: 131072
+evaluation:
+  pde_points: 256
+paths:
+  checkpoint: result/checkpoints/pinn-tc.pt
+  training_metrics: result/training/metrics.json
+  prediction: result/output/predictions.npz
+  evaluation_dir: result/evaluation
+paper_model:
+  2D: {hidden_layers: 8, width: 100}
+  3D: {hidden_layers: 4, width: 100}
+  real: {hidden_layers: 4, width: 50}
+  optimizer_sequence: [Adam, L-BFGS]

config.json

@@ -0,0 +1,55 @@
+{
+  "model_name": "PINN-TC",
+  "model_type": "pinn_tc",
+  "architectures": [
+    "PINNTC"
+  ],
+  "framework": "PyTorch",
+  "domain": "tropical cyclone dynamics",
+  "task": "continuous-domain tropical-cyclone field reconstruction",
+  "implementation": {
+    "entry_point": "model/pinn-tc.py",
+    "scope": "Core-method and full-data-dimension reduced-scale end-to-end engineering reproduction",
+    "train_script": "scripts/train.py",
+    "inference_script": "scripts/inference.py",
+    "evaluation_script": "scripts/result.py",
+    "synthetic_data_script": "scripts/fake_data.py"
+  },
+  "architecture": {
+    "inputs": ["y", "x", "t", "p"],
+    "outputs": ["u", "v", "geopotential", "omega"],
+    "supervised_outputs": ["u", "v", "geopotential"],
+    "pde_constrained_output": "omega",
+    "activation": "tanh",
+    "engineering_model": {"width": 16, "depth": 2},
+    "paper_models": {
+      "2D": {"hidden_layers": 8, "width": 100},
+      "3D": {"hidden_layers": 4, "width": 100},
+      "real": {"hidden_layers": 4, "width": 50}
+    },
+    "physics_residuals": ["zonal_momentum", "meridional_momentum", "pressure_continuity"]
+  },
+  "data": {
+    "dataset": "T-SHiELD Hurricane Ida case and NOAA/AOML observations",
+    "protocol": "Sparse storm-centered coordinate observations with collocation points over the full continuous domain",
+    "format_version": "pinn_tc_v1",
+    "input_shape": ["N", 4],
+    "output_shape": ["N", 4],
+    "dense_output_shape": [251, 251, 44, 3, 4],
+    "dense_output_layout": "YXPTV",
+    "horizontal_extent_m": [-400000, 400000],
+    "horizontal_resolution_km": 3.2,
+    "pressure_levels": 44,
+    "pressure_range_pa": [15000, 90000],
+    "observation_times_s": [-10800, 0, 10800],
+    "units": {"u_v": "m s-1", "geopotential": "m2 s-2", "omega": "Pa s-1"}
+  },
+  "configuration_sources": [
+    "conf/config.yaml",
+    "model/pinn-tc.py",
+    "scripts/fake_data.py",
+    "scripts/train.py",
+    "scripts/inference.py",
+    "scripts/result.py"
+  ]
+}

configuration.json

@@ -0,0 +1,12 @@
+{
+  "framework": "PyTorch",
+  "task": "continuous-domain tropical-cyclone field reconstruction",
+  "model": "PINN-TC",
+  "input_format": "coordinates [N,4] ordered as [y,x,t,p]",
+  "protocol": "Sparse-observation and PDE-constrained reconstruction on a 251x251x44x3 evaluation grid",
+  "default_config": "conf/config.yaml",
+  "training": "scripts/train.py",
+  "inference": "scripts/inference.py",
+  "evaluation": "scripts/result.py",
+  "visualization": "scripts/result.py"
+}

model/pinn-tc.py

@@ -0,0 +1,92 @@
+"""Continuous-coordinate physics-informed neural network for a tropical cyclone."""
+
+from __future__ import annotations
+
+import math
+
+import numpy as np
+import torch
+from torch import nn
+
+
+INPUT_ORDER = ("y", "x", "t", "p")
+OUTPUT_ORDER = ("u", "v", "geopotential", "omega")
+
+
+class PINNTC(nn.Module):
+    """Tanh MLP mapping physical [y, x, t, p] coordinates to four fields."""
+
+    def __init__(self, width: int = 16, depth: int = 2):
+        super().__init__()
+        if width < 1 or depth < 1:
+            raise ValueError("width and depth must be positive")
+        layers: list[nn.Module] = []
+        dimensions = [4] + [width] * depth + [4]
+        for index, (source, target) in enumerate(zip(dimensions[:-1], dimensions[1:])):
+            layer = nn.Linear(source, target)
+            nn.init.xavier_normal_(layer.weight, gain=nn.init.calculate_gain("tanh") if index < depth else 1.0)
+            nn.init.zeros_(layer.bias)
+            layers.append(layer)
+            if index < depth:
+                layers.append(nn.Tanh())
+        self.network = nn.Sequential(*layers)
+        self.register_buffer("coordinate_scale", torch.tensor([400_000.0, 400_000.0, 21_600.0, 75_000.0]))
+        self.register_buffer("coordinate_center", torch.tensor([0.0, 0.0, 0.0, 52_500.0]))
+        self.register_buffer("output_scale", torch.tensor([45.0, 45.0, 6_000.0, 0.5]))
+        self.register_buffer("output_center", torch.tensor([0.0, 0.0, 50_000.0, 0.0]))
+
+    def forward(self, coordinates: torch.Tensor) -> torch.Tensor:
+        if coordinates.ndim != 2 or coordinates.shape[1] != 4:
+            raise ValueError("coordinates must have shape [N,4] in [y,x,t,p] order")
+        normalized = (coordinates - self.coordinate_center) / self.coordinate_scale
+        return self.network(normalized) * self.output_scale + self.output_center
+
+
+def pde_residuals(
+    model: nn.Module,
+    coordinates: torch.Tensor,
+    f0: float = 5.0e-5,
+    beta: float = 2.0e-11,
+) -> tuple[torch.Tensor, torch.Tensor, torch.Tensor]:
+    """Return beta-plane pressure-coordinate momentum and continuity residuals."""
+    if not coordinates.requires_grad:
+        coordinates = coordinates.detach().requires_grad_(True)
+    fields = model(coordinates)
+    u, v, phi, omega = (fields[:, index] for index in range(4))
+
+    def gradient(value: torch.Tensor) -> torch.Tensor:
+        return torch.autograd.grad(value, coordinates, torch.ones_like(value), create_graph=True, retain_graph=True)[0]
+
+    grad_u, grad_v, grad_phi, grad_omega = map(gradient, (u, v, phi, omega))
+    y = coordinates[:, 0]
+    coriolis = f0 + beta * y
+    # Coordinate columns are [y, x, t, p]. Geopotential has units m2 s-2.
+    ru = grad_u[:, 2] + u * grad_u[:, 1] + v * grad_u[:, 0] + omega * grad_u[:, 3] - coriolis * v + grad_phi[:, 1]
+    rv = grad_v[:, 2] + u * grad_v[:, 1] + v * grad_v[:, 0] + omega * grad_v[:, 3] + coriolis * u + grad_phi[:, 0]
+    rc = grad_u[:, 1] + grad_v[:, 0] + grad_omega[:, 3]
+    return ru, rv, rc
+
+
+def analytic_vortex_numpy(coordinates: np.ndarray) -> np.ndarray:
+    """Structured, noise-free translating and pressure-decaying analytic cyclone."""
+    coordinates = np.asarray(coordinates, dtype=np.float64)
+    if coordinates.ndim != 2 or coordinates.shape[1] != 4:
+        raise ValueError("coordinates must have shape [N,4]")
+    y, x, time, pressure = coordinates.T
+    center_x = 45_000.0 * np.sin(time / 21_600.0 * np.pi / 2.0)
+    center_y = 30_000.0 * np.cos(time / 21_600.0 * np.pi / 2.0)
+    dx, dy = x - center_x, y - center_y
+    radius = np.sqrt(dx * dx + dy * dy + 4.0e6)
+    radial_scale = 95_000.0
+    vertical = np.exp(-((pressure - 55_000.0) / 27_000.0) ** 2)
+    tangential = 42.0 * (radius / radial_scale) * np.exp(1.0 - radius / radial_scale) * vertical
+    radial = -3.0 * np.exp(-radius / 180_000.0) * vertical
+    u = -tangential * dy / radius + radial * dx / radius + 2.5
+    v = tangential * dx / radius + radial * dy / radius + 1.0
+    geopotential = 50_000.0 - 5_500.0 * np.exp(-(radius / 130_000.0) ** 2) * vertical
+    omega = -0.32 * np.exp(-(radius / 110_000.0) ** 2) * np.sin(np.pi * (pressure - 15_000.0) / 75_000.0)
+    return np.stack((u, v, geopotential, omega), axis=1).astype(np.float32)
+
+
+def parameter_count(model: nn.Module) -> int:
+    return sum(parameter.numel() for parameter in model.parameters())

scripts/fake_data.py

@@ -0,0 +1,82 @@
+"""Generate sparse structured observations and collocation coordinates lazily."""
+
+import argparse
+import importlib.util
+from pathlib import Path
+
+import numpy as np
+import yaml
+
+
+ROOT = Path(__file__).resolve().parents[1]
+
+
+def load_model_module():
+    spec = importlib.util.spec_from_file_location("pinn_tc_model", ROOT / "model/pinn-tc.py")
+    module = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(module)
+    return module
+
+
+def structured_observations(config, mode):
+    data = config["data"]
+    extent = float(data["horizontal_extent_m"])
+    levels = np.linspace(float(data["pressure_min_pa"]), float(data["pressure_max_pa"]), int(data["pressure_levels"]))
+    times = np.asarray(data["observation_times_s"], dtype=np.float64)
+    if times.shape != (int(data["time_steps"]),) or not np.array_equal(times, [-10800.0, 0.0, 10800.0]):
+        raise ValueError("paper observation times must be [-3h, 0h, +3h]")
+    line = np.linspace(-extent, extent, int(data["observation_line_points"]))
+    coordinates = []
+    for time_index, time in enumerate(times):
+        active_mode = ("Cross" if time_index % 2 == 0 else "Plus") if mode == "Switch" else mode
+        for pressure in levels[:: int(data["observation_pressure_stride"])]:
+            if active_mode == "Plus":
+                coordinates.extend((y, x, time, pressure) for y, x in zip(np.zeros_like(line), line))
+                coordinates.extend((y, x, time, pressure) for y, x in zip(line, np.zeros_like(line)))
+            elif active_mode == "Cross":
+                coordinates.extend((y, x, time, pressure) for y, x in zip(line, line))
+                coordinates.extend((y, x, time, pressure) for y, x in zip(line, -line))
+            else:
+                raise ValueError("observation mode must be Cross, Plus, or Switch")
+            edge = np.linspace(-extent, extent, int(data["boundary_points_per_edge"]))
+            coordinates.extend((-extent, value, time, pressure) for value in edge)
+            coordinates.extend((extent, value, time, pressure) for value in edge)
+            coordinates.extend((value, -extent, time, pressure) for value in edge[1:-1])
+            coordinates.extend((value, extent, time, pressure) for value in edge[1:-1])
+    return np.unique(np.asarray(coordinates, dtype=np.float32), axis=0)
+
+
+def main():
+    parser = argparse.ArgumentParser()
+    parser.add_argument("--mode", choices=("Cross", "Plus", "Switch"))
+    parser.add_argument("--force", action="store_true")
+    args = parser.parse_args()
+    config = yaml.safe_load((ROOT / "conf/config.yaml").read_text())
+    mode = args.mode or config["data"]["observation_mode"]
+    output = ROOT / config["data"]["root"]
+    output.mkdir(parents=True, exist_ok=True)
+    path = output / "training_points.npz"
+    if path.exists() and not args.force:
+        print(f"exists={path.relative_to(ROOT)} (use --force to regenerate)")
+        return
+    module = load_model_module()
+    observations = structured_observations(config, mode)
+    targets = module.analytic_vortex_numpy(observations)[:, :3]
+    rng = np.random.default_rng(int(config["seed"]))
+    count = int(config["data"]["collocation_points"])
+    collocation = np.column_stack((
+        rng.uniform(-config["data"]["horizontal_extent_m"], config["data"]["horizontal_extent_m"], (count, 2)),
+        rng.uniform(-config["data"]["time_extent_s"], config["data"]["time_extent_s"], count),
+        rng.uniform(config["data"]["pressure_min_pa"], config["data"]["pressure_max_pa"], count),
+    )).astype(np.float32)
+    # The first two random columns already follow the required [y, x] order.
+    np.savez(path, observation_coordinates=observations, observation_targets=targets,
+             collocation_coordinates=collocation, input_order=np.asarray(module.INPUT_ORDER),
+             supervised_order=np.asarray(module.OUTPUT_ORDER[:3]), observation_mode=np.asarray(mode),
+             dense_shape=np.asarray([data_dim := int(config["data"]["grid_points"]), data_dim,
+                                     int(config["data"]["pressure_levels"]), int(config["data"]["time_steps"])]))
+    print(f"generated={path.relative_to(ROOT)} observations={len(observations)} collocation={len(collocation)} mode={mode}")
+
+
+if __name__ == "__main__":
+    main()

scripts/inference.py

@@ -0,0 +1,77 @@
+"""Run chunked inference on the full 251 x 251 x 44 x 3 paper-resolution domain."""
+
+import importlib.util
+from pathlib import Path
+
+import numpy as np
+import torch
+import yaml
+
+
+ROOT = Path(__file__).resolve().parents[1]
+
+
+def load_model_module():
+    spec = importlib.util.spec_from_file_location("pinn_tc_model", ROOT / "model/pinn-tc.py")
+    module = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(module)
+    return module
+
+
+def coordinates_for_indices(indices, y_axis, x_axis, times, pressures):
+    x_count, p_count, t_count = len(x_axis), len(pressures), len(times)
+    y_index = indices // (x_count * p_count * t_count)
+    remainder = indices % (x_count * p_count * t_count)
+    x_index = remainder // (p_count * t_count)
+    remainder %= p_count * t_count
+    p_index = remainder // t_count
+    t_index = remainder % t_count
+    return np.column_stack((y_axis[y_index], x_axis[x_index], times[t_index], pressures[p_index])).astype(np.float32)
+
+
+def main():
+    config = yaml.safe_load((ROOT / "conf/config.yaml").read_text())
+    module = load_model_module()
+    requested = config["runtime"]["device"]
+    device = torch.device("cuda" if requested == "auto" and torch.cuda.is_available() else ("cpu" if requested == "auto" else requested))
+    checkpoint = torch.load(ROOT / config["paths"]["checkpoint"], map_location=device, weights_only=False)
+    if checkpoint.get("format_version") != config["data"]["format_version"]:
+        raise ValueError("checkpoint format_version mismatch")
+    if checkpoint["input_order"] != list(module.INPUT_ORDER) or checkpoint["output_order"] != list(module.OUTPUT_ORDER):
+        raise ValueError("checkpoint coordinate or variable protocol mismatch")
+    model = module.PINNTC(**checkpoint["model_config"]).to(device)
+    model.load_state_dict(checkpoint["model_state_dict"]); model.eval()
+    data = config["data"]
+    y_axis = np.linspace(-data["horizontal_extent_m"], data["horizontal_extent_m"], data["grid_points"], dtype=np.float32)
+    x_axis = y_axis.copy()
+    times = np.asarray(data["observation_times_s"], dtype=np.float32)
+    pressures = np.linspace(data["pressure_min_pa"], data["pressure_max_pa"], data["pressure_levels"], dtype=np.float32)
+    shape = (len(y_axis), len(x_axis), len(pressures), len(times), 4)
+    expected = (251, 251, 44, 3, 4)
+    if shape != expected:
+        raise ValueError(f"dense grid must remain {expected}, got {shape}")
+    output = ROOT / config["paths"]["prediction"]
+    output.parent.mkdir(parents=True, exist_ok=True)
+    temporary = output.with_suffix(".work.npy")
+    prediction = np.lib.format.open_memmap(temporary, mode="w+", dtype=np.float32, shape=shape)
+    flat = prediction.reshape(-1, 4)
+    chunk_size = int(config["runtime"]["inference_chunk_size"])
+    with torch.inference_mode():
+        for start in range(0, len(flat), chunk_size):
+            stop = min(start + chunk_size, len(flat))
+            coordinates = coordinates_for_indices(np.arange(start, stop, dtype=np.int64), y_axis, x_axis, times, pressures)
+            values = model(torch.from_numpy(coordinates).to(device)).cpu().numpy().astype(np.float32)
+            if not np.isfinite(values).all():
+                raise FloatingPointError(f"non-finite prediction in rows {start}:{stop}")
+            flat[start:stop] = values
+    prediction.flush()
+    np.savez(output, predictions=prediction, y=y_axis, x=x_axis, time=times, pressure=pressures,
+             input_order=np.asarray(module.INPUT_ORDER), output_order=np.asarray(module.OUTPUT_ORDER),
+             format_version=np.asarray(config["data"]["format_version"]), layout=np.asarray("YXPTV"))
+    del prediction
+    temporary.unlink()
+    print(f"predictions={output.relative_to(ROOT)} shape={shape} bytes={output.stat().st_size}")
+
+
+if __name__ == "__main__":
+    main()

scripts/result.py

@@ -0,0 +1,108 @@
+"""Evaluate dense PINN-TC predictions and create sparse/target/prediction/error panels."""
+
+import importlib.util
+import json
+from pathlib import Path
+
+import matplotlib
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+import numpy as np
+import torch
+import yaml
+
+
+ROOT = Path(__file__).resolve().parents[1]
+
+
+def load_model_module():
+    spec = importlib.util.spec_from_file_location("pinn_tc_model", ROOT / "model/pinn-tc.py")
+    module = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(module)
+    return module
+
+
+def main():
+    config = yaml.safe_load((ROOT / "conf/config.yaml").read_text())
+    module = load_model_module()
+    archive = np.load(ROOT / config["paths"]["prediction"])
+    if str(archive["format_version"]) != config["data"]["format_version"]:
+        raise ValueError("prediction format_version mismatch")
+    prediction = archive["predictions"]
+    expected = (251, 251, 44, 3, 4)
+    if prediction.shape != expected or prediction.dtype != np.float32:
+        raise ValueError(f"prediction must be float32 with shape {expected}")
+    y_axis, x_axis, pressures, times = archive["y"], archive["x"], archive["pressure"], archive["time"]
+    sum_speed_error = sum_pred = sum_true = sum_pred2 = sum_true2 = sum_cross = 0.0
+    count = 0
+    radial_edges = np.linspace(0.0, np.sqrt(2.0) * config["data"]["horizontal_extent_m"], 13)
+    radial_sse = np.zeros(len(radial_edges) - 1); radial_count = np.zeros(len(radial_edges) - 1, dtype=np.int64)
+    yy, xx = np.meshgrid(y_axis, x_axis, indexing="ij")
+    radius = np.sqrt(xx * xx + yy * yy)
+    for pressure_index, pressure in enumerate(pressures):
+        for time_index, time in enumerate(times):
+            coordinates = np.column_stack((yy.ravel(), xx.ravel(), np.full(yy.size, time), np.full(yy.size, pressure))).astype(np.float32)
+            target = module.analytic_vortex_numpy(coordinates).reshape(251, 251, 4)
+            predicted_speed = np.linalg.norm(prediction[:, :, pressure_index, time_index, :2].astype(np.float64), axis=-1)
+            true_speed = np.linalg.norm(target[:, :, :2].astype(np.float64), axis=-1)
+            difference = predicted_speed - true_speed
+            sum_speed_error += np.square(difference).sum(); count += difference.size
+            sum_pred += predicted_speed.sum(); sum_true += true_speed.sum()
+            sum_pred2 += np.square(predicted_speed).sum(); sum_true2 += np.square(true_speed).sum()
+            sum_cross += (predicted_speed * true_speed).sum()
+            bins = np.clip(np.digitize(radius.ravel(), radial_edges) - 1, 0, len(radial_sse) - 1)
+            radial_sse += np.bincount(bins, weights=np.square(difference).ravel(), minlength=len(radial_sse))
+            radial_count += np.bincount(bins, minlength=len(radial_count))
+    covariance = sum_cross - sum_pred * sum_true / count
+    variance_pred = sum_pred2 - sum_pred * sum_pred / count
+    variance_true = sum_true2 - sum_true * sum_true / count
+    pearson = covariance / np.sqrt(max(variance_pred * variance_true, 1e-30))
+
+    checkpoint = torch.load(ROOT / config["paths"]["checkpoint"], map_location="cpu", weights_only=False)
+    if checkpoint.get("format_version") != config["data"]["format_version"]:
+        raise ValueError("checkpoint format_version mismatch")
+    model = module.PINNTC(**checkpoint["model_config"]); model.load_state_dict(checkpoint["model"]); model.eval()
+    rng = np.random.default_rng(int(config["seed"]) + 99)
+    residual_count = int(config["evaluation"]["pde_points"])
+    residual_coordinates = np.column_stack((rng.uniform(-400_000, 400_000, (residual_count, 2)),
+                                            rng.uniform(-21_600, 21_600, residual_count),
+                                            rng.uniform(15_000, 90_000, residual_count))).astype(np.float32)
+    tensor = torch.from_numpy(residual_coordinates).requires_grad_(True)
+    ru, rv, rc = module.pde_residuals(model, tensor, **config["physics"])
+    report = {"dense_shape": list(expected), "dense_points": int(np.prod(expected[:-1])),
+              "wind_speed_rmse_m_s": float(np.sqrt(sum_speed_error / count)), "wind_speed_pearson": float(pearson),
+              "pde_residual_rms": {"ru_m_s2": float(torch.sqrt(torch.mean(ru.square())).detach()),
+                                   "rv_m_s2": float(torch.sqrt(torch.mean(rv.square())).detach()),
+                                   "rc_s-1": float(torch.sqrt(torch.mean(rc.square())).detach())},
+              "radial_rmse": [{"radius_mid_km": float((left + right) / 2000.0), "rmse_m_s": float(np.sqrt(error / max(samples, 1))), "count": int(samples)}
+                              for left, right, error, samples in zip(radial_edges[:-1], radial_edges[1:], radial_sse, radial_count)],
+              "scope": "Structured analytic engineering validation, not a paper-performance claim."}
+    values = [report["wind_speed_rmse_m_s"], report["wind_speed_pearson"], *report["pde_residual_rms"].values()]
+    if not np.isfinite(values).all():
+        raise FloatingPointError("evaluation produced a non-finite metric")
+    output = ROOT / config["paths"]["evaluation_dir"]
+    output.mkdir(parents=True, exist_ok=True)
+    (output / "metrics.json").write_text(json.dumps(report, indent=2) + "\n")
+
+    pressure_index, time_index = int(np.argmin(abs(pressures - 50_000))), int(np.argmin(abs(times)))
+    coordinates = np.column_stack((yy.ravel(), xx.ravel(), np.full(yy.size, times[time_index]),
+                                   np.full(yy.size, pressures[pressure_index]))).astype(np.float32)
+    target_speed = np.linalg.norm(module.analytic_vortex_numpy(coordinates)[:, :2], axis=1).reshape(251, 251)
+    predicted_speed = np.linalg.norm(prediction[:, :, pressure_index, time_index, :2], axis=-1)
+    sparse = np.load(ROOT / config["data"]["root"] / "training_points.npz")["observation_coordinates"]
+    selected = sparse[(abs(sparse[:, 2] - times[time_index]) < 1) & (abs(sparse[:, 3] - pressures[pressure_index]) < 1000)]
+    figure, axes_plot = plt.subplots(1, 4, figsize=(16, 4), constrained_layout=True)
+    axes_plot[0].scatter(selected[:, 1] / 1000, selected[:, 0] / 1000, s=5, c="black"); axes_plot[0].set_title("Sparse observations")
+    panels = ((target_speed, "Target wind speed", "viridis"), (predicted_speed, "Predicted wind speed", "viridis"),
+              (predicted_speed - target_speed, "Prediction error", "coolwarm"))
+    for axis, (field, title, color_map) in zip(axes_plot[1:], panels):
+        image = axis.imshow(field, origin="lower", extent=(-400, 400, -400, 400), cmap=color_map)
+        axis.set_title(title); figure.colorbar(image, ax=axis, shrink=0.75)
+    for axis in axes_plot:
+        axis.set_xlabel("x (km)"); axis.set_ylabel("y (km)")
+    figure.savefig(output / "comparison.png", dpi=140); plt.close(figure)
+    print(f"evaluation={output.relative_to(ROOT)} wind_rmse={report['wind_speed_rmse_m_s']:.6g} pearson={pearson:.6g}")
+
+
+if __name__ == "__main__":
+    main()

scripts/train.py

@@ -0,0 +1,137 @@
+"""Train PINN-TC with Adam and an optional L-BFGS refinement under torchrun DDP."""
+
+import argparse
+import importlib.util
+import json
+import os
+import random
+from pathlib import Path
+
+import numpy as np
+import torch
+import yaml
+from torch.nn.parallel import DistributedDataParallel
+
+
+ROOT = Path(__file__).resolve().parents[1]
+
+
+def load_model_module():
+    spec = importlib.util.spec_from_file_location("pinn_tc_model", ROOT / "model/pinn-tc.py")
+    module = importlib.util.module_from_spec(spec)
+    spec.loader.exec_module(module)
+    return module
+
+
+def choose_device(config, local_rank):
+    requested = config["runtime"]["device"]
+    if requested == "auto":
+        return torch.device("cuda", local_rank) if torch.cuda.is_available() else torch.device("cpu")
+    return torch.device(requested)
+
+
+def main():
+    parser = argparse.ArgumentParser()
+    parser.add_argument("--resume", action="store_true", help="resume optimizer and model state from the configured checkpoint")
+    args = parser.parse_args()
+    config = yaml.safe_load((ROOT / "conf/config.yaml").read_text())
+    module = load_model_module()
+    seed = int(config["seed"])
+    random.seed(seed); np.random.seed(seed); torch.manual_seed(seed)
+    distributed = int(os.environ.get("WORLD_SIZE", "1")) > 1
+    local_rank = int(os.environ.get("LOCAL_RANK", "0"))
+    if distributed:
+        torch.distributed.init_process_group("nccl" if torch.cuda.is_available() else "gloo")
+    rank = torch.distributed.get_rank() if distributed else 0
+    device = choose_device(config, local_rank)
+    if device.type == "cuda":
+        torch.cuda.set_device(device); torch.cuda.manual_seed_all(seed)
+
+    source = np.load(ROOT / config["data"]["root"] / "training_points.npz")
+    observations = torch.from_numpy(source["observation_coordinates"]).to(device)
+    targets = torch.from_numpy(source["observation_targets"]).to(device)
+    collocation = torch.from_numpy(source["collocation_coordinates"]).to(device)
+    model = module.PINNTC(**config["model"]).to(device)
+    wrapped = DistributedDataParallel(model, device_ids=[local_rank] if device.type == "cuda" else None) if distributed else model
+    bare = wrapped.module if distributed else wrapped
+    gamma = float(config["train"]["gamma"])
+    batch_observations = int(config["train"]["observation_batch"])
+    batch_collocation = int(config["train"]["collocation_batch"])
+
+    def loss_for(obs_index, col_index):
+        prediction = wrapped(observations[obs_index])[:, :3]
+        speed = torch.sqrt(targets[obs_index, 0] ** 2 + targets[obs_index, 1] ** 2)
+        wind_weight = 1.0 + float(config["train"]["wind_speed_weight"]) * speed / 50.0
+        normalized_error = (prediction - targets[obs_index]) / torch.tensor([45.0, 45.0, 6000.0], device=device)
+        data_loss = (normalized_error.square() * wind_weight[:, None]).mean()
+        points = collocation[col_index].detach().requires_grad_(True)
+        ru, rv, rc = module.pde_residuals(wrapped, points, **config["physics"])
+        pde_loss = (ru.square().mean() + rv.square().mean()) / float(config["train"]["momentum_scale"]) ** 2
+        pde_loss = pde_loss + rc.square().mean() / float(config["train"]["continuity_scale"]) ** 2
+        return data_loss + gamma * pde_loss, data_loss, pde_loss
+
+    optimizer = torch.optim.Adam(wrapped.parameters(), lr=float(config["train"]["adam_lr"]))
+    history = []
+    checkpoint_path = ROOT / config["paths"]["checkpoint"]
+    completed_adam_steps = 0
+    checkpoint = None
+    if args.resume:
+        checkpoint = torch.load(checkpoint_path, map_location=device, weights_only=False)
+        if checkpoint["model_config"] != config["model"]:
+            raise ValueError("checkpoint model configuration mismatch")
+        bare.load_state_dict(checkpoint["model_state_dict"])
+        optimizer.load_state_dict(checkpoint["adam_optimizer_state_dict"])
+        history = list(checkpoint.get("history", []))
+        completed_adam_steps = int(checkpoint.get("completed_adam_steps", 0))
+    generator = torch.Generator(device="cpu").manual_seed(seed + rank)
+    for step in range(completed_adam_steps, completed_adam_steps + int(config["train"]["adam_steps"])):
+        obs_index = torch.randperm(len(observations), generator=generator)[:batch_observations].to(device)
+        col_index = torch.randperm(len(collocation), generator=generator)[:batch_collocation].to(device)
+        optimizer.zero_grad(set_to_none=True)
+        loss, data_loss, pde_loss = loss_for(obs_index, col_index)
+        loss.backward(); optimizer.step()
+        history.append({"stage": "Adam", "step": step + 1, "loss": float(loss.detach()),
+                        "data_loss": float(data_loss.detach()), "pde_loss": float(pde_loss.detach())})
+
+    lbfgs_steps = int(config["train"]["lbfgs_steps"])
+    optimizer_lbfgs = None
+    if lbfgs_steps > 0:
+        # Every DDP rank runs the same closure count and deterministic subset so gradient collectives remain aligned.
+        obs_index = torch.arange(min(batch_observations, len(observations)), device=device)
+        col_index = torch.arange(min(batch_collocation, len(collocation)), device=device)
+        optimizer_lbfgs = torch.optim.LBFGS(wrapped.parameters(), lr=float(config["train"]["lbfgs_lr"]),
+                                             max_iter=lbfgs_steps, history_size=int(config["train"]["lbfgs_history_size"]),
+                                             line_search_fn=None)
+        if checkpoint is not None and checkpoint.get("lbfgs_optimizer_state_dict") is not None:
+            optimizer_lbfgs.load_state_dict(checkpoint["lbfgs_optimizer_state_dict"])
+        latest = {}
+
+        def closure():
+            optimizer_lbfgs.zero_grad(set_to_none=True)
+            loss, data_loss, pde_loss = loss_for(obs_index, col_index)
+            loss.backward()
+            latest.update(loss=float(loss.detach()), data_loss=float(data_loss.detach()), pde_loss=float(pde_loss.detach()))
+            return loss
+
+        optimizer_lbfgs.step(closure)
+        history.append({"stage": "L-BFGS", "step": 1, **latest})
+
+    if rank == 0:
+        metrics_path = ROOT / config["paths"]["training_metrics"]
+        checkpoint_path.parent.mkdir(parents=True, exist_ok=True); metrics_path.parent.mkdir(parents=True, exist_ok=True)
+        torch.save({"model": bare.state_dict(), "model_state_dict": bare.state_dict(),
+                    "model_config": config["model"], "format_version": config["data"]["format_version"],
+                    "adam_optimizer_state_dict": optimizer.state_dict(),
+                    "lbfgs_optimizer_state_dict": optimizer_lbfgs.state_dict() if optimizer_lbfgs is not None else None,
+                    "completed_adam_steps": completed_adam_steps + int(config["train"]["adam_steps"]),
+                    "history": history,
+                    "input_order": list(module.INPUT_ORDER), "output_order": list(module.OUTPUT_ORDER),
+                    "domain": config["domain"], "physics": config["physics"], "seed": seed}, checkpoint_path)
+        metrics_path.write_text(json.dumps({"history": history, "distributed_world_size": int(os.environ.get("WORLD_SIZE", "1"))}, indent=2) + "\n")
+        print(f"checkpoint={checkpoint_path.relative_to(ROOT)} parameters={module.parameter_count(bare)} final_loss={history[-1]['loss']:.6g}")
+    if distributed:
+        torch.distributed.destroy_process_group()
+
+
+if __name__ == "__main__":
+    main()