Reaction mask implementation gpu

Source/PeleC.H

@@ -742,6 +742,10 @@ protected:
 
   static amrex::Vector<int> src_list;
 
+  static bool use_chem_mask;
+  static amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> lo_chem_mask_coordinate;
+  static amrex::GpuArray<amrex::Real, AMREX_SPACEDIM> hi_chem_mask_coordinate;
+
   static bool use_typical_vals_chem;
   static bool use_typical_vals_chem_usr;
   static amrex::Real typical_rhoY_val_min;

Source/PeleC.cpp

@@ -108,6 +108,12 @@ amrex::Vector<std::string> PeleC::m_diagVars;
 
 amrex::Vector<int> PeleC::src_list;
 
+bool PeleC::use_chem_mask = false; // Flag to check if mask is activated.
+amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>
+  PeleC::lo_chem_mask_coordinate; // Box coordinate low
+amrex::GpuArray<amrex::Real, AMREX_SPACEDIM>
+  PeleC::hi_chem_mask_coordinate; // Box coordinate high
+
 // this will be reset upon restart
 amrex::Real PeleC::previousCPUTimeUsed = 0.0;
 amrex::Real PeleC::startCPUTime = 0.0;
@@ -250,6 +256,16 @@ PeleC::read_params()
   pp.query("use_typ_vals_chem", use_typical_vals_chem);
   pp.query("use_typ_vals_chem_usr", use_typical_vals_chem_usr);
 
+  pp.query("chem_mask", use_chem_mask);
+
+  // Reading chemistry mask box coordinates
+  if (use_chem_mask) {
+    for (int n = 0; n < AMREX_SPACEDIM; ++n) {
+      pp.get("lo_chemmask", lo_chem_mask_coordinate[n], n);
+      pp.get("hi_chemmask", hi_chem_mask_coordinate[n], n);
+    }
+  }
+
   if (use_typical_vals_chem_usr) {
     use_typical_vals_chem = true;
   }

Source/React.cpp

@@ -38,6 +38,9 @@ PeleC::react_state(
 
   const amrex::Real strt_time = amrex::ParallelDescriptor::second();
 
+  const auto plo = geom.ProbLoArray();
+  const auto dx = geom.CellSizeArray();
+
   AMREX_ASSERT(do_react == 1);
 
   if ((verbose != 0) && amrex::ParallelDescriptor::IOProcessor()) {
@@ -126,6 +129,15 @@ PeleC::react_state(
 #pragma omp parallel if (amrex::Gpu::notInLaunchRegion())
 #endif
   {
+
+    amrex::RealBox Chem_Masked_Region(
+      AMREX_D_DECL(
+        lo_chem_mask_coordinate[0], lo_chem_mask_coordinate[1],
+        lo_chem_mask_coordinate[2]),
+      AMREX_D_DECL(
+        hi_chem_mask_coordinate[0], hi_chem_mask_coordinate[1],
+        hi_chem_mask_coordinate[2]));
+
     for (amrex::MFIter mfi(S_new, amrex::TilingIfNotGPU()); mfi.isValid();
          ++mfi) {
 
@@ -170,6 +182,20 @@ PeleC::react_state(
         auto const& mask = dummyMask.array(mfi);
         auto const& fc = fctCount.array(mfi);
 
+        if (use_chem_mask) {
+          amrex::ParallelFor(
+            bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
+              amrex::XDim3 point;
+              point.x = plo[0] + (i + 0.5) * dx[0];
+              point.y = plo[1] + (j + 0.5) * dx[1];
+              point.z = plo[2] + (k + 0.5) * dx[2];
+
+              if (Chem_Masked_Region.contains(point)) {
+                mask(i, j, k) = -1;
+              }
+            });
+        }
+
         amrex::ParallelFor(
           bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
             // work on old state
@@ -209,7 +235,7 @@ PeleC::react_state(
         );
 
         amrex::Gpu::Device::streamSynchronize();
-
+        bool use_chem_mask_d = use_chem_mask;
         // unpack data
         amrex::ParallelFor(
           bx, [=] AMREX_GPU_DEVICE(int i, int j, int k) noexcept {
@@ -236,6 +262,13 @@ PeleC::react_state(
                - rho_old * e_old) // old internal energy
               / dt;
 
+            if (use_chem_mask_d && mask(i, j, k) == -1) {
+              for (int nsp = 0; nsp < NUM_SPECIES; nsp++) {
+                rhoY(i, j, k, nsp) += nonrs_arr(i, j, k, UFS + nsp) * dt;
+              }
+              rhoY(i, j, k, NUM_SPECIES) += rhoedot_ext * dt;
+            }
+
             amrex::Real umnew =
               sold_arr(i, j, k, UMX) + dt * nonrs_arr(i, j, k, UMX);
             amrex::Real vmnew =