[WIP] cache Least Squares matrix for gradient computation

Common/include/CConfig.hpp

@@ -9733,6 +9733,14 @@ class CConfig {
    */
   bool GetIntegrated_HeatFlux() const { return Integrated_HeatFlux; }
 
+  /*!
+   * \brief Whether the factorized least-squares gradient metric terms are cached and reused
+   *        across evaluations. This is the default behavior, except for periodic boundaries
+   *        (their metric and RHS accumulations are fused in one exchange) and the discrete
+   *        adjoint (the coordinate dependence of the metrics must remain on the tape).
+   */
+  bool GetLSQMetricCaching() const { return (nMarker_PerBound == 0) && !DiscreteAdjoint; }
+
   /*!
    * \brief Get Compute Average.
    * \return YES if start computing averages

Common/include/geometry/CGeometry.hpp

@@ -217,6 +217,12 @@ class CGeometry {
   unsigned long edgeColorGroupSize{1};     /*!< \brief Size of the edge groups within each color. */
   unsigned long elemColorGroupSize{1};     /*!< \brief Size of the element groups within each color. */
 
+  /*--- Cached least-squares gradient metric terms (see computeGradientsLeastSquares.hpp). ---*/
+
+  su2activematrix LSQMetricCache[2];            /*!< \brief Cached LSQ metrics S = inv(A), upper triangle stored
+                                                            row-wise, [0] unweighted, [1] inverse-distance weighted. */
+  bool LSQMetricCacheValid[2] = {false, false}; /*!< \brief Validity of the cached LSQ metrics per weighting. */
+
   ColMajorMatrix<uint8_t> CoarseGridColor_; /*!< \brief Coarse grid levels, colorized. */
 
  public:
@@ -1919,6 +1925,28 @@ class CGeometry {
    */
   void SetNaturalEdgeColoring();
 
+  /*!
+   * \brief Get the cached least-squares metric terms (S = inv(A), upper triangle row-wise).
+   * \param[in] weighted - False for unweighted, true for inverse-distance weighting.
+   */
+  inline su2activematrix& GetLSQMetricCache(bool weighted) { return LSQMetricCache[weighted]; }
+  inline const su2activematrix& GetLSQMetricCache(bool weighted) const { return LSQMetricCache[weighted]; }
+
+  /*!
+   * \brief Check whether the cached least-squares metric terms are valid for a weighting.
+   */
+  inline bool LSQMetricCacheIsValid(bool weighted) const { return LSQMetricCacheValid[weighted]; }
+
+  /*!
+   * \brief Declare the cached least-squares metric terms valid for a weighting.
+   */
+  inline void SetLSQMetricCacheValid(bool weighted) { LSQMetricCacheValid[weighted] = true; }
+
+  /*!
+   * \brief Invalidate the cached least-squares metric terms (e.g. if node coordinates change).
+   */
+  inline void InvalidateLSQMetricCache() { LSQMetricCacheValid[0] = LSQMetricCacheValid[1] = false; }
+
   /*!
    * \brief Get the group size used in edge coloring.
    * \return Group size.

Common/src/CConfig.cpp

@@ -3886,6 +3886,7 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
       MUSCL_Flow = false;
     }
   }
+
   if (MUSCL_AdjFlow && (Kind_ConvNumScheme_AdjFlow == SPACE_CENTERED)) {
     if (OptionIsSet("MUSCL_ADJFLOW")) {
       SU2_MPI::Error("Centered schemes do not use MUSCL reconstruction (use MUSCL_ADJFLOW= NO).", CURRENT_FUNCTION);

Common/src/geometry/CMultiGridGeometry.cpp

@@ -1000,6 +1000,9 @@ void CMultiGridGeometry::MatchActuator_Disk(const CConfig* config) {
 
 void CMultiGridGeometry::SetControlVolume(const CGeometry* fine_grid, unsigned short action) {
   BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS {
+    /*--- Coarse coordinates change with the fine grid, any cached LSQ metrics are stale. ---*/
+    if (action != ALLOCATE) InvalidateLSQMetricCache();
+
     /*--- Compute the area of the coarse volume ---*/
     for (auto iCoarsePoint = 0u; iCoarsePoint < nPoint; iCoarsePoint++) {
       nodes->SetVolume(iCoarsePoint, 0.0);

Common/src/geometry/CPhysicalGeometry.cpp

@@ -6759,6 +6759,9 @@ void CPhysicalGeometry::SetControlVolume(CConfig* config, unsigned short action)
 
   BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS { /*--- The following is difficult to parallelize with threads. ---*/
 
+    /*--- Node coordinates changed, any cached LSQ gradient metrics are stale. ---*/
+    if (action != ALLOCATE) InvalidateLSQMetricCache();
+
     su2double my_DomainVolume = 0.0;
     for (auto iElem = 0ul; iElem < nElem; iElem++) {
       const auto nNodes = elem[iElem]->GetnNodes();