fix(phase_change): bound Newton relaxation loops by max_iter

[sbryngelson]

Description

The three Newton solvers in src/common/m_phase_change.fpp — s_infinite_pt_relaxation_k, s_infinite_ptg_relaxation_k, and s_TSat — looped on a convergence-only do while condition with no iteration limit. The module parameter max_iter was declared but never referenced anywhere (grep-confirmed): it has been dead since the original phase-change PR (#179), which declared it — and even added it to the !$acc declare list — but never wired it into the loop conditions. A later cleanup (#999) removed the stale GPU declaration while migrating to the macro, leaving a fully orphaned declaration. So it is a latent birth defect, not a regression.

When a cell fails to converge, the loop never exits and the run hangs indefinitely (on GPU the kernel never returns and the host blocks in cuStreamSynchronize). Reproduced with the stock, unmodified examples/2D_phasechange_bubble, which hangs at t_step = 8 independent of compiler, MPI, and grid resolution.

Root cause of the non-convergence itself (tracked as a separate follow-up): s_TSat and the pTg solver compare a dimensional Gibbs-free-energy residual (natural magnitude ~10²–10³) against the fixed tolerance ptgalpha_eps = 1e-2 under heavy Om = 1e-3 underrelaxation, so for low-pressure vapor cells the exit condition is effectively unreachable. This PR does not change that convergence behavior — it only removes the ability to hang.

This change:

• Bounds each of the three loops with if (ns >= max_iter) exit, accepting the last iterate on exhaustion. An unbounded hang becomes a bounded, finite step.
• Reduces the misleading real-literal initialization max_iter = 1e8_wp (an integer initialized from a real literal) to a clean integer 100000, comfortably above the observed legitimate convergence ceiling.

Closes #1646.

Type of change

• Bug fix

Testing

• All 18 phase-change regression tests (Phase Change model 5/6, 1D/2D/3D, 2/3 fluids, model_eqns 2/3) pass byte-identical against the existing golden files. This confirms the cap sits above the legitimate convergence ceiling and truncates no converging cell — one case (2D, model 6) legitimately iterates deep (~145 s vs ~7 s for its siblings) yet still matches bit-for-bit, so the cap is neither too low (clipping valid slow cells) nor merely cosmetic.
• Reproduced the original hang on examples/2D_phasechange_bubble (freezes at t_step = 8); with this fix the same case advances past it and continues stepping steadily (ran to t_step = 31).
• ./mfc.sh precheck passes (formatting, spell, lint, docs, params, case validation).

CPU-only change in src/common/; no GPU code paths altered.

Checklist

• I added or updated tests for new behavior
• I updated documentation if user-facing behavior changed

The change is golden-safe (no behavior change for converging cells), so it adds no new golden test; a dedicated "does not hang" regression would require a case that loops unboundedly without the fix. Happy to add one if reviewers prefer.

[Copilot]

Pull request overview

This PR addresses a solver hang in the phase-change module by bounding previously unbounded Newton iteration loops, ensuring non-convergent cells cannot trap the run in an infinite loop (including on GPU where the kernel may never return).

Changes:

• Replace the unused/incorrect max_iter initialization (integer from a real literal) with a clean integer constant (100000).
• Add max_iter exhaustion guards to the three Newton solvers (s_infinite_pt_relaxation_k, s_infinite_ptg_relaxation_k, s_TSat) so they exit after too many iterations and accept the last iterate.

[sbryngelson]

@copilot good catch on the (pS - pO)/pO division — you're right that it's evaluated on the first pass (Fortran .or. isn't short-circuit) with pO = 0. It's benign in this PR (the numerator is nonzero, so it's +Inf not NaN, and the .or. (ns == 0) guard forces the first iteration regardless), and it's pre-existing rather than introduced here — but it would trip an FPE-trapping build, so it's worth removing.

Fixed in the follow-up #1654 (which reworks these solvers): the relative criterion is rewritten in multiply form, abs(pS - pO) > (palpha_eps/1.e4_wp)*abs(pO), which is mathematically identical for the exit decision and never divides by pO.

[codecov]

Codecov Report

❌ Patch coverage is 0% with 3 lines in your changes missing coverage. Please review.
✅ Project coverage is 59.54%. Comparing base (e090d56) to head (208e9c3).
⚠️ Report is 5 commits behind head on master.

Files with missing lines	Patch %	Lines
src/common/m_phase_change.fpp	0.00%	0 Missing and 3 partials ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #1653      +/-   ##
==========================================
- Coverage   60.95%   59.54%   -1.42%     
==========================================
  Files          83       83              
  Lines       20003    21114    +1111     
  Branches     2983     3131     +148     
==========================================
+ Hits        12193    12572     +379     
- Misses       5782     6437     +655     
- Partials     2028     2105      +77     

☔ 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.