Put loop over iter in MLMG smooth (#4352)

AlexanderSinn · web-flow · commit eae3e36137ad · 2025-03-25T09:52:19.000-07:00
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H b/Src/LinearSolvers/MLMG/AMReX_MLCellLinOp.H
@@ -83,7 +83,7 @@ public:
     void apply (int amrlev, int mglev, MF& out, MF& in, BCMode bc_mode,
                         StateMode s_mode, const MLMGBndryT<MF>* bndry=nullptr) const override;
     void smooth (int amrlev, int mglev, MF& sol, const MF& rhs,
-                         bool skip_fillboundary=false) const final;
+                 bool skip_fillboundary, int niter) const final;
 
     void solutionResidual (int amrlev, MF& resid, MF& x, const MF& b,
                                    const MF* crse_bcdata=nullptr) override;
@@ -1211,15 +1211,17 @@ MLCellLinOpT<MF>::apply (int amrlev, int mglev, MF& out, MF& in, BCMode bc_mode,
 template <typename MF>
 void
 MLCellLinOpT<MF>::smooth (int amrlev, int mglev, MF& sol, const MF& rhs,
-                          bool skip_fillboundary) const
+                          bool skip_fillboundary, int niter) const
 {
     BL_PROFILE("MLCellLinOp::smooth()");
-    for (int redblack = 0; redblack < 2; ++redblack)
-    {
-        applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Solution,
-                nullptr, skip_fillboundary);
-        Fsmooth(amrlev, mglev, sol, rhs, redblack);
-        skip_fillboundary = false;
+    for (int i = 0; i < niter; ++i) {
+        for (int redblack = 0; redblack < 2; ++redblack)
+        {
+            applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Solution,
+                    nullptr, skip_fillboundary);
+            Fsmooth(amrlev, mglev, sol, rhs, redblack);
+            skip_fillboundary = false;
+        }
     }
 }
 
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLCurlCurl.H b/Src/LinearSolvers/MLMG/AMReX_MLCurlCurl.H
@@ -73,7 +73,7 @@ public:
                 StateMode s_mode, const MLMGBndryT<MF>* bndry=nullptr) const override;
 
     void smooth (int amrlev, int mglev, MF& sol, const MF& rhs,
-                 bool skip_fillboundary=false) const override;
+                 bool skip_fillboundary, int niter) const override;
 
     void solutionResidual (int amrlev, MF& resid, MF& x, const MF& b,
                            const MF* crse_bcdata=nullptr) override;
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLCurlCurl.cpp b/Src/LinearSolvers/MLMG/AMReX_MLCurlCurl.cpp
@@ -319,28 +319,29 @@ MLCurlCurl::apply (int amrlev, int mglev, MF& out, MF& in, BCMode /*bc_mode*/,
 }
 
 void MLCurlCurl::smooth (int amrlev, int mglev, MF& sol, const MF& rhs,
-                         bool skip_fillboundary) const
+                         bool skip_fillboundary, int niter) const
 {
     AMREX_ASSERT(rhs[0].nGrowVect().allGE(1));
 
     applyBC(amrlev, mglev, const_cast<MF&>(rhs), CurlCurlStateType::b);
-
 #if (AMREX_SPACEDIM == 1)
     int ncolors = 2;
 #else
     int ncolors = 4;
 #endif
 
-    for (int color = 0; color < ncolors; ++color) {
-        if (!skip_fillboundary) {
-            applyBC(amrlev, mglev, sol, CurlCurlStateType::x);
-        }
-        skip_fillboundary = false;
+    for (int i = 0; i < niter; ++i) {
+        for (int color = 0; color < ncolors; ++color) {
+            if (!skip_fillboundary) {
+                applyBC(amrlev, mglev, sol, CurlCurlStateType::x);
+            }
+            skip_fillboundary = false;
 #if (AMREX_SPACEDIM == 1)
-        smooth1D(amrlev, mglev, sol, rhs, color);
+            smooth1D(amrlev, mglev, sol, rhs, color);
 #else
-        smooth4(amrlev, mglev, sol, rhs, color);
+            smooth4(amrlev, mglev, sol, rhs, color);
 #endif
+        }
     }
 }
 
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLLinOp.H b/Src/LinearSolvers/MLMG/AMReX_MLLinOp.H
@@ -355,9 +355,10 @@ public:
      * \param sol               unknowns
      * \param rhs               RHS
      * \param skip_fillboundary flag controlling whether ghost cell filling can be skipped.
+     * \param niter number of (multicolor) smoothing iterations
      */
     virtual void smooth (int amrlev, int mglev, MF& sol, const MF& rhs,
-                         bool skip_fillboundary=false) const = 0;
+                         bool skip_fillboundary, int niter) const = 0;
 
     //! Divide mf by the diagonal component of the operator. Used by bicgstab.
     virtual void normalize (int amrlev, int mglev, MF& mf) const {
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLMG.H b/Src/LinearSolvers/MLMG/AMReX_MLMG.H
@@ -1317,10 +1317,7 @@ MLMGT<MF>::mgVcycle (int amrlev, int mglev_top)
 
         setVal(cor[amrlev][mglev], RT(0.0));
         bool skip_fillboundary = true;
-        for (int i = 0; i < nu1; ++i) {
-            linop.smooth(amrlev, mglev, cor[amrlev][mglev], res[amrlev][mglev], skip_fillboundary);
-            skip_fillboundary = false;
-        }
+        linop.smooth(amrlev, mglev, cor[amrlev][mglev], res[amrlev][mglev], skip_fillboundary, nu1);
 
         // rescor = res - L(cor)
         computeResOfCorrection(amrlev, mglev);
@@ -1364,11 +1361,8 @@ MLMGT<MF>::mgVcycle (int amrlev, int mglev_top)
         }
         setVal(cor[amrlev][mglev_bottom], RT(0.0));
         bool skip_fillboundary = true;
