EB Updates: (1) TerrainIF improvement, (2) SRD for RhoTheta_comp in erf_slow_rhs_post.

CMakeLists.txt

@@ -238,7 +238,7 @@ add_library(${PROJECT_NAME}::erf_api ALIAS erf_srclib)
 
 # Collect all headers and make them installable with the target
 
-set(ERF_INCLUDES "Source/ERF.H;Source/ERF_Constants.H;Source/WindFarmParametrization/SimpleActuatorDisk/ERF_SimpleAD.H;Source/WindFarmParametrization/EWP/ERF_EWP.H;Source/WindFarmParametrization/Null/ERF_NullWindFarm.H;Source/WindFarmParametrization/ERF_WindFarm.H;Source/WindFarmParametrization/Fitch/ERF_Fitch.H;Source/BoundaryConditions/ERF_PhysBCFunct.H;Source/BoundaryConditions/ERF_MOSTAverage.H;Source/BoundaryConditions/ERF_MOSTRoughness.H;Source/BoundaryConditions/ERF_ABLMost.H;Source/BoundaryConditions/ERF_FillPatcher.H;Source/BoundaryConditions/ERF_MOSTStress.H;Source/BoundaryConditions/ERF_TimeInterpolatedData.H;Source/Utils/ERF_Interpolation.H;Source/Utils/ERF_TileNoZ.H;Source/Utils/ERF_PlaneAverage.H;Source/Utils/ERF_Interpolation_WENO.H;Source/Utils/ERF_DirectionSelector.H;Source/Utils/ERF_ParFunctions.H;Source/Utils/ERF_Wstar.H;Source/Utils/ERF_MicrophysicsUtils.H;Source/Utils/ERF_SatMethods.H;Source/Utils/ERF_Interpolation_1D.H;Source/Utils/ERF_Interpolation_UPW.H;Source/Utils/ERF_TerrainMetrics.H;Source/Utils/ERF_Interpolation_WENO_Z.H;Source/Utils/ERF_Thetav.H;Source/Utils/ERF_WaterVaporSaturation.H;Source/Utils/ERF_Utils.H;Source/Utils/ERF_Orbit.H;Source/Utils/ERF_EOS.H;Source/Utils/ERF_HSEUtils.H;Source/Utils/ERF_Container.H;Source/EB/ERF_EB.H;Source/EB/ERF_EBAux.H;Source/EB/ERF_EBCutCell.H;Source/EB/ERF_EBIF.H;Source/EB/ERF_EBPolygon.H;Source/EB/ERF_EBUtils.H;Source/Particles/ERFPC.H;Source/Particles/ERF_ParticleData.H;Source/Prob/ERF_InitDensityHSEDry.H;Source/Prob/ERF_InitRayleighDamping.H;Source/Prob/ERF_InitConstantDensityHSE.H;Source/ERF_ProbCommon.H;Source/ERF_Derive.H;Source/Radiation/ERF_RRTMGP_Interface.H;Source/Radiation/ERF_RRTMGP_Utils.H;Source/Radiation/ERF_Radiation.H;Source/Radiation/ERF_OrbCosZenith.H;Source/SourceTerms/ERF_NumericalDiffusion.H;Source/SourceTerms/ERF_SrcHeaders.H;Source/SourceTerms/ERF_ForestDrag.H;Source/IO/ERF_SampleData.H;Source/IO/ERF_NCInterface.H;Source/IO/ERF_NCWpsFile.H;Source/IO/ERF_NCPlotFile.H;Source/IO/ERF_ReadBndryPlanes.H;Source/IO/ERF_WriteBndryPlanes.H;Source/PBL/ERF_MYNNStruct.H;Source/PBL/ERF_PBLModels.H;Source/PBL/ERF_PBLHeight.H;Source/TimeIntegration/ERF_TI_substep_fun.H;Source/TimeIntegration/ERF_TI_slow_headers.H;Source/TimeIntegration/ERF_TI_slow_rhs_fun.H;Source/TimeIntegration/ERF_TI_fast_headers.H;Source/TimeIntegration/ERF_TI_utils.H;Source/TimeIntegration/ERF_MRI.H;Source/TimeIntegration/ERF_TI_no_substep_fun.H;Source/LandSurfaceModel/Null/ERF_NullSurf.H;Source/LandSurfaceModel/ERF_LandSurface.H;Source/LandSurfaceModel/MM5/ERF_MM5.H;Source/LandSurfaceModel/SLM/ERF_SLM.H;Source/ERF_IndexDefines.H;Source/Advection/ERF_AdvectionSrcForMom_N.H;Source/Advection/ERF_AdvectionSrcForScalars.H;Source/Advection/ERF_AdvectionSrcForMom_T.H;Source/Advection/ERF_Advection.H;Source/Initialization/ERF_MetgridUtils.H;Source/Diffusion/ERF_EddyViscosity.H;Source/Diffusion/ERF_Diffusion.H;Source/Microphysics/Null/ERF_NullMoistLagrangian.H;Source/Microphysics/Null/ERF_NullMoist.H;Source/Microphysics/ERF_Microphysics.H;Source/Microphysics/ERF_LagrangianMicrophysics.H;Source/Microphysics/ERF_EulerianMicrophysics.H;Source/Microphysics/Kessler/ERF_Kessler.H;Source/Microphysics/Morrison/ERF_Morrison.H;Source/Microphysics/SAM/ERF_SAM.H;Source/Microphysics/SatAdj/ERF_SatAdj.H;Source/DataStructs/ERF_InputSpongeData.H;Source/DataStructs/ERF_TurbPertStruct.H;Source/DataStructs/ERF_SpongeStruct.H;Source/DataStructs/ERF_AdvStruct.H;Source/DataStructs/ERF_DataStruct.H;Source/DataStructs/ERF_InputSoundingData.H;Source/DataStructs/ERF_DiffStruct.H;Source/DataStructs/ERF_TurbStruct.H;Source/LinearSolvers/ERF_TerrainPoisson.H;Source/LinearSolvers/ERF_FFTUtils.H")
+set(ERF_INCLUDES "Source/ERF.H;Source/ERF_Constants.H;Source/WindFarmParametrization/SimpleActuatorDisk/ERF_SimpleAD.H;Source/WindFarmParametrization/EWP/ERF_EWP.H;Source/WindFarmParametrization/Null/ERF_NullWindFarm.H;Source/WindFarmParametrization/ERF_WindFarm.H;Source/WindFarmParametrization/Fitch/ERF_Fitch.H;Source/BoundaryConditions/ERF_PhysBCFunct.H;Source/BoundaryConditions/ERF_MOSTAverage.H;Source/BoundaryConditions/ERF_MOSTRoughness.H;Source/BoundaryConditions/ERF_ABLMost.H;Source/BoundaryConditions/ERF_FillPatcher.H;Source/BoundaryConditions/ERF_MOSTStress.H;Source/BoundaryConditions/ERF_TimeInterpolatedData.H;Source/Utils/ERF_Interpolation.H;Source/Utils/ERF_TileNoZ.H;Source/Utils/ERF_PlaneAverage.H;Source/Utils/ERF_Interpolation_WENO.H;Source/Utils/ERF_DirectionSelector.H;Source/Utils/ERF_ParFunctions.H;Source/Utils/ERF_Wstar.H;Source/Utils/ERF_MicrophysicsUtils.H;Source/Utils/ERF_SatMethods.H;Source/Utils/ERF_Interpolation_1D.H;Source/Utils/ERF_Interpolation_UPW.H;Source/Utils/ERF_TerrainMetrics.H;Source/Utils/ERF_Interpolation_WENO_Z.H;Source/Utils/ERF_Thetav.H;Source/Utils/ERF_WaterVaporSaturation.H;Source/Utils/ERF_Utils.H;Source/Utils/ERF_Orbit.H;Source/Utils/ERF_EOS.H;Source/Utils/ERF_HSEUtils.H;Source/Utils/ERF_Container.H;Source/EB/ERF_EB.H;Source/EB/ERF_EBAdvection.H;Source/EB/ERF_EBAdvectionSrcForScalars.H;Source/EB/ERF_EBAux.H;Source/EB/ERF_EBCutCell.H;Source/EB/ERF_EBIF.H;Source/EB/ERF_EBIFTerrain.H;Source/EB/ERF_EBPolygon.H;Source/EB/ERF_EBRedistribute.H;Source/EB/ERF_EBToPVD.H;Source/EB/ERF_EBUtils.H;Source/Particles/ERFPC.H;Source/Particles/ERF_ParticleData.H;Source/Prob/ERF_InitDensityHSEDry.H;Source/Prob/ERF_InitRayleighDamping.H;Source/Prob/ERF_InitConstantDensityHSE.H;Source/ERF_ProbCommon.H;Source/ERF_Derive.H;Source/Radiation/ERF_RRTMGP_Interface.H;Source/Radiation/ERF_RRTMGP_Utils.H;Source/Radiation/ERF_Radiation.H;Source/Radiation/ERF_OrbCosZenith.H;Source/SourceTerms/ERF_NumericalDiffusion.H;Source/SourceTerms/ERF_SrcHeaders.H;Source/SourceTerms/ERF_ForestDrag.H;Source/IO/ERF_SampleData.H;Source/IO/ERF_NCInterface.H;Source/IO/ERF_NCWpsFile.H;Source/IO/ERF_NCPlotFile.H;Source/IO/ERF_ReadBndryPlanes.H;Source/IO/ERF_WriteBndryPlanes.H;Source/PBL/ERF_MYNNStruct.H;Source/PBL/ERF_PBLModels.H;Source/PBL/ERF_PBLHeight.H;Source/TimeIntegration/ERF_TI_substep_fun.H;Source/TimeIntegration/ERF_TI_slow_headers.H;Source/TimeIntegration/ERF_TI_slow_rhs_fun.H;Source/TimeIntegration/ERF_TI_fast_headers.H;Source/TimeIntegration/ERF_TI_utils.H;Source/TimeIntegration/ERF_MRI.H;Source/TimeIntegration/ERF_TI_no_substep_fun.H;Source/LandSurfaceModel/Null/ERF_NullSurf.H;Source/LandSurfaceModel/ERF_LandSurface.H;Source/LandSurfaceModel/MM5/ERF_MM5.H;Source/LandSurfaceModel/SLM/ERF_SLM.H;Source/ERF_IndexDefines.H;Source/Advection/ERF_AdvectionSrcForMom_N.H;Source/Advection/ERF_AdvectionSrcForScalars.H;Source/Advection/ERF_AdvectionSrcForMom_T.H;Source/Advection/ERF_Advection.H;Source/Initialization/ERF_MetgridUtils.H;Source/Diffusion/ERF_EddyViscosity.H;Source/Diffusion/ERF_Diffusion.H;Source/Microphysics/Null/ERF_NullMoistLagrangian.H;Source/Microphysics/Null/ERF_NullMoist.H;Source/Microphysics/ERF_Microphysics.H;Source/Microphysics/ERF_LagrangianMicrophysics.H;Source/Microphysics/ERF_EulerianMicrophysics.H;Source/Microphysics/Kessler/ERF_Kessler.H;Source/Microphysics/Morrison/ERF_Morrison.H;Source/Microphysics/SAM/ERF_SAM.H;Source/Microphysics/SatAdj/ERF_SatAdj.H;Source/DataStructs/ERF_InputSpongeData.H;Source/DataStructs/ERF_TurbPertStruct.H;Source/DataStructs/ERF_SpongeStruct.H;Source/DataStructs/ERF_AdvStruct.H;Source/DataStructs/ERF_DataStruct.H;Source/DataStructs/ERF_InputSoundingData.H;Source/DataStructs/ERF_DiffStruct.H;Source/DataStructs/ERF_TurbStruct.H;Source/LinearSolvers/ERF_TerrainPoisson.H;Source/LinearSolvers/ERF_FFTUtils.H")
 
 
 set_target_properties(

Exec/DevTests/EB_Test/inputs

@@ -66,7 +66,8 @@ erf.alpha_C = 0.0
 
 erf.init_type = "uniform"
 
-erf.terrain_type = EB
+erf.terrain_type = "EB"
+erf.write_eb_surface = TRUE
 
 # PROBLEM PARAMETERS
 prob.rho_0 = 1.0

Source/EB/ERF_EBIFTerrain.H

@@ -23,45 +23,91 @@ public:
           dz(a_geom.CellSize(2)),
           prob_lo_x(a_geom.ProbLo(0)),
           prob_lo_y(a_geom.ProbLo(1)),
-          prob_lo_z(a_geom.ProbLo(2))
+          prob_lo_z(a_geom.ProbLo(2)),
+          prob_hi_x(a_geom.ProbHi(0)),
+          prob_hi_y(a_geom.ProbHi(1)),
+          i_hi(static_cast<int>(std::round((a_geom.ProbHi(0)-a_geom.ProbLo(0))/a_geom.CellSize(0)))),
+          j_hi(static_cast<int>(std::round((a_geom.ProbHi(1)-a_geom.ProbLo(1))/a_geom.CellSize(1))))
         {}
 
     AMREX_GPU_HOST_DEVICE inline
     amrex::Real operator() (AMREX_D_DECL(amrex::Real x, amrex::Real y, amrex::Real z))
         const noexcept
     {
-        const int i1 = static_cast<int>(std::floor((x-prob_lo_x) / dx));
-        const int j1 = static_cast<int>(std::floor((y-prob_lo_y) / dy));
+        int i1{};
+        int j1{};
+        int i2{};
+        int j2{};
+
+        amrex::Real x1{};
+        amrex::Real x2{};
+        amrex::Real y1{};
+        amrex::Real y2{};
+        amrex::Real w1{};
+        amrex::Real w2{};
+        amrex::Real terr_z{0.};
+
+        // z_stretched
         const int k1 = static_cast<int>(std::floor((z-prob_lo_z) / dz));
-
-        const int i2 = i1+1;
-        const int j2 = j1+1;
-
-        amrex::Real x1 = prob_lo_x + i1*dx;
-        amrex::Real y1 = prob_lo_y + j1*dy;
-        amrex::Real z1 = prob_lo_z + k1*dz;
-
-        amrex::Real x2 = prob_lo_x + i2*dx;
-        amrex::Real y2 = prob_lo_y + j2*dy;
-
-        amrex::Real remainder_z = (z - z1)/dz;
-
+        const amrex::Real z1 = prob_lo_z + k1*dz;
+        const amrex::Real remainder_z = (z - z1)/dz;
         amrex::Real z_stretched = prob_lo_z;
-
         for (int kk = 0; kk < k1; ++kk) {
             z_stretched += dz_s[kk];
         }
         z_stretched += remainder_z * dz_s[k1];
 
-        // Do bilinear interpolation of the terrain surface
-        amrex::Real denom = (x2-x1)*(y2-y1);
-        amrex::Real  w_11 = (x2-x)*(y2-y)/denom;
-        amrex::Real  w_12 = (x2-x)*(y-y1)/denom;
-        amrex::Real  w_21 = (x-x1)*(y2-y)/denom;
-        amrex::Real  w_22 = (x-x1)*(y-y1)/denom;
-
-        amrex::Real terr_z = w_11*terr_arr(i1,j1,0) + w_12*terr_arr(i1,j2,0) + w_21*terr_arr(i2,j1,0) + w_22*terr_arr(i2,j2,0);
-
+        // Interpolation (nine subregions of x-y plane)
+        if        (x <= prob_lo_x && y <= prob_lo_y) {
+            terr_z = terr_arr(0   ,0   ,0);
+        } else if (x >= prob_hi_x && y <= prob_lo_y) {
+            terr_z = terr_arr(i_hi,0   ,0);
+        } else if (x >= prob_hi_x && y >= prob_hi_y) {
+            terr_z = terr_arr(i_hi,j_hi,0);
+        } else if (x <= prob_lo_x && y >= prob_hi_y) {
+            terr_z = terr_arr(0   ,j_hi,0);
+        } else if (x > prob_lo_x && x < prob_hi_x && (y <= prob_lo_y || y >= prob_hi_y) ) {
+            i1 = static_cast<int>(std::floor((x-prob_lo_x) / dx));
+            i2 = i1+1;
+            x1 = prob_lo_x + i1*dx;
+            x2 = x1 + dx;
+            w1 = (x2-x)/dx;
+            w2 = (x-x1)/dx;
+            if        (y <= prob_lo_y){
+                terr_z = w1*terr_arr(i1,0   ,0) + w2*terr_arr(i2,0   ,0);
+            } else if (y >= prob_hi_y) {
+                terr_z = w1*terr_arr(i1,j_hi,0) + w2*terr_arr(i2,j_hi,0);
+            }
+        } else if (y > prob_lo_y && y < prob_hi_y && (x <= prob_lo_x || x >= prob_hi_x) ) {
+            j1 = static_cast<int>(std::floor((y-prob_lo_y) / dy));
+            j2 = j1+1;
+            y1 = prob_lo_y + j1*dy;
+            y2 = y1 + dy;
+            w1 = (y2-y)/dy;
+            w2 = (y-y1)/dy;
+            if        (x <= prob_lo_x){
+                terr_z = w1*terr_arr(0   ,j1,0) + w2*terr_arr(0   ,j2,0);
+            } else if (x >= prob_hi_x) {
+                terr_z = w1*terr_arr(i_hi,j1,0) + w2*terr_arr(i_hi,j2,0);
+            }
+        } else {
+            // Do bilinear interpolation of the terrain surface
+            i1 = static_cast<int>(std::floor((x-prob_lo_x) / dx));
+            i2 = i1+1;
+            j1 = static_cast<int>(std::floor((y-prob_lo_y) / dy));
+            j2 = j1+1;
+            x1 = prob_lo_x + i1*dx;
+            x2 = x1 + dx;
+            y1 = prob_lo_y + j1*dy;
+            y2 = y1 + dy;
+
+            const amrex::Real denom = dx*dy;
+            const amrex::Real w11 = (x2-x)*(y2-y)/denom;
+            const amrex::Real w12 = (x2-x)*(y-y1)/denom;
+            const amrex::Real w21 = (x-x1)*(y2-y)/denom;
+            const amrex::Real w22 = (x-x1)*(y-y1)/denom;
+            terr_z = w11*terr_arr(i1,j1,0) + w12*terr_arr(i1,j2,0) + w21*terr_arr(i2,j1,0) + w22*terr_arr(i2,j2,0);
+        }
         return -(z_stretched - terr_z);
     }
 
@@ -76,6 +122,8 @@ protected:
     amrex::Real const* dz_s;
     amrex::Real dx, dy, dz;
     amrex::Real prob_lo_x, prob_lo_y, prob_lo_z;
+    amrex::Real prob_hi_x, prob_hi_y;
+    int i_hi, j_hi;
 };
 
