[sbryngelson]

What

Three changes to the mean-source machinery of #2177 and #2182 (m_body_forces.fpp).

1. mean_temp_wall, a self-similar thermostat for thermally developed flow past isothermal walls. With mean_temp_wall = T_s, the mean gas temperature source becomes the local source of the periodic part psi = T - T_s of a thermally fully developed field T = T_s + psi exp(-beta x) (Sun, Tenneti & Subramaniam, IJHMT 86, 2015): S = beta [rho cp u psi - d(lam psi)/dx - lam dpsi/dx + beta lam psi]. The mean_temp controller still sets the source total qdot V_f, and in this mode holds the flux-weighted mixing-cup temperature T_m instead of the mass mean T_b. beta is the root of qdot V_f = beta (A1 + beta A2), with A1 the sum of the same discrete terms the source applies (so the balance holds next to bodies and for varying properties) and A2 = sum(lam psi dV); beta, A1 and A2 are columns 7 to 9 of mean_temp.dat. The uniform source is unchanged when the option is off. The thermostat now accepts no-slip and slip walls in y and z (x stays periodic). The validator requires, with mean_temp_wall, chem_params%diffusion, mean_mflux(1) > 0 (the mixing cup needs a held flow) and every wall and immersed boundary isothermal at mean_temp_wall.
2. Immersed-boundary support for mean_dpdx/mean_mflux and mean_qdot/mean_temp, which the validator rejected: the controller sums (fluid volume, mass flux, work, temperature averages) run over fluid cells only (zero volume inside a body), and the force and the heat source are not applied inside bodies. The two validator rejections are removed.
3. examples/2D_thermal_channel: a plane channel between isothermal walls, the Nu = 7.541 acceptance case (the Nu is computed from mean_temp.dat).

Tests: IBM -> Periodic Circle -> Mean Pressure Gradient (single and two fluid) and IBM Reacting Surface -> Thermostat and -> Wall Temperature. The channel has no golden test: see the end of Verification.

Verification

Results change only with mean_temp_wall set or with ib and a mean source together; the existing tests of #2177 and #2182 (Bodyforces, MeanDpdx, Thermostat) are unchanged.

Nu = 7.541 (plane channel, constant properties), Tioga. Hot wall (330 K) and cold wall (270 K) pair, whose viscous heating cancels in the average: Nu = 7.7725 (hot) and 7.3224 (cold) at 64 cells across, average 7.548, +0.09% from 7.541; at 32 cells 7.7885 and 7.3238, average 7.556 (+0.20%). Nu is -qdot H^2 / (k (T_w - T_0)) from mean_temp.dat. These ran before #2154 merged; on the final head the 32-cell hot wall gives 7.7821 (7.7885 before the merge).

Immersed boundary controller (IBM -> Periodic Circle -> Mean Pressure Gradient run for 800 steps): the logged flux settles at 0.04999999 against the target 0.05. An independent mean of rho u over the fluid cells only, from the saved fields, is 0.0499 (0.24% low, from my geometric fluid mask), while the mean over all cells is 0.0484, 3% low: what a volume that included the body would give.

