Fix multigrid early exit functionality (#2831)

Base early exit on proper residual
fix some implicit line agglomeration
add more debugging info
---------

Co-authored-by: Pedro Gomes <pcarruscag@gmail.com>
Co-authored-by: Copilot Autofix powered by AI <62310815+github-advanced-security[bot]@users.noreply.github.com>

Author: 3 people

Common/include/CConfig.hpp

@@ -1147,6 +1147,8 @@ class CConfig {
   /*--- Multigrid options  ---*/
   unsigned short nMG_PreSmooth_p{0}, nMG_PostSmooth_p{0}, nMG_CorrecSmooth_p{0};
   unsigned short *MG_PreSmooth_p{nullptr}, *MG_PostSmooth_p{nullptr}, *MG_CorrecSmooth_p{nullptr};
+  unsigned short nMG_CflScaling_p{0};
+  su2double *MG_CflScaling_p{nullptr};
 
   ENUM_STREAMWISE_PERIODIC Kind_Streamwise_Periodic; /*!< \brief Kind of Streamwise periodic flow (pressure drop or massflow) */
   bool Streamwise_Periodic_Temperature;              /*!< \brief Use real periodicity for Energy equation or otherwise outlet source term. */

Common/include/option_structure.hpp

@@ -1117,9 +1117,12 @@ struct CMGOptions {
   std::vector<unsigned short> MG_PreSmooth;    /*!< \brief Multigrid pre-smoothing iterations per level. */
   std::vector<unsigned short> MG_PostSmooth;   /*!< \brief Multigrid post-smoothing iterations per level. */
   std::vector<unsigned short> MG_CorrecSmooth; /*!< \brief Multigrid Jacobi correction-smoothing per level. */
+  std::vector<su2double> MG_CflScaling;        /*!< \brief Per-level CFL scaling factors relative to the previous (finer) level. Entry [i] scales level i+1 from level i. Size = nMGLevels. */
   bool MG_Smooth_EarlyExit{false};        /*!< \brief Enable early exit for MG smoothing iterations. */
   bool MG_Smooth_Output{false};           /*!< \brief Output compact per-cycle smoothing summary. */
+  su2double MG_Smooth_StagnationTol{0.0}; /*!< \brief Stagnation early exit: stop if current_rms >= prev_rms * tol. 0 = disabled. */
   bool MG_Implicit_Lines{false};          /*!< \brief Enable implicit-lines agglomeration from walls. */
+  unsigned long MG_Implicit_Lines_MaxLength{20}; /*!< \brief Maximum nodes on a wall-normal implicit line (including wall seed). */
 };
 
 /*!

Common/src/CConfig.cpp

@@ -1985,16 +1985,22 @@ void CConfig::SetConfig_Options() {
   addDoubleOption("MG_DAMP_PROLONGATION", Damp_Correc_Prolong, 0.5);
   /*!\brief MG_SMOOTH_EARLY_EXIT\n DESCRIPTION: Enable early exit for MG smoothing when RMS drops below threshold. DEFAULT: NO \ingroup Config*/
   addBoolOption("MG_SMOOTH_EARLY_EXIT", MGOptions.MG_Smooth_EarlyExit, true);
-  /*!\brief MG_SMOOTH_RES_THRESHOLD\n DESCRIPTION: Smoothing stops when current_rms < threshold * initial_rms. DEFAULT: 0.1 \ingroup Config*/
-  addDoubleOption("MG_SMOOTH_RES_THRESHOLD", MGOptions.MG_Smooth_Res_Threshold, 0.5);
+  /*!\brief MG_SMOOTH_RES_THRESHOLD\n DESCRIPTION: Early exit smoothing when current_rms drops below threshold * initial_rms. DEFAULT: 0.9 \ingroup Config*/
+  addDoubleOption("MG_SMOOTH_RES_THRESHOLD", MGOptions.MG_Smooth_Res_Threshold, 0.9);
   /*!\brief MG_SMOOTH_OUTPUT\n DESCRIPTION: Print compact per-cycle smoothing iteration summary. DEFAULT: NO \ingroup Config*/
   addBoolOption("MG_SMOOTH_OUTPUT", MGOptions.MG_Smooth_Output, false);
+  /*!\brief MG_SMOOTH_STAGNATION_TOL\n DESCRIPTION: Stop smoothing if current_rms >= previous_rms * this value. Values < 1.0 enable early exit on stagnation, 1.0 only exits on defect growth. DEFAULT: 0.99 \ingroup Config*/
+  addDoubleOption("MG_SMOOTH_STAGNATION_TOL", MGOptions.MG_Smooth_StagnationTol, 0.99);
   /*!\brief MG_SMOOTH_COEFF\n DESCRIPTION: Smoothing coefficient for the correction prolongation Jacobi smoother. DEFAULT: 1.25 \ingroup Config*/
   addDoubleOption("MG_SMOOTH_COEFF", MGOptions.MG_Smooth_Coeff, 1.25);
   /*!\brief MG_MIN_MESHSIZE\n DESCRIPTION: Minimum number of CVs on the coarsest multigrid level. Levels that would produce fewer CVs are not created. DEFAULT: 50 \ingroup Config*/
   addUnsignedLongOption("MG_MIN_MESHSIZE", MGOptions.MG_Min_MeshSize, 500);
   /*!\brief MG_IMPLICIT_LINES\n DESCRIPTION: Enable agglomeration along implicit lines from wall seeds. DEFAULT: NO \ingroup Config*/
   addBoolOption("MG_IMPLICIT_LINES", MGOptions.MG_Implicit_Lines, false);
+  /*!\brief MG_IMPLICIT_LINES_MAX_LENGTH\n DESCRIPTION: Maximum number of nodes on a wall-normal implicit agglomeration line (including the wall seed node). DEFAULT: 20 \ingroup Config*/
+  addUnsignedLongOption("MG_IMPLICIT_LINES_MAX_LENGTH", MGOptions.MG_Implicit_Lines_MaxLength, 20);
+  /*!\brief MG_CFL_SCALING\n DESCRIPTION: Per-level CFL scaling factors for coarse MG levels. Entry i is the ratio CFL(i+1)/CFL(i). If fewer values than nMGLevels are given, the last value is repeated. DEFAULT: 0.25 (i.e., 1/4 per level) \ingroup Config*/
+  addDoubleListOption("MG_CFL_SCALING", nMG_CflScaling_p, MG_CflScaling_p);
 
   /*!\par CONFIG_CATEGORY: Spatial Discretization \ingroup Config*/
   /*--- Options related to the spatial discretization ---*/
@@ -4789,6 +4795,16 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
                [](unsigned short  ) { return (unsigned short)0; });
     fillSmooth(nMG_CorrecSmooth_p, MG_CorrecSmooth_p, MGOptions.MG_CorrecSmooth,
                [](unsigned short  ) { return (unsigned short)0; });
+
+    /*--- Fill MG_CflScaling to size nMGLevels (one entry per coarse level transition). ---*/
+    MGOptions.MG_CflScaling.resize(nMGLevels);
+    if (nMG_CflScaling_p != 0) {
+      for (unsigned short i = 0; i < nMGLevels; ++i)
+        MGOptions.MG_CflScaling[i] = (i < nMG_CflScaling_p) ? MG_CflScaling_p[i] : MG_CflScaling_p[nMG_CflScaling_p - 1];
+    } else {
+      for (unsigned short i = 0; i < nMGLevels; ++i)
+        MGOptions.MG_CflScaling[i] = 0.25;
+    }
   }
 
   /*--- Override MG Smooth parameters ---*/
