closed_form.hpp divides by the air deficit, not the leaf-to-air one (#93 drift)

[dfalster]

Found by a staleness sweep in #115, which records the divergence at the site but deliberately does not repair it.

The divergence

Every vapour-pressure deficit in inst/include/phylloptim/closed_form.hpp is atm_vpd_ — the deficit at air temperature:

closed_form.hpp:82    l.H2O_CO2_stom_diff_ratio_ * 1e-3 * assim * l.atm_vpd_ / ...
closed_form.hpp:137   l.atm_kpa_ * E * kg_to_mol_h2o / l.atm_vpd_ / ...
closed_form.hpp:172   l.H2O_CO2_stom_diff_ratio_ * 1e-3 * l.atm_vpd_ / ...
closed_form.hpp:203   l.atm_kpa_ * E * kg_to_mol_h2o / l.atm_vpd_ / ...

That is what Fick's law divided by when the header was written. #93 moved the live model to vpd_leaf_, the leaf-to-air deficit, which on the energy-balance path is 3–4× the air's over Tair 25–45. So on that path these formulae approximate a model that no longer exists.

Off the energy-balance path the two are equal by construction — set_leaf_vpd returns atm_vpd_ exactly, which is the arithmetic #93 arranged to keep the golden file bit-identical — so the prescribed-temperature path, which is the default and the only one this header has ever been benchmarked on, is unaffected.

Why this is not a rename

The obvious repair is atm_vpd_ → vpd_leaf_. It may not be valid, and that is the substance of this issue rather than a caveat on it.

vpd_leaf_ is a function of Tleaf; Tleaf is a function of E; and E is what these expressions solve for. The header's whole claim is that given λ the optimum collapses to the Medlyn form

ci/ca = xi/(xi + sqrt(D)),   xi = sqrt(Q/lambda)

with D a known constant. If D becomes D(Tleaf(E)) then D is no longer known when the inversion is performed, and the closed form is not obviously still closed. At beta2 = 1/stem_c, where the header currently reports 47× because "the leaf solve itself has essentially vanished", there would now be something left to solve.

So the possibilities are, at least:

1. the substitution works with one extra Picard pass over D, at some cost to the 47×;
2. the explicit form survives for a linearised esat(T) (the same linearisation plantecophys uses for radiation conductance — see Temperature-dependent outgoing longwave in the Penman-Monteith Rn (from plant #581 / #567 review) #28, Fixed longwave offset is what stops an energy-balance alignment being exact #97);
3. the closed form is simply a prescribed-temperature-path optimisation and should say so and be gated on !use_energy_balance_.

Option 3 is cheap and honest and may be the right answer, given the header's own STATUS says it is not wired in and nothing calls it by default.

Why it has not bitten

• Not wired into Leaf; nothing calls it by default (its own STATUS note).
• within_guard requires establishing accuracy per use before it can replace the exact solve on a production path, so any adoption is already gated on a comparison that would surface this.
• No test exercises it on the PM path, which is also why nothing failed when Leaf-to-air VPD, and Sperry's ProfitMax #93 landed.

What would settle it

Benchmark and accuracy-check the closed form against optimise_psi_stem_TF with the gate on, over the Tair 25–45 range where the deficits diverge, and see whether option 1's error and cost are acceptable. If not, take option 3 and gate it.

⚠️ Do not measure this by A/B-ing two builds run in sequence — hazard 5, and the header's own numbers were taken in the companion analysis rather than here.