[WIP] cache Least Squares matrix for gradient computation

[tbellosta]

Proposed Changes

This PR adds an opt-in option, CACHE_LSQ_METRICS (default NO), caches the factorized metric terms of the (weighted or unweighted) least-squares gradients. The LSQ normal matrix A depends only on the node coordinates and on the weighting, yet it is currently re-assembled and re-factorized on every gradient evaluation, every nonlinear iteration. With the option enabled:

• On the first evaluation the normal matrix is assembled and factorized as usual, and its inverse S = A^-1 (6 unique entries per node in 3D, 3 in 2D) is stored in CGeometry.
• Every subsequent evaluation only accumulates the right-hand side b = sum_k wd(u_k - u_i) in a loop over the edges (using the edge coloring; the contribution of an edge is identical for both of its end points, so each edge is visited once) followed by one small symmetric matrix-vector product S*b per node.
• The cache is stored per grid level with one slot per weighting, shared by all solvers (flow primitive and reconstruction gradients, turbulence/species/scalar solution gradients, auxiliary variables). The common practice of unweighted LSQ for MUSCL reconstruction plus inverse-distance-weighted LSQ for viscous gradients is fully supported: both slots are built once and reused; e.g. the turbulence solver directly reuses the metrics built by the flow solver.
• Moving/deforming grids are supported: the cache is invalidated in SetControlVolume(..., UPDATE) (the common point of all mesh motion/deformation paths, on fine and coarse MG levels) and rebuilt on the next evaluation, i.e. once per mesh update, so the savings scale with the number of inner iterations per time step.
• The option is automatically disabled (with a warning) for periodic boundaries (the periodic LSQ communication fuses the matrix and RHS accumulations; an RHS-only exchange would be needed) and for the discrete adjoint (the coordinate dependence of the metrics must remain on the tape).

Saved operations (per node per gradient evaluation, 3D, k edge-neighbors, N variables)

Eliminated from every evaluation after the first:

• normal-matrix assembly: ~21k flops (~300 flops per node on tetrahedra, k=14),
• per-node Cholesky factorization + inversion: ~35 flops and 3 sqrt (no sqrt/div remain in the iteration path),
• half of the geometric edge work (distance vectors and solution differences are computed once per edge instead of once per node-neighbor pair).

What remains is the RHS accumulation (~4N flops per node per edge) and one S*b product (18N flops per node). For the N=6 compressible primitives this roughly halves the gradient-kernel flops; for scalar solvers (N=1-2, e.g. SA) where the assembly dominates it is a ~4-5x reduction. Measured end-to-end (serial, linear-solver work held fixed): ~3.3% of total wall time on an inviscid ONERA M6 (57.5k nodes, tets, one gradient set per iteration) and ~6% on a 2D RANS-SA case (three LSQ gradient sets per iteration). Cases dominated by the linear solver will see proportionally less.

Memory footprint

The cache adds nDim*(nDim+1)/2 su2doubles per node per weighting used (48 B/node per weighting in 3D, 24 B/node in 2D; at most two weightings), allocated lazily per grid level and shared by all solvers. For reference, every solver already allocates an nDim*nDim Rmatrix scratch (72 B/node in 3D), so for a typical RANS case the addition is small compared to the existing gradient machinery and negligible compared to overall solver memory.

Validation

With the option enabled, results are identical to develop (to output precision) for:

• inviscid NACA0012, unweighted LSQ reconstruction, 200 iterations (serial and 2 MPI ranks),
• RANS-SA NACA0012 with NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES and NUM_METHOD_GRAD_RECON= LEAST_SQUARES (all residuals including the turbulence variable),
• an unsteady pitching NACA0012 with GRID_MOVEMENT= RIGID_MOTION (dual time stepping), exercising the invalidation/rebuild path,
• hybrid-parallel runs (OpenMP edge-colored loop, bit-reproducible run-to-run; falls back to a race-free node loop if the natural coloring is in use).

config_template.cfg documents the new option. Default behavior (CACHE_LSQ_METRICS= NO) is bit-identical to current develop.

• I am submitting my contribution to the develop branch.
• My contribution generates no new compiler warnings (try with --warnlevel=3 when using meson).
• My contribution is commented and consistent with SU2 style (https://su2code.github.io/docs_v7/Style-Guide/).
• I used the pre-commit hook to prevent dirty commits and used pre-commit run --all to format old commits.
• I have added a test case that demonstrates my contribution, if necessary.
• I have updated appropriate documentation (Tutorials, Docs Page, config_template.cfg), if necessary.

[pcarruscag]

Nice. Make it the default, compute in geometry preprocessing when necessary, no opt-in