-        for (int i = 0; i < nu1; ++i) {
-            linop.smooth(amrlev, mglev_bottom, cor[amrlev][mglev_bottom],
-                         res[amrlev][mglev_bottom], skip_fillboundary);
-            skip_fillboundary = false;
-        }
+        linop.smooth(amrlev, mglev_bottom, cor[amrlev][mglev_bottom],
+            res[amrlev][mglev_bottom], skip_fillboundary, nu1);
         if (verbose >= 4)
         {
             computeResOfCorrection(amrlev, mglev_bottom);
@@ -1391,9 +1385,7 @@ MLMGT<MF>::mgVcycle (int amrlev, int mglev_top)
             amrex::Print() << print_ident << "AT LEVEL "  << amrlev << " " << mglev
                            << "   UP: Norm before smooth " << norm << "\n";
         }
-        for (int i = 0; i < nu2; ++i) {
-            linop.smooth(amrlev, mglev, cor[amrlev][mglev], res[amrlev][mglev]);
-        }
+        linop.smooth(amrlev, mglev, cor[amrlev][mglev], res[amrlev][mglev], false, nu2);
 
         if (cf_strategy == CFStrategy::ghostnodes) { computeResOfCorrection(amrlev, mglev); }
 
@@ -1512,10 +1504,7 @@ MLMGT<MF>::actualBottomSolve ()
     if (bottom_solver == BottomSolver::smoother)
     {
         bool skip_fillboundary = true;
-        for (int i = 0; i < nuf; ++i) {
-            linop.smooth(amrlev, mglev, x, b, skip_fillboundary);
-            skip_fillboundary = false;
-        }
+        linop.smooth(amrlev, mglev, x, b, skip_fillboundary, nuf);
     }
     else
     {
@@ -1589,9 +1578,7 @@ MLMGT<MF>::actualBottomSolve ()
                 setVal(cor[amrlev][mglev], RT(0.0));
             }
             const int n = (ret==0) ? nub : nuf;
-            for (int i = 0; i < n; ++i) {
-                linop.smooth(amrlev, mglev, x, b);
-            }
+            linop.smooth(amrlev, mglev, x, b, false, n);
         }
     }
 
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLNodeLinOp.H b/Src/LinearSolvers/MLMG/AMReX_MLNodeLinOp.H
@@ -45,7 +45,7 @@ public:
                         StateMode s_mode, const MLMGBndry* bndry=nullptr) const final;
 
     void smooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs,
