Add Mercier gradient with respect to the flux slopes

[krystophny]

Summary

Completes the analytic Mercier scalar-slope gradient by adding ∂d_mercier/∂toroidal_flux_slope and ∂d_mercier/∂poloidal_flux_slope (pressure_slope and iota_slope were already done). Both flux slopes are genuine paths in the objective gradient — a boundary/profile change moves the equilibrium's flux derivatives.

• toroidal_flux_slope (=flux_slope) couples through nearly every term: psi_slope (→ shear, iota_psi, dp_dpsi, d2v_dpsi2, current_slope_ratio), grad_psi = |psi_slope|·area/|J| (in every integral, ∂grad_psi/∂flux = grad_psi/flux), AND the β denominator. poloidal_flux_slope enters only the β denominator.
• The β-mediated part reuses the solved β harmonics with no extra solve: ∂β_mn/∂flux = −β_mn·(∂denom_mn/∂flux)/denom_mn, with ∂denom/∂φ_t = −2πn, ∂denom/∂φ_p = +2πm, zeroing the same resonant modes.
• New result fields d_mercier_d_toroidal_flux_slope(:), d_mercier_d_poloidal_flux_slope(:); existing fields and output kept.

Verification

Before

The flux-slope gradients and their test do not exist on main.

After

FD (re-running grad_psi, the β solve, and the integrals at the perturbed slope) matches the analytic full gradient:

analytic d/d toroidal   =   1.4719140325233543E+01   fd =   1.4719140336072467E+01   rel err = 7.4E-10
analytic d/d poloidal   =   2.7511213297157305E-02   fd =   2.7511213073005823E-02   rel err = 8.1E-09
PASS

ctest -R mercier → 4/4 (mercier_diagnostic, mercier_gradient, mercier_full_gradient, mercier_fluxslope_gradient); full suite 38/38. Behaviour-preserving: gliss_mercier on W7X_gst_ns64_MN12.nc byte-identical to pristine main (md5 unchanged). -Werror, fo lint, fo fmt all clean.