Generalize coarsening for MR (#945)

* Move interpolation functions for MR to new folder

* Preparatory clean-up of old namespace Average for future MR functions

* Add interpolation for MR in new namespace Coarsen

* Change file names Average.H/.cpp to Coarsen.H/.cpp

* Remove Source/Parallelization/WarpXComm.H (not necessary anymore)

* Coarsening for IO and MR in separate files

* Clean up IO and MR Coarsen namespaces

* Remove old interpolation functions (charge and current)

* Void commit: trigger Travis CI build

* Fix GPU build

* Void commit: trigger Travis CI build

* Add Python script to test interpolation points/weights in 1D

* Move using-directives inside namespaces

* Add Doxygen documentation and comments

* Minor style changes

* Improve new Python script

1 files changed, 9 insertions, 66 deletions

diff --git a/Source/Parallelization/WarpXComm.cpp b/Source/Parallelization/WarpXComm.cpp
index f55a9ef9a..e1aa6c316 100644
--- a/Source/Parallelization/WarpXComm.cpp
+++ b/Source/Parallelization/WarpXComm.cpp
@@ -6,12 +6,10 @@
  *
  * License: BSD-3-Clause-LBNL
  */
-#include "WarpXComm.H"
 #include "WarpXComm_K.H"
 #include "WarpX.H"
 #include "WarpXSumGuardCells.H"
-#include "InterpolateCurrentFineToCoarse.H"
-#include "InterpolateDensityFineToCoarse.H"
+#include "Utils/CoarsenMR.H"
 
 #include <algorithm>
 #include <cstdlib>
@@ -568,7 +566,9 @@ WarpX::SyncCurrent ()
         std::array<      MultiFab*,3> crse { current_cp[lev][0].get(),
                                              current_cp[lev][1].get(),
                                              current_cp[lev][2].get() };
-        interpolateCurrentFineToCoarse(fine, crse, refinement_ratio[0]);
+        CoarsenMR::Coarsen( *crse[0], *fine[0], refinement_ratio );
+        CoarsenMR::Coarsen( *crse[1], *fine[1], refinement_ratio );
+        CoarsenMR::Coarsen( *crse[2], *fine[2], refinement_ratio );
     }
 
     // For each level
@@ -581,40 +581,6 @@ WarpX::SyncCurrent ()
 }
 
 void
-interpolateCurrentFineToCoarse ( std::array< amrex::MultiFab const *, 3 > const & fine,
-                                 std::array< amrex::MultiFab       *, 3 > const & coarse,
-                                 int const refinement_ratio)
-{
-    WARPX_PROFILE("interpolateCurrentFineToCoarse()");
-    BL_ASSERT(refinement_ratio == 2);
-    const IntVect& ng = (fine[0]->nGrowVect() + 1) / refinement_ratio; // add equivalent no. of guards to coarse patch
-
-#ifdef _OPENMP
-#pragma omp parallel if (Gpu::notInLaunchRegion())
-#endif
-    {
-        for (int idim = 0; idim < fine.size(); ++idim)  // j-field components
-        {
-            // OMP in-box decomposition of coarse into tilebox
-            for (MFIter mfi(*coarse[idim], TilingIfNotGPU()); mfi.isValid(); ++mfi)
-            {
-                const Box& bx = mfi.growntilebox(ng); // only grow to outer directions of tileboxes for filling guards
-
-                auto const & arrFine = fine[idim]->const_array(mfi);
-                auto const & arrCoarse = coarse[idim]->array(mfi);
-
-                if( idim == 0 )
-                    amrex::ParallelFor( bx, InterpolateCurrentFineToCoarse<0>(arrFine, arrCoarse, refinement_ratio) );
-                else if( idim == 1 )
-                    amrex::ParallelFor( bx, InterpolateCurrentFineToCoarse<1>(arrFine, arrCoarse, refinement_ratio) );
-                else if( idim == 2 )
-                    amrex::ParallelFor( bx, InterpolateCurrentFineToCoarse<2>(arrFine, arrCoarse, refinement_ratio) );
-            }
-        }
-    }
-}
-
-void
 WarpX::SyncRho ()
 {
     WARPX_PROFILE("SyncRho()");
@@ -628,7 +594,7 @@ WarpX::SyncRho ()
     {
         rho_cp[lev]->setVal(0.0);
         const IntVect& refinement_ratio = refRatio(lev-1);
-        interpolateDensityFineToCoarse(*rho_fp[lev], *rho_cp[lev], refinement_ratio[0]);
+        CoarsenMR::Coarsen( *rho_cp[lev], *rho_fp[lev], refinement_ratio );
     }
 
     // For each level
@@ -640,31 +606,6 @@ WarpX::SyncRho ()
     }
 }
 
-void
-interpolateDensityFineToCoarse (const MultiFab& fine, MultiFab& coarse, int const refinement_ratio)
-{
-    WARPX_PROFILE("interpolateDensityFineToCoarse()");
-    BL_ASSERT(refinement_ratio == 2);
-    const IntVect& ng = (fine.nGrowVect() + 1) / refinement_ratio;  // add equivalent no. of guards to coarse patch
-    const int nc = fine.nComp();
-
-#ifdef _OPENMP
-#pragma omp parallel if (Gpu::notInLaunchRegion())
-#endif
-    {
-        // OMP in-box decomposition of coarse into tilebox
-        for (MFIter mfi(coarse, TilingIfNotGPU()); mfi.isValid(); ++mfi)
-        {
-            const Box& bx = mfi.growntilebox(ng); // only grow to outer directions of tileboxes for filling guards
-
-            amrex::ParallelFor(
-                bx,
-                InterpolateDensityFineToCoarse(fine.const_array(mfi), coarse.array(mfi), refinement_ratio, nc)
-            );
-        }
-    }
-}
-
 /** \brief Fills the values of the current on the coarse patch by
  *  averaging the values of the current of the fine patch (on the same level).
  */
@@ -683,7 +624,9 @@ WarpX::RestrictCurrentFromFineToCoarsePatch (int lev)
     std::array<      MultiFab*,3> crse { current_cp[lev][0].get(),
                                          current_cp[lev][1].get(),
                                          current_cp[lev][2].get() };
-    interpolateCurrentFineToCoarse(fine, crse, refinement_ratio[0]);
+    CoarsenMR::Coarsen( *crse[0], *fine[0], refinement_ratio );
+    CoarsenMR::Coarsen( *crse[1], *fine[1], refinement_ratio );
+    CoarsenMR::Coarsen( *crse[2], *fine[2], refinement_ratio );
 }
 
 void
@@ -792,7 +735,7 @@ WarpX::RestrictRhoFromFineToCoarsePatch (int lev)
     if (rho_fp[lev]) {
         rho_cp[lev]->setVal(0.0);
         const IntVect& refinement_ratio = refRatio(lev-1);
-        interpolateDensityFineToCoarse(*rho_fp[lev], *rho_cp[lev], refinement_ratio[0]);
+        CoarsenMR::Coarsen( *rho_cp[lev], *rho_fp[lev], refinement_ratio );
     }
 }