@@ -7534,10 +7550,12 @@ void CConfig::SetOutput(SU2_COMPONENT val_software, unsigned short val_izone) {
         MGTable.AddColumn("Presmooth",     10);
         MGTable.AddColumn("PostSmooth",    10);
         MGTable.AddColumn("CorrectSmooth", 10);
+        MGTable.AddColumn("CFL Scaling",   12);
         MGTable.SetAlign(PrintingToolbox::CTablePrinter::RIGHT);
         MGTable.PrintHeader();
         for (unsigned short iLevel = 0; iLevel < nMGLevels+1; iLevel++) {
-          MGTable << iLevel << MGOptions.MG_PreSmooth[iLevel] << MGOptions.MG_PostSmooth[iLevel] << MGOptions.MG_CorrecSmooth[iLevel];
+          const string cflStr = (iLevel == 0) ? "-" : std::to_string(MGOptions.MG_CflScaling[iLevel-1]).substr(0,6);
+          MGTable << iLevel << MGOptions.MG_PreSmooth[iLevel] << MGOptions.MG_PostSmooth[iLevel] << MGOptions.MG_CorrecSmooth[iLevel] << cflStr;
         }
         MGTable.PrintFooter();
       }

Common/src/geometry/CMultiGridGeometry.cpp

@@ -1280,24 +1280,28 @@ su2double CMultiGridGeometry::ComputeLocalCurvature(const CGeometry* fine_grid,
 void CMultiGridGeometry::AgglomerateImplicitLines(unsigned long& Index_CoarseCV, const CGeometry* fine_grid,
                                                   const CConfig* config, CMultiGridQueue& MGQueue_InnerCV) {
   /*--- Parameters ---*/
-  const su2double ANGLE_THRESHOLD_DEG = 20.0;   /*!< Stop line if direction deviates more than this. */
-  constexpr unsigned long MAX_LINE_LENGTH = 20; /*!< Max nodes on implicit line (including wall). */
+  const su2double ANGLE_THRESHOLD_DEG = 20.0; /*!< Stop line if direction deviates more than this. */
+  const unsigned long MAX_LINE_LENGTH = config->GetMGOptions().MG_Implicit_Lines_MaxLength;
   const su2double cos_threshold = cos(ANGLE_THRESHOLD_DEG * PI_NUMBER / 180.0);
 
   const unsigned long nPointFine = fine_grid->GetnPoint();
 
-  /*--- Collect implicit lines starting at wall vertices.
+  /*--- Collect implicit lines starting at viscous (no-slip) wall vertices only.
+   *    Seeding from non-wall boundaries (farfield, inlet, outlet, symmetry) would
+   *    claim interior BL cells before the wall lines can reach them, leaving
+   *    wall-seeded lines with length < 3 (discarded).  Restricting to viscous walls
+   *    ensures the boundary-layer cells are agglomerated wall-first.
    *    Each line: [wall_node, interior_1, interior_2, ...].
    *    The wall node (index 0) is already agglomerated by boundary agglomeration;
    *    only interior nodes (index >= 1) are paired into coarse CVs. ---*/
   vector<vector<unsigned long>> lines;
 
   for (auto iMarker = 0u; iMarker < fine_grid->GetnMarker(); iMarker++) {
-    /*--- Skip non-physical markers ---*/
-    if (config->GetMarker_All_KindBC(iMarker) == SEND_RECEIVE ||
-        config->GetMarker_All_KindBC(iMarker) == INTERNAL_BOUNDARY ||
-        config->GetMarker_All_KindBC(iMarker) == NEARFIELD_BOUNDARY)
-      continue;
+    /*--- Only seed lines from viscous (no-slip) wall markers.
+     *    Non-wall boundaries (farfield, inlet, outlet, symmetry) must NOT seed
+     *    lines because they would prematurely claim boundary-layer interior nodes. ---*/
+    const auto bc = config->GetMarker_All_KindBC(iMarker);
+    if (bc != HEAT_FLUX && bc != ISOTHERMAL && bc != CHT_WALL_INTERFACE && bc != SMOLUCHOWSKI_MAXWELL) continue;
 
     for (auto iVertex = 0ul; iVertex < fine_grid->GetnVertex(iMarker); iVertex++) {
       const auto iPoint = fine_grid->vertex[iMarker][iVertex]->GetNode();

SU2_CFD/include/integration/CMultiGridIntegration.hpp

@@ -217,90 +217,64 @@ class CMultiGridIntegration final : public CIntegration {
                      unsigned short RunTime_EqSystem, unsigned long Iteration, unsigned short iZone);
 
   /*!
-   * \brief Compute adaptive CFL for multigrid coarse levels.
-   * \param[in] config - Problem configuration.
-   * \param[in] solver_coarse - Coarse grid solver.
-   * \param[in] geometry_coarse - Coarse grid geometry.
-   * \param[in] iMesh - Current multigrid level.
-   * \param[in] CFL_fine - Fine grid CFL value (passive).
-   * \param[in] CFL_coarse_current - Current coarse grid CFL value (passive).
-   * \param[in] rms_res_coarse - Coarse-grid RMS residual (already MPI-reduced, from lastPreSmoothRMS).
-   * \return New CFL value for the coarse grid.
-   */
-  passivedouble computeMultigridCFL(CConfig* config, unsigned short iMesh,
-                                     passivedouble CFL_fine, passivedouble CFL_coarse_current,
-                                     passivedouble rms_res_coarse);
-
-  /*!
-   * \brief Adapt the residual restriction damping factor.
+   * \brief Adapt both restriction and prolongation damping factors from the global-trend signal.
    *
-   * Uses \c lastPreSmoothIters[] (filled by the previous multigrid cycle) to assess
-   * whether the pre-smoother is converging fast or slow on coarse levels, then adjusts
-   * \c Damp_Res_Restric in \p config accordingly.
+   * Uses the cross-cycle EMA ratio (crossCycleRatio = fine_d0 / EMA(fine_d0)) to detect
+   * long-term convergence or divergence, then adjusts both \c Damp_Res_Restric and
+   * \c Damp_Correc_Prolong with a single shared signal.  The EMA filters per-cycle noise;
+   * no per-level aggregation or floor counter is needed.
    *
-   * Signal logic:
-   *  - any coarse level ran its full configured iterations: reduce damping
-   *  - all coarse levels exited early: increase damping
-   *  - mixed (some full, some partial): no change
-   *
-   * \param[in,out] config - Problem configuration.
+   * \param[in,out] config          - Problem configuration.
+   * \param[in]     crossCycleRatio - Current fine_d0 divided by the EMA of fine_d0.
    */
-  void adaptRestrictionDamping(CConfig* config);
+  void adaptDampingFactors(CConfig* config, passivedouble crossCycleRatio);
 
   /*!
-   * \brief Adapt the correction prolongation damping factor.
-   *
-   * Uses \c lastCorrecSmoothIters[] (filled by the previous multigrid cycle) to assess
-   * whether the correction smoother is struggling or converging fast,
-   * then adjusts \c Damp_Correc_Prolong in \p config accordingly.
-   *
-   * Signal logic:
-   *  - any level ran its full correction-smooth iterations: reduce damping
-   *  - all levels exited early: increase damping
-   *  - mixed: no change
-   *
-   * \param[in,out] config - Problem configuration; \c SetDamp_Correc_Prolong is called to persist the result.
+   * \brief Helper function for early-exit logic during pre/post-smoothing.
+   * \param[in] iSmooth - Current smoothing iteration index.
+   * \param[in] iMesh - Index of the mesh in multigrid computations.
+   * \param[in] defect - Current RMS defect value.
+   * \param[in] mgOpts - Reference to multigrid options.
+   * \param[in] stag_tol - Stagnation tolerance value.
+   * \param[in] early_exit - Whether early exit is enabled.
+   * \param[out] lastRMS - Array to store RMS values [start, end].
+   * \param[out] exitReason - Character for early exit reason ('T', 'S', 'A', or ' ').
+   * \param[out] worstStepRatio - Worst step-to-step ratio seen.
+   * \param[out] worstStep - Iteration number of worst step.
    */
-  void adaptProlongationDamping(CConfig* config);
+  void prePostEarlyExit(unsigned short iSmooth, unsigned short iMesh,
+                        passivedouble defect, const CMGOptions& mgOpts,
+                        passivedouble stag_tol, bool early_exit,
+                        passivedouble lastRMS[2], char& exitReason,
+                        passivedouble& worstStepRatio, unsigned short& worstStep);
 
-  /*--- CFL adaptation state variables.
-   *    These must be passivedouble: AD::Reset() clears the tape between adjoint recordings,
-   *    but class members survive. If these were su2double their stale AD indices would
-   *    reference the cleared tape, causing invalid memory access during the backward pass. ---*/
   static constexpr int MAX_MG_LEVELS = 10;
-  passivedouble current_avg[MAX_MG_LEVELS] = {};
-  passivedouble prev_avg[MAX_MG_LEVELS] = {};
-  passivedouble last_res[MAX_MG_LEVELS] = {};
-  bool last_was_increase[MAX_MG_LEVELS] = {};
-  int oscillation_count[MAX_MG_LEVELS] = {};
-  unsigned long last_check_iter[MAX_MG_LEVELS] = {};
-  unsigned long last_update_iter[MAX_MG_LEVELS] = {};
-  unsigned long last_reset_iter = std::numeric_limits<unsigned long>::max();
 
   /*--- Early-exit smoothing state (shared across OMP threads via master write + barrier). ---*/
-  bool mg_early_exit_flag = false;             /*!< \brief Shared flag for early exit across OMP threads. */
-  passivedouble mg_initial_smooth_rms = 0.0;  /*!< \brief Initial RMS before current smoothing phase. */
-  passivedouble mg_last_smooth_rms = 0.0;     /*!< \brief Last computed RMS; cached to avoid redundant Allreduce. */
+  bool mg_early_exit_flag = false;              /*!< \brief Shared flag for early exit across OMP threads. */
+  passivedouble mg_initial_smooth_rms = 0.0; /*!< \brief Initial RMS residual before current smoothing phase (FAS). */
+  passivedouble mg_prev_smooth_rms = 0.0;    /*!< \brief RMS residual from previous smoothing step; used for stagnation detection. */
+  passivedouble mg_fine_rms_ema = 0.0;      /*!< \brief EMA of fine-grid pre-smooth RMS across cycles; cross-cycle trend signal. */
+  passivedouble last_crossCycleRatio = 1.0; /*!< \brief crossCycleRatio from the most recent cycle; stored for display only. */
 
   /*--- Actual iteration counts per MG level, filled each cycle for the compact output summary. ---*/
   unsigned short lastPreSmoothIters[MAX_MG_LEVELS+1] = {};
   unsigned short lastPostSmoothIters[MAX_MG_LEVELS+1] = {};
   unsigned short lastCorrecSmoothIters[MAX_MG_LEVELS+1] = {};
+  /*--- Early-exit reason per level: 'T'=threshold, 'S'=stagnation, ' '=ran to completion. ---*/
+  char lastPreSmoothExitReason[MAX_MG_LEVELS+1]  = {};
+  char lastPostSmoothExitReason[MAX_MG_LEVELS+1] = {};
 
-  /*--- Per-level residual progress flags: true if the final RMS after that phase was lower
-   *    than the initial RMS.  Used by the adaptive damping routines to distinguish
-   *    "hit max iters but still converging" from "hit max iters and stagnated". ---*/
-  bool lastPreSmoothProgress[MAX_MG_LEVELS+1] = {};
-  bool lastPostSmoothProgress[MAX_MG_LEVELS+1] = {};
-  bool lastCorrecSmoothProgress[MAX_MG_LEVELS+1] = {};
-
-  /*--- Per-level start/end RMS for the compact output summary.
-   *    [0] = initial RMS before smoothing, [1] = final RMS after smoothing.
-   *    Filled unconditionally (early-exit path and exhaustion path).
-   *    Must be passivedouble: class members survive tape resets; su2double would
-   *    carry stale AD indices referencing a cleared tape. ---*/
+  /*--- Per-level start/end RMS residual for adaptive damping. ---*/
   passivedouble lastPreSmoothRMS[MAX_MG_LEVELS+1][2] = {};
   passivedouble lastPostSmoothRMS[MAX_MG_LEVELS+1][2] = {};
   passivedouble lastCorrecSmoothRMS[MAX_MG_LEVELS+1][2] = {};
 
+  /*--- Per-level worst step-to-step amplification seen inside a smoothing phase.
+   *    step==0 means no intra-smoother ratio was available (fewer than 2 sweeps). ---*/
+  passivedouble lastPreSmoothWorstStepRatio[MAX_MG_LEVELS+1] = {};
+  passivedouble lastPostSmoothWorstStepRatio[MAX_MG_LEVELS+1] = {};
+  unsigned short lastPreSmoothWorstStep[MAX_MG_LEVELS+1] = {};
+  unsigned short lastPostSmoothWorstStep[MAX_MG_LEVELS+1] = {};
+
 };