Johns Hopkins University — Semester 2 Graduate Research Advisors: Dr. Ján Drgoňa · Dr. Liang Wu
This project develops and benchmarks surrogate-assisted Differentiable Predictive Control (DPC) for a 2D phase-change material (PCM) thermal energy storage system governed by the Stefan enthalpy PDE.
The core challenge: the full PDE simulator is too slow for real-time optimization. We train neural operator surrogates (DeepONet, FNO, CNO) to approximate the PDE, then differentiate through them to solve the control problem via gradient descent — following the DPC framework of Drgoňa et al. (2022). The final deployable artifact is a behaviour-cloned offline policy (single neural network forward pass, ~2 ms/decision) that matches the quality of online gradient-based solvers at zero optimization cost.
Control objective: Regulate State of Charge (SoC = 0.50) in a 0.1 m × 0.1 m PCM slab under a sinusoidal inlet-temperature disturbance, by modulating HTF mass flow rate ṁ ∈ [0.001, 0.05] kg/s.
| Stage | Status |
|---|---|
| PDE simulator (NumPy FDM) | Complete — FEM-validated (4-level cascade, SoC RMSE 0.0008 vs FEniCS) |
| Dataset generation | Complete — 1000 trajectories, trajectory-level train/test split |
| Surrogate training | Complete — DeepONet2D (MSE 1.47e-4), FNO2D, CNO2D |
| Controller benchmark (6-way) | Complete — see Key Results |
| Offline BC policy | Complete — 200 epochs, val MSE 2.22e-2, ~2 ms inference |
| Multi-step surrogate rollout training | Planned (addresses compound prediction error in Surrogate-DPC) |
| Paper draft | In progress (docs/ACC_PCM_Surrogate_DPC.pdf) |
| Controller | SoC RMSE | Pk-Pk SoC | Mean ṁ (g/s) | Replan (s) | Wall time |
|---|---|---|---|---|---|
| Open-loop (ṁ = const) | 0.0583 | 0.0168 | 25.00 | — | <1 s |
| DPC-DeepONet2D | 0.0519 | 0.0050 | 12.92 | 2.5 | 164 s |
| DPC-FNO2D | 0.0502 | 0.0051 | 26.46 | 8.7 | 518 s |
| MPPI-FNO (K=200) | 0.0516 | 0.0041 | 16.96 | 10.8 | 658 s |
| BC-Offline (neural policy) | 0.0549 | 0.0100 | 25.26 | 0.002 | 0.2 s |
All feedback controllers outperform open-loop. Best SoC tracking: DPC-FNO2D (−13.8%). Best flow smoothness: MPPI-FNO. BC-Offline matches online DPC quality within 9% at 0.002 s/decision vs 2.5–10.8 s for online solvers — making it the only controller deployable at the re-plan frequency used here.
- Controllers (DPC-DeepONet2D, DPC-FNO2D, MPPI-FNO): Online optimisation oracles — upper-bound benchmarks for each surrogate architecture, not intended for real-time deployment as-is.
- BC-Offline: Deployable policy — trained once via behaviour cloning on DPC decisions, then deployed as a single forward pass. Matches Drgoňa et al. (2022) offline-policy deployment mode.
| R2 | SoC RMSE | TV(ṁ) change |
|---|---|---|
| 0.0 | 0.05384 | baseline |
| 0.5 | 0.05386 | −15.9% |
| 2.0 | 0.05385 | −20.7% |
R2 reduces flow-rate chatter by ~16–21% with negligible SoC cost (flat to 4 decimal places).
Stefan phase-change problem using the apparent heat capacity (enthalpy) formulation (Voller, 1990):
ρ · Ceff(T) · ∂T/∂t = k · ∇²T
Ceff(T) = Cp + (ΔH / (σ√(2π))) · exp(-(T - Tm)² / (2σ²))
| Parameter | Value | Units |
|---|---|---|
| Melting temperature Tm | 310.15 (37°C) | K |
| Latent heat ΔH | 210,000 | J/kg |
| Thermal conductivity k | 0.22 | W/(m·K) |
| Density ρ | 790 | kg/m³ |
| Specific heat Cp | 2,000 | J/(kg·K) |
| Mushy-zone half-width σ | 0.5 | K |
| Domain Lx = Ly | 0.1 | m |
| Grid Nx × Ny | 15 × 15 | cells |
| Time step dt | 30 | s |
- Left wall (x=0): Robin convective —
-k ∂T/∂x = h_c(ṁ)(Tin − T),h_c = h₀(ṁ/ṁ_ref) - All other walls: Insulated Neumann —
∂T/∂n = 0
h(T) = Cp·(T − Tm) + (ΔH/2)·(1 + erf((T − Tm) / (σ√2)))
E = ρ · ∫∫ h(T) dx dy
SoC = clip((E − E_min) / E_range, 0, 1)
E_min and E_max anchored at Tm ± 15 K to avoid SoC saturation near phase-transition endpoints.
Following Drgoňa et al. (2022):
min Σ_{k=0}^{H-1} Q·(SoC_k − SoC_r)² + R·ṁ_k² + R2·(ṁ_k − ṁ_{k−1})²
Parameters: Q=100, R=0.1, R2=0.5, H=10, resolve_every=4 timesteps.
The NumPy FDM simulator (pcm_datagen/simulator_2d.py) was validated through a four-level cascade
before use as a data generator:
| Level | Test | Result |
|---|---|---|
| 1 | Analytical BC energy balance | PASS (< 1e-10 K error) |
| 2a | Pure-liquid cross-simulator (NumPy vs PhiFlow) | PASS (0.001 K) |
| 2b | Mushy-zone cross-simulator | PASS (0.133 K — expected: different mushy kernels) |
| 2c | Left-wall BC formula | PASS (0.02% error) |
| 3 | SoC enthalpy self-consistency | PASS (< 0.001) |
| 4 | FEM ground-truth (FEniCS 40×40, Crank-Nicolson) | PASS — SoC RMSE = 0.0008 |
Full reports: docs/validation/pcm_validation_report.pdf, docs/validation/numpy_datagen_justification.pdf
| ID | Status | Issue | Resolution |
|---|---|---|---|
| D1 | Fixed | SoC normalisation at Tm±5 K caused saturation in ~40% of trajectories | Widened to ±15 K in simulator_2d.py:141 |
| D2 | Documented | Gaussian σ=0.5 K narrower than grid spacing; Voller (1990) warns of latent-heat jump-over risk | Explicit Euler CFL check present; noted in simulator |
| D3 | Documented | h_conv linear in ṁ is a simplification (Dittus-Boelter exponent: 0.8) | Modelling assumption noted in simulator_2d.py |
| D5 | Complete | DPC should deploy offline neural policy, not online gradient descent | BC policy trained (200 epochs, val MSE 2.22e-2); deployed as BC-Offline |
| D6 | Fixed | Timestep-level train/test split caused data leakage | Changed to trajectory-level split in dataset_2d.py:311 |
| D8 | Fixed | PCM-37 parameters uncited | Rubitherm RT37HC datasheet citation added to simulator_2d.py:45 |
Research-Work/
│
├── src/ # Controllers and training scripts
│ ├── neural_operator_2d.py # DeepONet2D (branch/trunk, augmented trunk)
│ ├── neural_operator_benchmarks.py # FNOSurrogate2D, CNOSurrogate2D
│ ├── dpc_controller_2d.py # Online DPC via DeepONet (oracle baseline)
│ ├── dpc_controller_2d_fno.py # Online DPC via FNO (oracle baseline)
│ ├── dpc_controller_2d_mppi.py # MPPI controller (H=20, K=200)
│ ├── dpc_controller_2d_phiflow.py # DPC via differentiable PhiFlow PDE (no surrogate error)
│ ├── dpc_controller_2d_offline.py # BC-Offline policy inference (single forward pass)
│ ├── train_deeponet_2d.py # DeepONet2D training
│ ├── train_bc_policy.py # Behaviour-cloning policy training
│ ├── collect_bc_data.py # Collect (state → m_opt) pairs from online DPC
│ ├── run_comparison_5way_2d.py # 6-controller comparison with timing breakdown
│ ├── run_r2_sweep.py # R2 smoothness sweep
│ ├── run_horizon_sweep.py # Prediction horizon H sweep
│ ├── run_diurnal_comparison.py # 24-hour diurnal disturbance comparison
│ └── train_surrogates_benchmark.py # DeepONet vs FNO vs CNO benchmark
│
├── scripts/ # Post-processing and utility tools
│ ├── analyze_comparison.py # Post-run diagnostics and result parsing
│ ├── make_animation_2d.py # 2D trajectory animation
│ ├── make_animation_24h.py # 24-hour scenario animation
│ └── consolidate_to_hdf5.py # Convert NPZ trajectories → HDF5
│
├── pcm_datagen/ # Data generation engine (FEM-validated)
│ ├── simulator_2d.py # Stefan PDE solver (FDM, NumPy, cell-centred)
│ ├── simulator_2d_phiflow.py # Differentiable PhiFlow step (for PhiFlow-DPC)
│ ├── dataset_2d.py # Trajectory generator; trajectory-level split
│ └── signals.py # Randomised ṁ and Tin signal generators
│
├── fem_benchmark/
│ └── fenics_pcm_2d.py # FEniCS benchmark (40×40, Crank-Nicolson) — validation oracle
│
├── docs/ # All documentation
│ ├── ACC_PCM_Surrogate_DPC.tex # Conference paper (ACC submission)
│ ├── ACC_PCM_Surrogate_DPC.pdf
│ ├── IEEEtran.cls
│ ├── PCM_DPC_Report.tex # Semester research report
│ ├── PCM_DPC_Report.pdf
│ ├── PCM_DPC_Report.docx
│ ├── PCM_DPC_Author_Response.md # Reviewer response
│ ├── REVIEW_LOG.md # Living peer-review log
│ ├── CONSISTENCY_LOG.md # Cross-file consistency audit
│ ├── NEXT_STEPS.md
│ └── validation/ # Validation reports and figures
│ ├── pcm_validation_report.pdf
│ ├── numpy_datagen_justification.pdf
│ └── level{1,2,3,4}_*.png
│
├── models/ # Trained weights — git-ignored; available on request
├── data/ # Datasets — git-ignored
├── results/ # Experiment outputs and plots — git-ignored
├── output/ # Script output directory — git-ignored
│
├── generate_dataset.py # Dataset generation entry point
├── run_validation.py # Run full 4-level validation cascade
└── requirements.txt
Note:
models/andoutput/dataset_2d/are git-ignored (large binaries). Pre-trained weights and datasets available on request.
pip install torch numpy scipy matplotlib phiflow neuraloperatorGPU (CUDA) recommended. Tested on RTX 4070 Laptop (CUDA 12).
python run_validation.py
# Runs Levels 1–3; Level 4 requires FEniCS Docker (see below)python -c "
from pcm_datagen.dataset_2d import generate_dataset_2d
generate_dataset_2d(N=1000, output_dir='output/dataset_2d', n_workers=6, consolidate=True)
"python -u src/train_deeponet_2d.py \
--data_dir output/dataset_2d \
--epochs 100 --batch 512 --p 128 --subsample 4 \
--save models/deeponet_2d.ptpython src/run_comparison_5way_2d.py --skip_phiflow
# Output: output/plots/comparison_5way_2d.png
# output/comparison_5way_summary.txt# Collect DPC decisions (~46 min, 150 episodes)
python -u src/collect_bc_data.py --n_episodes 150 --H 5 --n_iter 20 --seed 42
# Train via behaviour cloning (~25 sec, 200 epochs)
python -u src/train_bc_policy.py --epochs 200
# Output: models/dpc_policy_2d.pt# R2 smoothness penalty sweep
python src/run_r2_sweep.py --r2_values 0 0.01 0.05 0.1 0.2 0.5 1.0 2.0
# Prediction horizon sweep
python src/run_horizon_sweep.py --device cuda --H_values 5 10 20 40
# 24-hour diurnal scenario
python src/run_diurnal_comparison.py --device cudaLevel 4 validation requires FEniCS via Docker:
docker cp fem_benchmark/fenics_pcm_2d.py <container_id>:/tmp/
docker exec <container_id> python3 /tmp/fenics_pcm_2d.py
docker cp <container_id>:/tmp/fenics_pcm_2d_output.json validation_results/Container: quay.io/fenicsproject/stable — FEniCS 2019.2.0.dev0, Python 3.10.12.
| Component | Tool |
|---|---|
| Deep learning | PyTorch (CUDA) |
| Neural operators | neuraloperator (FNO/CNO), custom DeepONet |
| Differentiable PDE | PhiFlow ≥ 3.4 |
| FEM oracle | FEniCS 2019 (Docker) |
| Numerics | NumPy, SciPy |
- Voller, V.R. (1990). Fixed grid techniques for phase change problems. Int. J. Num. Methods Eng.
- Drgoňa, J. et al. (2022). Differentiable predictive control. Neural Networks.
- Sarkar, S., Drgoňa, J., Goswami, D. (2025). Learning to Control PDEs with DPC and TI-Neural Operators. (preprint)
- Nayak, P., Goswami, D. (2025). TI-DeepONet: Learnable Time Integration for Stable Neural Operators.
- Lu, L. et al. (2021). Learning nonlinear operators via DeepONet. Nature Machine Intelligence.
- Rubitherm GmbH. RT37HC datasheet. rubitherm.eu.