-#endif
+#endif

Source/EB/ERF_EBRedistribute.H

@@ -0,0 +1,29 @@
+#ifndef ERF_EBRedistribute_H_
+#define ERF_EBRedistribute_H_
+
+#include <AMReX.H>
+#include <AMReX_MultiFab.H>
+#include <AMReX_Geometry.H>
+#include <AMReX_BCRec.H>
+#include <AMReX_MFIter.H>
+
+void redistribute_term (int ncomp,
+                        const amrex::Geometry& geom,
+                        amrex::MultiFab& result,
+                        amrex::MultiFab& result_tmp,
+                        amrex::MultiFab const& state,
+                        amrex::EBFArrayBoxFactory const& ebfact,
+                        amrex::BCRec const* bc,
+                        amrex::Real const dt);
+
+void redistribute_term (amrex::MFIter const& mfi,
+                        int ncomp,
+                        const amrex::Geometry& geom,
+                        amrex::MultiFab& result,
+                        amrex::MultiFab& result_tmp,
+                        amrex::MultiFab const& state,
+                        amrex::EBFArrayBoxFactory const& ebfact,
+                        amrex::BCRec const* bc,
+                        amrex::Real const dt);
+
+#endif

Source/EB/ERF_EBRedistribute.cpp

@@ -1,55 +1,62 @@
 #include <AMReX_Config.H>
