diff --git a/benchmarks/README.md b/benchmarks/README.md
new file mode 100644
index 0000000..b92d0e7
--- /dev/null
+++ b/benchmarks/README.md
@@ -0,0 +1,89 @@
+# LyoPRONTO Benchmarks
+
+Tools to compare Pyomo (finite differences + orthogonal collocation) vs Scipy optimizers across parameter grids.
+
+## Current Workflow (Jan 2025)
+
+### 1. Generate Grid Data
+Use `grid_cli.py` to run Cartesian product benchmarks across N parameters:
+
+```bash
+# Generate 3×3 grid: scipy baseline + FD + collocation
+python benchmarks/grid_cli.py generate \
+ --task Tsh --scenario baseline \
+ --vary product.A1=16,18,20 \
+ --vary ht.KC=2.75e-4,3.3e-4,4.0e-4 \
+ --methods scipy,fd,colloc \
+ --n-elements 24 --n-collocation 3 \
+ --out benchmarks/results/grid_A1_KC.jsonl
+```
+
+**Options:**
+- `--vary key=val1,val2,...` (repeatable) — parameter sweeps
+- `--methods scipy,fd,colloc` — which solvers to run
+- `--n-elements N` — finite elements (both FD and collocation base)
+- `--n-collocation NCP` — collocation points per element
+- `--warmstart` — enable staged solve with scipy trajectory (off by default for robustness)
+- `--raw-colloc` — disable effective-nfe parity reporting
+- `--force` — regenerate even if output exists (reuse-first by default)
+
+**Make shortcut:**
+```bash
+make bench VARY='product.A1=16,18,20 ht.KC=2.75e-4,3.3e-4,4.0e-4' METHODS=fd,colloc
+```
+
+### 2. Analyze Results
+Open `benchmarks/grid_analysis.ipynb` (read-only; expects pre-generated JSONL):
+
+```bash
+JSONL_PATH=benchmarks/results/grid_A1_KC.jsonl jupyter notebook benchmarks/grid_analysis.ipynb
+```
+
+Notebook cells produce:
+- CSV pivot tables (objectives, ratios, speedups)
+- Heatmaps (objective difference, speedup, parity)
+- Scatter plots and histograms
+- Summary interpretation
+
+**Headless execution (optional):**
+```bash
+make analyze OUT=benchmarks/results/grid_A1_KC.jsonl METRIC=ratio.pyomo_over_scipy
+```
+
+## Core Modules
+
+- **`grid_cli.py`** — CLI for N-dimensional grid generation (replaces legacy `run_grid*.py`)
+- **`adapters.py`** — Normalized scipy/Pyomo runners with discretization metadata
+- **`scenarios.py`** — Scenario definitions (vial, product, ht, eq_cap, nVial)
+- **`schema.py`** — Versioned record serialization (v2: trajectories + hashing)
+- **`validate.py`** — Physics checks (mass balance, dryness, constraint violations)
+- **`grid_analysis.ipynb`** — Analysis notebook (read-only data loader)
+- **`results/`** — Default output directory
+
+## Tasks
+
+- `Tsh` — optimize shelf temperature only (pressure fixed)
+- `Pch` — optimize chamber pressure only (shelf fixed)
+- `both` — joint optimization (pressure + temperature)
+
+## Output Schema (v2)
+
+Each JSONL record includes:
+- `version`: schema version (currently 2)
+- `hash.inputs`, `hash.record`: SHA-256 hashes for deduplication
+- `environment`: Python, Pyomo, Ipopt versions, OS, hostname, timestamp
+- `task`, `scenario`: optimization variant and scenario name
+- `grid.param1`, `grid.param2`, ... : swept parameters with paths and values
+- `scipy`: `{success, wall_time_s, objective_time_hr, solver, metrics}`
+- `pyomo`: same as scipy, plus:
+ - `discretization`: `{method, n_elements_requested, n_elements_applied, n_collocation, effective_nfe, total_mesh_points}`
+ - `warmstart_used`: bool
+- `failed`: overall failure flag (any solver failed or dryness unmet)
+
+## Notes
+
+- **Warmstart disabled by default** for robustness testing; enable with `--warmstart`.
+- **Effective-nfe true by default** for collocation parity with FD mesh density.
+- **Reuse-first**: if JSONL exists, generation skipped unless `--force` supplied.
+- **Trajectories embedded** in records (numpy arrays → lists during serialization).
+- **Hashing** prevents duplicate runs (schema v2 `hash.inputs` field).
diff --git a/benchmarks/adapters.py b/benchmarks/adapters.py
new file mode 100644
index 0000000..c8f476e
--- /dev/null
+++ b/benchmarks/adapters.py
@@ -0,0 +1,244 @@
+"""Adapters normalizing scipy and Pyomo optimizer outputs.
+
+Each adapter returns a dictionary with standardized keys:
+- trajectory: np.ndarray (time, Tsub, Tbot, Tsh, Pch_mTorr, flux, percent_dried)
+- success: bool
+- message: str
+- raw: original solver output or model reference
+- solver_stats: dict (iterations, evals, etc.)
+"""
+
+from __future__ import annotations
+
+import time
+from typing import Any
+
+import numpy as np
+from lyopronto import opt_Pch, opt_Pch_Tsh, opt_Tsh
+from lyopronto.pyomo_models import optimizers as pyomo_opt
+
+DRYNESS_TARGET = 98.9 # Percentage (0-100)
+
+# Scipy adapters -------------------------------------------------------------
+
+
+def scipy_adapter(
+ task: str,
+ vial: dict[str, float],
+ product: dict[str, float],
+ ht: dict[str, float],
+ eq_cap: dict[str, float],
+ nVial: int,
+ scenario: dict[str, Any],
+ dt: float = 0.01,
+) -> dict[str, Any]:
+ """Run scipy baseline optimizer variant for specified task.
+
+ task 'Tsh': optimize shelf temperature only (pressure schedule fixed)
+ task 'Pch': optimize chamber pressure only (shelf temperature schedule fixed)
+ task 'both': optimize both pressure and temperature concurrently
+ """
+ if task == "Tsh":
+ # Pressure schedule fixed at 0.1 Torr with long hold allowing completion
+ Pchamber = {"setpt": [0.1], "dt_setpt": [1800.0], "ramp_rate": 0.5}
+ # Shelf temperature optimization bounds
+ Tshelf = {"min": -45.0, "max": 120.0}
+ runner = opt_Tsh.dry
+ args = (vial, product, ht, Pchamber, Tshelf, dt, eq_cap, nVial)
+ elif task == "Pch":
+ # Pressure optimization lower bound
+ Pchamber = {"min": 0.05}
+ # Shelf multi-step schedule with sufficient time for drying completion
+ # NOTE: dt_setpt in MINUTES (opt_Pch expects minutes, converts internally)
+ # High resistance products (A1=20) need ~86 hours, use 100 hours for margin
+ Tshelf = {
+ "init": -35.0,
+ "setpt": [-20.0, 120.0],
+ "dt_setpt": [300.0, 5700.0],
+ "ramp_rate": 1.0,
+ }
+ runner = opt_Pch.dry
+ args = (vial, product, ht, Pchamber, Tshelf, dt, eq_cap, nVial)
+ elif task == "both":
+ Pchamber = {"min": 0.05}
+ Tshelf = {"min": -45.0, "max": 120.0}
+ runner = opt_Pch_Tsh.dry
+ args = (vial, product, ht, Pchamber, Tshelf, dt, eq_cap, nVial)
+ else:
+ raise ValueError(f"Unknown task '{task}'")
+
+ start = time.perf_counter()
+ out = runner(*args)
+ wall = time.perf_counter() - start
+
+ traj = out
+ success = traj.size > 0
+ message = "scipy run completed"
+ objective_time_hr = float(traj[-1, 0]) if success else None
+
+ return {
+ "trajectory": traj,
+ "success": success,
+ "message": message,
+ "wall_time_s": wall,
+ "objective_time_hr": objective_time_hr,
+ "solver": {"status": "n/a", "termination_condition": "n/a"},
+ "solver_stats": {},
+ "raw": out,
+ }
+
+
+# Pyomo adapters -------------------------------------------------------------
+
+
+def pyomo_adapter(
+ task: str,
+ vial: dict[str, float],
+ product: dict[str, float],
+ ht: dict[str, float],
+ eq_cap: dict[str, float],
+ nVial: int,
+ scenario: dict[str, Any],
+ dt: float = 0.01,
+ warmstart: bool = False,
+ method: str = "fd", # 'fd' or 'colloc'
+ n_elements: int = 24,
+ n_collocation: int = 3,
+ effective_nfe: bool = True,
+ tsh_ramp_rate: float | None = None, # Max Tsh ramp rate [°C/hr]
+ pch_ramp_rate: float | None = None, # Max Pch ramp rate [Torr/hr]
+ use_secant_ramp_constraints: bool = True, # Add secant slope constraints for collocation
+) -> dict[str, Any]:
+ """Run Pyomo optimizer counterpart for specified task with discretization controls.
+
+ Parameters
+ ----------
+ method : str
+ 'fd' for finite differences, 'colloc' for orthogonal collocation.
+ n_elements : int
+ Base number of finite elements requested.
+ n_collocation : int
+ Number of collocation points per element (only if method == 'colloc').
+ effective_nfe : bool
+ If True (collocation only), treat n_elements as effective (parity with FD) for reporting.
+ warmstart : bool
+ Use staged solve with scipy warmstart trajectory. Default False for robustness benchmarking.
+ use_secant_ramp_constraints : bool
+ When using collocation with ramp constraints, add explicit secant slope constraints.
+ Default True. If False, only polynomial derivative constraints are used, which may
+ allow numerical derivatives between mesh points to exceed the ramp limit (by ~15%).
+ """
+ if task == "Tsh":
+ Pchamber = {"setpt": [0.1], "dt_setpt": [180.0], "ramp_rate": 0.5}
+ Tshelf = {"min": -45.0, "max": 120.0}
+ if tsh_ramp_rate is not None:
+ Tshelf["max_ramp_rate"] = tsh_ramp_rate
+ runner = pyomo_opt.optimize_Tsh_pyomo
+ args = (vial, product, ht, Pchamber, Tshelf, dt, eq_cap, nVial)
+ elif task == "Pch":
+ Pchamber = {"min": 0.05}
+ if pch_ramp_rate is not None:
+ Pchamber["max_ramp_rate"] = pch_ramp_rate
+ # Shelf multi-step schedule with sufficient time for drying completion
+ # NOTE: dt_setpt in MINUTES (Pyomo internally uses same convention as scipy)
+ # High resistance products (A1=20) need ~86 hours, use 100 hours for margin
+ Tshelf = {
+ "init": -35.0,
+ "setpt": [-20.0, 120.0],
+ "dt_setpt": [300.0, 5700.0],
+ "ramp_rate": 1.0,
+ }
+ runner = pyomo_opt.optimize_Pch_pyomo
+ args = (vial, product, ht, Pchamber, Tshelf, dt, eq_cap, nVial)
+ elif task == "both":
+ Pchamber = {"min": 0.05, "max": 0.5}
+ if pch_ramp_rate is not None:
+ Pchamber["max_ramp_rate"] = pch_ramp_rate
+ Tshelf = {"min": -45.0, "max": 120.0, "init": -35.0}
+ if tsh_ramp_rate is not None:
+ Tshelf["max_ramp_rate"] = tsh_ramp_rate
+ runner = pyomo_opt.optimize_Pch_Tsh_pyomo
+ args = (vial, product, ht, Pchamber, Tshelf, dt, eq_cap, nVial)
+ else:
+ raise ValueError(f"Unknown task '{task}'")
+
+ use_fd = method.lower() == "fd"
+ treat_eff = (not use_fd) and bool(effective_nfe)
+
+ start = time.perf_counter()
+ meta: dict[str, Any] = {}
+ try:
+ res = runner(
+ *args,
+ n_elements=int(n_elements),
+ n_collocation=int(n_collocation),
+ use_finite_differences=use_fd,
+ treat_n_elements_as_effective=treat_eff,
+ warmstart_scipy=warmstart,
+ return_metadata=True,
+ tee=False,
+ use_secant_ramp_constraints=use_secant_ramp_constraints,
+ )
+ if isinstance(res, dict) and "output" in res:
+ traj = res["output"]
+ meta = res.get("metadata", {})
+ else:
+ traj = res
+ success = True
+ message = "pyomo run completed"
+ except Exception as e:
+ traj = np.empty((0, 7))
+ success = False
+ message = f"pyomo failure: {e.__class__.__name__}: {e}"[:300]
+ wall = time.perf_counter() - start
+
+ if success and isinstance(traj, np.ndarray) and traj.size:
+ default_t = float(traj[-1, 0])
+ else:
+ default_t = None
+ objective_time_hr = (
+ float(meta.get("objective_time_hr", default_t)) if success else None
+ )
+
+ solver_info = {
+ "status": meta.get("status"),
+ "termination_condition": meta.get("termination_condition"),
+ "ipopt_iterations": meta.get("ipopt_iterations"),
+ "n_points": meta.get("n_points"),
+ "staged_solve_success": meta.get("staged_solve_success"),
+ }
+
+ # Discretization reporting
+ if use_fd:
+ total_mesh_points = int(n_elements) + 1 # simple FD time mesh
+ else:
+ # For collocation, total interior points per element is n_collocation
+ # Effective parity: treat effective_nfe True => report n_elements as comparable to FD
+ total_mesh_points = int(n_elements) * int(n_collocation) + 1
+ discretization = {
+ "method": "fd" if use_fd else "colloc",
+ "n_elements_requested": int(n_elements),
+ "n_elements_applied": int(n_elements), # could differ if transformation adjusts
+ "n_collocation": int(n_collocation) if not use_fd else None,
+ "effective_nfe": bool(treat_eff) if not use_fd else False,
+ "total_mesh_points": total_mesh_points,
+ }
+ # Merge any mesh_info from metadata if present
+ mesh_info = meta.get("mesh_info")
+ if isinstance(mesh_info, dict):
+ discretization.update(
+ {k: v for k, v in mesh_info.items() if k not in discretization}
+ )
+
+ return {
+ "trajectory": traj,
+ "success": success,
+ "message": message,
+ "wall_time_s": wall,
+ "objective_time_hr": objective_time_hr,
+ "solver": solver_info,
+ "solver_stats": {},
+ "raw": traj,
+ "warmstart_used": warmstart,
+ "discretization": discretization,
+ }
diff --git a/benchmarks/grid_analysis.ipynb b/benchmarks/grid_analysis.ipynb
new file mode 100644
index 0000000..50fbc62
--- /dev/null
+++ b/benchmarks/grid_analysis.ipynb
@@ -0,0 +1,1918 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "7d8b0d88",
+ "metadata": {},
+ "source": [
+ "# Benchmark Analysis: Pyomo vs Scipy\n",
+ "\n",
+ "This notebook analyzes pre-generated benchmark data comparing Pyomo optimization methods (Finite Differences and Orthogonal Collocation) against the Scipy baseline.\n",
+ "\n",
+ "## Data Files\n",
+ "\n",
+ "**Full benchmark files** (tracked in git):\n",
+ "- `baseline_Tsh_3x3.jsonl` - Shelf temperature optimization, 3x3 parameter grid\n",
+ "- `baseline_Pch_3x3.jsonl` - Chamber pressure optimization, 3x3 parameter grid\n",
+ "\n",
+ "**Additional benchmark files** (gitignored, regenerate locally):\n",
+ "- `both_no_ramp_3x3.jsonl` - Both Tsh and Pch optimization\n",
+ "- `both_ramp_constrained.jsonl` - Ramp-constrained optimization (Tsh ≤10°C/hr, Pch ≤0.1 Torr/hr)\n",
+ "- `high_resistance_*.jsonl` - Stress tests with high product resistance\n",
+ "- `sanity_check.jsonl` - Quick 2x2 validation grid\n",
+ "\n",
+ "## Data Generation (CLI Reference)\n",
+ "\n",
+ "```bash\n",
+ "# Shelf temperature optimization (3x3 grid)\n",
+ "python benchmarks/grid_cli.py generate --task Tsh --scenario baseline \\\n",
+ " --vary product.A1=16,18,20 --vary ht.KC=2.75e-4,3.5e-4,4e-4 \\\n",
+ " --methods scipy,fd,colloc --n-elements 12 \\\n",
+ " --out benchmarks/results/baseline_Tsh_3x3.jsonl --force\n",
+ "\n",
+ "# Chamber pressure optimization (3x3 grid)\n",
+ "python benchmarks/grid_cli.py generate --task Pch --scenario baseline \\\n",
+ " --vary product.A1=16,18,20 --vary ht.KC=2.75e-4,3.5e-4,4e-4 \\\n",
+ " --methods scipy,fd,colloc --n-elements 12 \\\n",
+ " --out benchmarks/results/baseline_Pch_3x3.jsonl --force\n",
+ "\n",
+ "# Both Tsh and Pch optimization (3x3 grid)\n",
+ "python benchmarks/grid_cli.py generate --task both --scenario baseline \\\n",
+ " --vary product.A1=16,18,20 --vary ht.KC=2.75e-4,3.5e-4,4e-4 \\\n",
+ " --methods scipy,fd,colloc --n-elements 12 \\\n",
+ " --out benchmarks/results/both_no_ramp_3x3.jsonl --force\n",
+ "\n",
+ "# With trajectories for visualization (adds --save-trajectories flag)\n",
+ "python benchmarks/grid_cli.py generate --task both --scenario baseline \\\n",
+ " --vary product.A1=18 --vary ht.KC=3.3e-4 \\\n",
+ " --methods scipy,fd,colloc --n-elements 12 --save-trajectories \\\n",
+ " --out benchmarks/results/both_trajectories.jsonl --force\n",
+ "\n",
+ "# Ramp-constrained verification (exercising Pyomo derivative constraints)\n",
+ "python benchmarks/grid_cli.py generate --task both --scenario baseline \\\n",
+ " --vary product.A1=18 --vary ht.KC=3.3e-4 \\\n",
+ " --methods scipy,fd,colloc --save-trajectories \\\n",
+ " --tsh-ramp 10.0 --pch-ramp 0.1 \\\n",
+ " --out benchmarks/results/both_ramp_constrained.jsonl --force\n",
+ "```\n",
+ "\n",
+ "## Key Metrics\n",
+ "\n",
+ "- **Objective parity**: % difference in optimized drying time (hr)\n",
+ "- **Speedup**: wall-clock time ratio (scipy/pyomo)\n",
+ "- **Success rate**: solver convergence across parameter grid"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4daae57e",
+ "metadata": {},
+ "source": [
+ "# Grid Analysis: FD vs Collocation vs Scipy\n",
+ "\n",
+ "Load and analyze benchmark results across a parameter grid comparing scipy baseline with Pyomo finite differences (FD) and orthogonal collocation methods.\n",
+ "\n",
+ "**Note on First-Run Overhead**: The first Pyomo run may be slower (~1s for FD) due to one-time initialization (Pyomo DAE transformation, IPOPT library loading). Subsequent runs are faster (~0.03-0.05s). This is expected behavior and does not affect objective parity or speedup analysis.\n",
+ "\n",
+ "**Parameter Grid**: The benchmarks use a 3×3 grid varying `product.A1` (product resistance coefficient) and `ht.KC` (heat transfer coefficient) to test optimizer robustness across different product and vial characteristics."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 1,
+ "id": "4d3637ff",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2025-11-15T01:58:21.819127Z",
+ "iopub.status.busy": "2025-11-15T01:58:21.818909Z",
+ "iopub.status.idle": "2025-11-15T01:58:22.015994Z",
+ "shell.execute_reply": "2025-11-15T01:58:22.014659Z"
+ }
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "Loaded 3 records from /home/bernalde/repos/LyoPRONTO/benchmarks/results/both_ramp_constrained.jsonl\n",
+ "Task type: both\n",
+ "Output folder: /home/bernalde/repos/LyoPRONTO/benchmarks/results/both_ramp_constrained\n",
+ "Scipy: 1, FD: 1, Collocation: 1\n",
+ "Trajectories stored: True\n",
+ "Ramp constraints: Tsh=10.0 C/hr, Pch=0.1 Torr/hr\n",
+ "\n",
+ "Sample parameters: A1=18.0, KC=3.30e-04\n"
+ ]
+ }
+ ],
+ "source": [
+ "# Load benchmark grid data\n",
+ "import json\n",
+ "import sys\n",
+ "from pathlib import Path\n",
+ "\n",
+ "import numpy as np\n",
+ "import pandas as pd\n",
+ "\n",
+ "# Find repository root\n",
+ "repo_root = Path.cwd()\n",
+ "while not (repo_root / \".git\").exists() and repo_root != repo_root.parent:\n",
+ " repo_root = repo_root.parent\n",
+ "\n",
+ "# =============================================================================\n",
+ "# PATH CONFIGURATION - All paths defined here, reused throughout notebook\n",
+ "# =============================================================================\n",
+ "RESULTS_DIR = repo_root / \"benchmarks\" / \"results\"\n",
+ "\n",
+ "# Input data file for grid analysis - choose one:\n",
+ "# Baseline benchmarks (tracked in git):\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'baseline_Tsh_3x3.jsonl'\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'baseline_Pch_3x3.jsonl'\n",
+ "\n",
+ "# Additional benchmarks (gitignored, regenerate locally):\n",
+ "GRID_DATA_PATH = RESULTS_DIR / \"both_ramp_constrained.jsonl\"\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'both_no_ramp_3x3.jsonl'\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'high_resistance_Tsh_3x3.jsonl'\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'high_resistance_Pch_3x3.jsonl'\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'high_resistance_both_3x3.jsonl'\n",
+ "# GRID_DATA_PATH = RESULTS_DIR / 'sanity_check.jsonl'\n",
+ "\n",
+ "# Trajectory data file (can be same as GRID_DATA_PATH if it has trajectories,\n",
+ "# or a separate single-point file for faster trajectory visualization)\n",
+ "TRAJ_DATA_PATH = GRID_DATA_PATH # Use same file by default\n",
+ "# TRAJ_DATA_PATH = RESULTS_DIR / 'baseline_Tsh_3x3_trajectories.jsonl' # Single-point alternative\n",
+ "\n",
+ "# Output folder: derived from input filename (e.g., baseline_Tsh_3x3.jsonl -> baseline_Tsh_3x3/)\n",
+ "EXPERIMENT_NAME = GRID_DATA_PATH.stem.replace(\n",
+ " \"_summary\", \"\"\n",
+ ") # Strip _summary suffix if present\n",
+ "OUTPUT_DIR = RESULTS_DIR / EXPERIMENT_NAME\n",
+ "OUTPUT_DIR.mkdir(exist_ok=True)\n",
+ "\n",
+ "# Output files (plots and exports) - saved to experiment-specific folder\n",
+ "GRID_OBJECTIVE_DIFF_PNG = OUTPUT_DIR / \"grid_objective_diff.png\"\n",
+ "GRID_SPEEDUP_PNG = OUTPUT_DIR / \"grid_speedup.png\"\n",
+ "GRID_COMPARISON_CSV = OUTPUT_DIR / \"grid_comparison.csv\"\n",
+ "TRAJ_COMPARISON_PNG = OUTPUT_DIR / \"trajectory_comparison.png\"\n",
+ "FLUX_ANALYSIS_PNG = OUTPUT_DIR / \"flux_analysis.png\"\n",
+ "CONSTRAINT_ANALYSIS_PNG = OUTPUT_DIR / \"constraint_analysis.png\"\n",
+ "RAMP_ANALYSIS_PNG = OUTPUT_DIR / \"ramp_rate_analysis.png\"\n",
+ "# =============================================================================\n",
+ "\n",
+ "if not GRID_DATA_PATH.exists():\n",
+ " raise FileNotFoundError(\n",
+ " f\"Data file not found: {GRID_DATA_PATH}\\n\\n\"\n",
+ " f\"Generate benchmark data using grid_cli.py first.\\n\"\n",
+ " f\"See the CLI reference in the notebook header for example commands.\"\n",
+ " )\n",
+ "\n",
+ "grid_recs = []\n",
+ "with GRID_DATA_PATH.open(\"r\") as f:\n",
+ " for line in f:\n",
+ " line = line.strip()\n",
+ " if line:\n",
+ " grid_recs.append(json.loads(line))\n",
+ "\n",
+ "# Detect task type from data\n",
+ "task_type = grid_recs[0].get(\"task\", \"unknown\") if grid_recs else \"unknown\"\n",
+ "print(f\"Loaded {len(grid_recs)} records from {GRID_DATA_PATH}\")\n",
+ "print(f\"Task type: {task_type}\")\n",
+ "print(f\"Output folder: {OUTPUT_DIR}\")\n",
+ "\n",
+ "# Separate by method\n",
+ "scipy_recs = [r for r in grid_recs if r.get(\"pyomo\") is None]\n",
+ "fd_recs = [\n",
+ " r\n",
+ " for r in grid_recs\n",
+ " if r.get(\"pyomo\") and r[\"pyomo\"].get(\"discretization\", {}).get(\"method\") == \"fd\"\n",
+ "]\n",
+ "colloc_recs = [\n",
+ " r\n",
+ " for r in grid_recs\n",
+ " if r.get(\"pyomo\") and r[\"pyomo\"].get(\"discretization\", {}).get(\"method\") == \"colloc\"\n",
+ "]\n",
+ "\n",
+ "print(f\"Scipy: {len(scipy_recs)}, FD: {len(fd_recs)}, Collocation: {len(colloc_recs)}\")\n",
+ "\n",
+ "# Check if trajectories are available\n",
+ "has_trajectories = any(\n",
+ " r.get(\"scipy\", {}).get(\"trajectory\") is not None for r in grid_recs\n",
+ ")\n",
+ "print(f\"Trajectories stored: {has_trajectories}\")\n",
+ "if not has_trajectories:\n",
+ " print(\n",
+ " \" -> No trajectories. For trajectory plots, regenerate with --save-trajectories flag\"\n",
+ " )\n",
+ "\n",
+ "# Check for ramp constraints\n",
+ "ramp_constraints = None\n",
+ "for r in grid_recs:\n",
+ " if r.get(\"pyomo\") and r[\"pyomo\"].get(\"ramp_constraints\"):\n",
+ " ramp_constraints = r[\"pyomo\"][\"ramp_constraints\"]\n",
+ " break\n",
+ "if ramp_constraints:\n",
+ " print(\n",
+ " f\"Ramp constraints: Tsh={ramp_constraints.get('tsh_ramp_rate')} C/hr, Pch={ramp_constraints.get('pch_ramp_rate')} Torr/hr\"\n",
+ " )\n",
+ "else:\n",
+ " print(\"Ramp constraints: None\")\n",
+ "\n",
+ "# Show sample record structure\n",
+ "if scipy_recs:\n",
+ " print(\n",
+ " f\"\\nSample parameters: A1={scipy_recs[0]['grid']['param1']['value']}, KC={scipy_recs[0]['grid']['param2']['value']:.2e}\"\n",
+ " )"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "id": "32752773",
+ "metadata": {
+ "execution": {
+ "iopub.execute_input": "2025-11-15T01:58:22.018022Z",
+ "iopub.status.busy": "2025-11-15T01:58:22.017719Z",
+ "iopub.status.idle": "2025-11-15T01:58:22.033999Z",
+ "shell.execute_reply": "2025-11-15T01:58:22.032886Z"
+ }
+ },
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "\n",
+ "\n",
+ "
\n",
+ " \n",
+ " \n",
+ " | \n",
+ " A1 | \n",
+ " KC | \n",
+ " scipy_obj | \n",
+ " fd_obj | \n",
+ " colloc_obj | \n",
+ " fd_diff_pct | \n",
+ " colloc_diff_pct | \n",
+ " scipy_wall | \n",
+ " fd_wall | \n",
+ " colloc_wall | \n",
+ " fd_speedup | \n",
+ " colloc_speedup | \n",
+ "
\n",
+ " \n",
+ " \n",
+ " \n",
+ " | 0 | \n",
+ " 18.0 | \n",
+ " 3.30e-04 | \n",
+ " 11.756321 | \n",
+ " 11.266449 | \n",
+ " 10.913175 | \n",
+ " -4.166888 | \n",
+ " -7.171859 | \n",
+ " 13.527728 | \n",
+ " 0.720486 | \n",
+ " 0.060607 | \n",
+ " 18.775847 | \n",
+ " 223.205661 | \n",
+ "
\n",
+ " \n",
+ "
\n",
+ "