GPU tests, at the final head 014f4fb95 (master with #2193 to #2196 merged): AMD MI250X (Tioga, amdflang, OpenMP offload), goldens generated there: the 4 new tests pass; Bodyforces 6, MeanDpdx 1 and Thermostat pass. NVIDIA + NVHPC (A100, OpenACC): the 4 new tests pass against the Tioga goldens at the default tolerance for the immersed-boundary mean-flux tests and at 1e-6 for the two reacting-surface thermostat tests; Bodyforces, MeanDpdx and Thermostat (4 tests) pass. ./mfc.sh precheck and the Python validator tests (179) pass.

Merging #2194 exposed one conflict, fixed here: its rule that the temperature controller needs a constant time step rejected the two reacting-surface thermostat tests (the example runs an adaptive step), so they now use cfl_const_dt and their goldens were regenerated.

An independent review led to: the applied source and the beta balance now share one set of discrete terms (they differed next to bodies); the validator requirements listed above; a one-time message when the beta quadratic has no real root (it then uses the vertex 2 src/A1); and one macro for the repeated "zero inside bodies" pattern.

Why the channel has no golden test. A golden of the channel differed between AMD and NVIDIA by up to 7e-3 after five steps even at a controller gain of 0.001, so it would not hold on other CI platforms. The cause is not isolated: the thermochemistry temperature solve stops at a relative tolerance of 1e-6, and both the isothermal wall and the controller (which differences the temperature and divides by dt) go through it. The channel is therefore checked by the Nu runs above, and the example is skipped in the Example suite for the same reason. The four immersed-boundary tests, which use a prescribed-temperature body and a small gain (0.001), agree between vendors at the default tolerance.

Not shown: only the plane channel is verified (the docs say so). The hot and cold wall cases in flows with streamlines of very different speed relax at different rates and beta is not converged there. The IBM Reacting Surface -> Thermostat -> Wall Temperature test exercises that mode as a regression golden only (its temperature gap is far above the few percent the docs advise).

Name	CPU	GPU	Compiler / offload	Test flags
wingtip (Georgia Tech)	Intel Xeon Gold 6338 (Ice Lake)	NVIDIA A100 80GB PCIe (cc 8.0)	NVHPC nvfortran 25.11, OpenACC	--gpu acc -g 0 1 2 3 -j 4
Tioga (LLNL)	AMD EPYC 7A53 (64 cores)	AMD Instinct MI250X (gfx90a)	AMD flang (ROCm AFAR 24.3.0), OpenMP offload	--gpu mp --mpi -j 4

This PR was made with Claude Code (AI).

[github-actions]

Claude Code Review

Head SHA: 21fbc82

Files changed:

• 11
• docs/documentation/case.md
• docs/documentation/equations.md
• examples/2D_thermal_channel/case.py
• src/simulation/m_body_forces.fpp
• src/simulation/m_global_parameters.fpp
• src/simulation/m_time_steppers.fpp
• toolchain/mfc/case_validator.py
• toolchain/mfc/params/definitions.py
• toolchain/mfc/params/descriptions.py
• toolchain/mfc/test/cases.py

Findings:

• medium src/simulation/m_body_forces.fpp: In s_compute_mts_terms, Ys is built only from cell j, but the loop do o = 1, 3 calls get_mixture_thermal_conductivity_mixavg(T_o, Ys, lam_o) with the temperature of cell j+o-2 and cell j's composition. The neighbour conductivities lam(1) and lam(3) therefore pair the wrong Y with T. This is consistent between the applied source and the controller's A1/A2, because both use mts_g1 and mts_g2. It is still not the conductivity of the neighbour cell, so the discrete d(lam psi)/dx is first-order inconsistent where composition varies in x. It is invisible in the constant-composition channel that the PR verifies. Either compute Ys from the neighbour cell's conserved variables, or document the approximation. I am moderately sure of this; whether it was intended is unclear.
• medium src/simulation/m_body_forces.fpp: I could not confirm that get_mixture_thermal_conductivity_mixavg carries GPU_ROUTINE(seq) and is safe to call in a device loop. It is not in the prepared context, and the call is new to this module. get_mixture_specific_heat_cp_mass is already used in a GPU loop here. Please verify the routine, since GPU-only failures are not caught by the CPU tests. 2D_thermal_channel is also excluded from the golden-file tests because it differed between AMD and NVIDIA GPUs.

[sbryngelson]

What this fixes / adds. mean_temp_wall: a mean gas temperature source for thermally developed flow (the periodic part psi = T - T_s of T_s + psi exp(-beta x)), with immersed-boundary support (isothermal IB walls at the same temperature), and the periodic-circle mean-pressure-gradient tests.

How it was tested

• wingtip, --gpu acc (NVHPC 25.11, 4x A100 80GB) and a CPU build (--no-gpu), both on head 21fbc82fe plus its patch: the four named Thermostat/wall-temperature tests 4 passed, 0 failed; Bodyforces 6/0; MeanDpdx 1/0; Thermostat 4/0, on both builds, against the existing goldens.
• Tioga (MI250X, AFAR amdflang, OpenMP offload), head 644a5918b: the two new thermostat tests passed.

What we found

• CI went red on GNU, Intel and NVHPC CPU when I replaced the Fypp code macros with plain functions: that edit also moved the source terms (g1, g2) from "recomputed every Runge-Kutta stage" to "stored once per step", which shifted the wall-temperature thermostat result by about 1e-4 against the golden. 21fbc82fe restores the per-stage numerics (one shared kernel s_compute_mts_terms) and keeps the macro removal; the goldens did not change.
• The thermostat tests use a constant time step and override_tol 1e-6: the Newton temperature solve (tol 1e-6) makes them differ across vendors, and the controller magnifies it. The thermal-channel golden was dropped for the same reason.
• 64-cell channel run, hot and cold wall: Nu 7.7699 and 7.3220, mean about 7.546, +0.07% from 7.541, stopped at t of about 2.62e-5 rather than run to steady state. The README for that example is not part of this PR yet.

CI on the current head was still running when this was written. (Hand tests run with Claude Code.)

[sbryngelson]

Reply to the Claude Code Review above (two medium findings):

1. Composition at the neighbour cells. Correct as read: in s_compute_mts_terms, Ys is built from cell j only, and lam(1), lam(3) use the neighbours' temperatures with that composition. This is what the original code did (the MTS_TERMS macro this replaced), and the applied source and the controller's A1/A2 share the same values, so they stay consistent. It is exact when composition is uniform in x, which is the thermally developed state this model assumes, and it is not exercised otherwise (the verified channel has constant composition). It is first-order inconsistent where composition varies in x. I did not change it: using the neighbour cell's composition needs species values at j-1 and j+1 in the halo, and I have not confirmed those are valid where s_update_mean_sources is called. I will add a comment stating the approximation with the next push, and I am happy to switch to neighbour-cell composition if a maintainer prefers.
2. get_mixture_thermal_conductivity_mixavg on the device. It is declared $:GPU_ROUTINE(function_name='get_mixture_thermal_conductivity_mixavg', parallelism='[seq]') in toolchain/mfc/thermochem/module.fpp.mako, and it is already called from device loops in m_chemistry.fpp (the diffusion kernel) and m_ibm.fpp. This kernel ran on an A100 (NVHPC 25.11, --gpu acc) with the Thermostat tests passing, and on Tioga's MI250X (AFAR amdflang) at the earlier commit with the same call. It has not run on Cray CCE yet; the Frontier lanes will be the first.

(Reply written with Claude Code.)

[codecov]

Codecov Report

❌ Patch coverage is 70.00000% with 21 lines in your changes missing coverage. Please review.
✅ Project coverage is 64.91%. Comparing base (9649d50) to head (ab9d69a).
⚠️ Report is 7 commits behind head on master.

Files with missing lines	Patch %	Lines
src/simulation/m_body_forces.fpp	69.11%	7 Missing and 14 partials ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #2200      +/-   ##
==========================================
+ Coverage   63.71%   64.91%   +1.20%     
==========================================
  Files          90      110      +20     
  Lines       24569    29538    +4969     
  Branches     3584     4481     +897     
==========================================
+ Hits        15653    19175    +3522     
- Misses       6135     6812     +677     
- Partials     2781     3551     +770     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[sbryngelson]

Hand test on NVIDIA after the merge with master (#2198 and #2204). wingtip, NVHPC 25.11 (nvfortran 25.11-0), 4x A100 80GB; a GPU build (--gpu acc) and a CPU build (--no-gpu) of the same commit, 191b51cac, against the existing goldens (none regenerated):

tests	GPU (acc)	CPU
the four Thermostat / wall-temperature tests (88ADD5C1 D7D6DEEB F8D466BD C02C399E)	4 passed, 0 failed	4 passed, 0 failed
Bodyforces	6 / 0	6 / 0
MeanDpdx	1 / 0	1 / 0
Thermostat	4 / 0	4 / 0
IBM Reacting Surface (8 tests, including #2198's lumped-wall cases next to this PR's thermostat cases)	8 / 0	not run

The conflict in toolchain/mfc/test/cases.py with #2198 was resolved by keeping both blocks as separate dicts; the 8 passing tests confirm that. The current head fee66d26d adds a merge of #2169 (ghost pressure for spinning and moving bodies), which touches no code path used by these fixed-body, periodic-box tests; I did not rerun on it by hand and rely on CI for the IBM set. Not yet run on Cray CCE.

(Hand test run with Claude Code.)

[github-actions]

Lines of Code

File	Lines	Diff
src/simulation/m_body_forces.fpp	606	+90
src/simulation/m_global_parameters.fpp	946	+1

Directory	Lines	Diff
simulation	39320	+91
total	59267	+91