+#include <AMReX_Geometry.H>
 
 #include <ERF.H>
+#include <ERF_EBRedistribute.H>
 #include "AMReX_EB_Redistribution.H"
 
 using namespace amrex;
 
 void
-ERF::redistribute_term (int lev, int ncomp,
-                        MultiFab& result,
-                        MultiFab& result_tmp,
-                        MultiFab const& state,
-                        BCRec const* bc, // this is bc for the state (needed for SRD slopes)
-                        Real const local_dt)
+redistribute_term ( int ncomp,
+                    const Geometry& geom,
+                    MultiFab& result,
+                    MultiFab& result_tmp,
+                    MultiFab const& state,
+                    EBFArrayBoxFactory const& ebfact,
+                    BCRec const* bc, // this is bc for the state (needed for SRD slopes)
+                    Real const local_dt)
 {
     // ************************************************************************
     // Redistribute result_tmp and pass out result
     // ************************************************************************
     AMREX_ASSERT(result.nComp() == state.nComp());
 
-    result_tmp.FillBoundary(geom[lev].periodicity());
+    result_tmp.FillBoundary(geom.periodicity());
 
 #ifdef _OPENMP
 #pragma omp parallel if (Gpu::notInLaunchRegion())
 #endif
     for (MFIter mfi(state,TilingIfNotGPU()); mfi.isValid(); ++mfi)
     {
-        redistribute_term(mfi, lev, ncomp, result, result_tmp, state, bc, local_dt);
+        redistribute_term(mfi, ncomp, geom, result, result_tmp, state, ebfact, bc, local_dt);
     }
 }
 
 void