-                         bool skip_fillboundary=false) const override;
+                 bool skip_fillboundary, int niter) const override;
 
     void solutionResidual (int amrlev, MultiFab& resid, MultiFab& x, const MultiFab& b,
                                    const MultiFab* crse_bcdata=nullptr) override;
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLNodeLinOp.cpp b/Src/LinearSolvers/MLMG/AMReX_MLNodeLinOp.cpp
@@ -167,12 +167,15 @@ MLNodeLinOp::apply (int amrlev, int mglev, MultiFab& out, MultiFab& in, BCMode b
 
 void
 MLNodeLinOp::smooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs,
-                     bool skip_fillboundary) const
+                     bool skip_fillboundary, int niter) const
 {
-    if (!skip_fillboundary) {
-        applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Correction);
+    for (int i = 0; i < niter; ++i) {
+        if (!skip_fillboundary) {
+            applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Correction);
+        }
+        Fsmooth(amrlev, mglev, sol, rhs);
+        skip_fillboundary = false;
     }
-    Fsmooth(amrlev, mglev, sol, rhs);
 }
 
 Real
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLNodeTensorLaplacian.H b/Src/LinearSolvers/MLMG/AMReX_MLNodeTensorLaplacian.H
@@ -58,7 +58,7 @@ public:
                  MultiFab& fine_res, MultiFab& fine_sol, const MultiFab& fine_rhs) const final;
 
     void smooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs,
-                 bool skip_fillboundary=false) const final;
+                 bool skip_fillboundary, int niter) const final;
 
     void prepareForSolve () final;
     void Fapply (int amrlev, int mglev, MultiFab& out, const MultiFab& in) const final;
diff --git a/Src/LinearSolvers/MLMG/AMReX_MLNodeTensorLaplacian.cpp b/Src/LinearSolvers/MLMG/AMReX_MLNodeTensorLaplacian.cpp
@@ -225,18 +225,20 @@ MLNodeTensorLaplacian::Fapply (int amrlev, int mglev, MultiFab& out, const Multi
 
 void
 MLNodeTensorLaplacian::smooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs,
-                               bool skip_fillboundary) const
+                               bool skip_fillboundary, int niter) const
 {
     BL_PROFILE("MLNodeTensorLaplacian::smooth()");
-    for (int redblack = 0; redblack < 4; ++redblack) {
-        if (!skip_fillboundary) {
-            applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Correction);
+    for (int i = 0; i < niter; ++i) {
+        for (int redblack = 0; redblack < 4; ++redblack) {
+            if (!skip_fillboundary) {
+                applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Correction);
+            }
+            m_redblack = redblack;
+            Fsmooth(amrlev, mglev, sol, rhs);
+            skip_fillboundary = false;
         }
-        m_redblack = redblack;
-        Fsmooth(amrlev, mglev, sol, rhs);
-        skip_fillboundary = false;
+        nodalSync(amrlev, mglev, sol);
     }
-    nodalSync(amrlev, mglev, sol);
 }
 
 void

Original file line number	Diff line number	Diff line change
`@@ -167,12 +167,15 @@ MLNodeLinOp::apply (int amrlev, int mglev, MultiFab& out, MultiFab& in, BCMode b`
`167`	`167`
`168`	`168`	`void`
`169`	`169`	`MLNodeLinOp::smooth (int amrlev, int mglev, MultiFab& sol, const MultiFab& rhs,`
`170`		`- bool skip_fillboundary) const`
	`170`	`+ bool skip_fillboundary, int niter) const`
`171`	`171`	`{`
`172`		`- if (!skip_fillboundary) {`
`173`		`- applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Correction);`
	`172`	`+ for (int i = 0; i < niter; ++i) {`
	`173`	`+ if (!skip_fillboundary) {`
	`174`	`+ applyBC(amrlev, mglev, sol, BCMode::Homogeneous, StateMode::Correction);`
	`175`	`+ }`
	`176`	`+ Fsmooth(amrlev, mglev, sol, rhs);`
	`177`	`+ skip_fillboundary = false;`
`174`	`178`	`}`
`175`		`- Fsmooth(amrlev, mglev, sol, rhs);`
`176`	`179`	`}`
`177`	`180`
`178`	`181`	`Real`