-ERF::redistribute_term (MFIter const& mfi, int lev, int ncomp,
-                        MultiFab& result,
-                        MultiFab& result_tmp,
-                        MultiFab const& state,
-                        BCRec const* bc, // this is bc for the state (needed for SRD slopes)
-                        Real const local_dt)
+redistribute_term ( MFIter const& mfi,
+                    int ncomp,
+                    const Geometry& geom,
+                    MultiFab& result,
+                    MultiFab& result_tmp,
+                    MultiFab const& state,
+                    EBFArrayBoxFactory const& ebfact,
+                    BCRec const* bc, // this is bc for the state (needed for SRD slopes)
+                    Real const local_dt)
 {
     AMREX_ASSERT(result.nComp() == state.nComp());
 
     Box const& bx = mfi.tilebox();
 
-    EBFArrayBoxFactory const& ebfact = EBFactory(lev);
+    // EBFArrayBoxFactory const& ebfact = EBFactory(lev);
     EBCellFlagFab const& flagfab = ebfact.getMultiEBCellFlagFab()[mfi];
     Array4<EBCellFlag const> const& flag = flagfab.const_array();
 
     bool regular = (flagfab.getType(amrex::grow(bx,4)) == FabType::regular);
     bool covered = (flagfab.getType(bx) == FabType::covered);
 
-    Array4<Real      > out = result.array(mfi);
-    Array4<Real      > in  = result_tmp.array(mfi);
+    Array4<Real> out = result.array(mfi);
+    Array4<Real> in  = result_tmp.array(mfi);
 
     if (!regular && !covered)
     {
@@ -84,7 +91,7 @@ ERF::redistribute_term (MFIter const& mfi, int lev, int ncomp,
                              scratch, flag,
                              apx, apy, apz, vfrac,
                              fcx, fcy, fcz, ccc,
-                             bc, geom[lev], local_dt, redistribution_type);
+                             bc, geom, local_dt, redistribution_type);
     }
     else
     {

Source/EB/Make.package

@@ -18,6 +18,7 @@ CEXE_headers += ERF_EBAdvection.H
 CEXE_headers += ERF_EBAdvectionSrcForScalars.H
 
 CEXE_sources += ERF_EBBox.cpp
+CEXE_headers += ERF_EBRedistribute.H
 CEXE_sources += ERF_EBRedistribute.cpp
 CEXE_sources += ERF_EBWriteSurface.cpp
 

Source/ERF.H

@@ -477,18 +477,6 @@ public:
     void MakeEBGeometry ();
     void make_eb_box ();
     void make_eb_regular ();
-    void redistribute_term (int lev, int ncomp,
-                            amrex::MultiFab& result,
-                            amrex::MultiFab& result_tmp,
-                            amrex::MultiFab const& state,
-                            amrex::BCRec const* bc,
-                            amrex::Real const dt);
-    void redistribute_term (amrex::MFIter const& mfi, int lev, int ncomp,
-                            amrex::MultiFab& result,
-                            amrex::MultiFab& result_tmp,
-                            amrex::MultiFab const& state,
-                            amrex::BCRec const* bc,
-                            amrex::Real const dt);
 
     // more flexible version of AverageDown() that lets you average down across multiple levels
     void AverageDownTo (int crse_lev, int scomp, int ncomp); // NOLINT

Source/TimeIntegration/ERF_AdvanceDycore.cpp

@@ -9,6 +9,7 @@
 #include <ERF_Diffusion.H>
 #include <ERF_TileNoZ.H>
 #include <ERF_Utils.H>
+#include <ERF_EBRedistribute.H>
 
 using namespace amrex;