Merge branch 'development' into updated_scaling

AMReX-Astro · Jul 10, 2024 · 42c9441 · 42c9441
2 parents 7e17909 + 5590794
commit 42c9441
Show file tree

Hide file tree

Showing 24 changed files with 170 additions and 349 deletions.
diff --git a/Exec/mhd_tests/LoopAdvection/problem_initialize.H b/Exec/mhd_tests/LoopAdvection/problem_initialize.H
@@ -16,8 +16,8 @@ void problem_initialize ()
       eos_state.rho = problem::rho_0;
       eos_state.p = problem::p_0;
       eos_state.T = 100000.0_rt;  // initial guess
-      for (int n = 0; n < NumSpec; n++) {
-          eos_state.xn[n] = 0.0_rt;
+      for (auto & X : eos_state.xn) {
+          X = 0.0_rt;
       }
       eos_state.xn[0] = 1.0_rt;
 

diff --git a/Exec/mhd_tests/LoopAdvection/problem_initialize_mhd_data.H b/Exec/mhd_tests/LoopAdvection/problem_initialize_mhd_data.H
@@ -7,6 +7,8 @@ Real A_z (int i, int j, int k,
           const GeometryData& geomdata)
 {
 
+    amrex::ignore_unused(k);
+
     // Compute A_z.  This lives on edges, e.g., {A_z}_{i-1/2,j-1/2,k}
     // so it is centered only in the z-direction.
 

diff --git a/Exec/mhd_tests/LoopAdvection/problem_initialize_state_data.H b/Exec/mhd_tests/LoopAdvection/problem_initialize_state_data.H
@@ -10,6 +10,8 @@ void problem_initialize_state_data (int i, int j, int k,
                                     const GeometryData& geomdata)
 {
 
+    amrex::ignore_unused(geomdata);
+
     state(i,j,k,URHO) = problem::rho_0;
     state(i,j,k,UMX) = problem::rho_0 * problem::u_x;
     state(i,j,k,UMY) = problem::rho_0 * problem::u_y;

diff --git a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_simp_sdc b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_simp_sdc
@@ -1,5 +1,5 @@
 ## Latest inputs file being used to reproduce initial conditions from Pakmor et al. 2022
-## with 12 km resolution using simplified sdc
+## with 12.5 km resolution using simplified sdc
 
 ############################## CASTRO INPUTS ###############################################
 

diff --git a/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_strang b/Exec/science/wdmerger/tests/he_double_det/inputs_pakmor_strang
@@ -1,5 +1,5 @@
 ## Latest inputs file being used to reproduce initial conditions from Pakmor et al. 2022
-## with 12 km resolution using default strang-splitting
+## with 12.5 km resolution using default strang-splitting
 
 ############################## CASTRO INPUTS ###############################################
 

diff --git a/Source/hydro/Castro_hydro.H b/Source/hydro/Castro_hydro.H
@@ -118,19 +118,6 @@
                  amrex::Array4<amrex::Real> const& q_arr,
                  amrex::Array4<amrex::Real> const& qaux_arr);
 
-///
-/// compute the flattening coefficient.  This is 0 if we are in a shock and
-/// 1 if we are not applying any flattening
-///
-/// @param bx         the box to operate over
-/// @param q_arr      the primitive variable state
-/// @param flatn      the flattening coefficient
-/// @param pres_comp  index into q_arr of the pressure component
-///
-    static void uflatten(const amrex::Box& bx,
-                  amrex::Array4<amrex::Real const> const& q_arr,
-                  amrex::Array4<amrex::Real> const& flatn,
-                  const int pres_comp);
 
 ///
 /// A multidimensional shock detection algorithm

diff --git a/Source/hydro/Castro_mol_hydro.cpp b/Source/hydro/Castro_mol_hydro.cpp
@@ -8,6 +8,7 @@
 #include <advection_util.H>
 
 #include <fourth_center_average.H>
+#include <flatten.H>
 
 using namespace amrex;
 
@@ -117,19 +118,23 @@ Castro::construct_mol_hydro_source(Real time, Real dt, MultiFab& A_update)
         Array4<Real> const flatn_arr = flatn.array();
 
         if (first_order_hydro == 1) {
-          amrex::ParallelFor(obx,
-          [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
-          {
-            flatn_arr(i,j,k) = 0.0;
-          });
+            amrex::ParallelFor(obx,
+            [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
+            {
+                flatn_arr(i,j,k) = 0.0;
+            });
         } else if (use_flattening == 1) {
-          uflatten(obx, q_arr, flatn_arr, QPRES);
+            amrex::ParallelFor(obx,
+            [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
+            {
+                flatn_arr(i,j,k) = hydro::flatten(i, j, k, q_arr, QPRES);
+            });
         } else {
-          amrex::ParallelFor(obx,
-          [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
-          {
-            flatn_arr(i,j,k) = 1.0;
-          });
+            amrex::ParallelFor(obx,
+            [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
+            {
+                flatn_arr(i,j,k) = 1.0;
+            });
         }
 
 

diff --git a/Source/hydro/Make.package b/Source/hydro/Make.package
@@ -10,7 +10,6 @@ CEXE_sources += Castro_mol.cpp
 CEXE_headers += advection_util.H
 CEXE_sources += advection_util.cpp
 CEXE_headers += flatten.H
-CEXE_sources += flatten.cpp
 
 ifeq ($(USE_TRUE_SDC),TRUE)
   CEXE_sources += Castro_mol_hydro.cpp

diff --git a/Source/hydro/flatten.cpp b/Source/hydro/flatten.cpp
diff --git a/Source/hydro/riemann_constants.H b/Source/hydro/riemann_constants.H
@@ -9,6 +9,9 @@ namespace riemann_constants {
     constexpr amrex::Real smlp1 = 1.e-10_rt;
     constexpr amrex::Real small = 1.e-8_rt;
     constexpr amrex::Real smallu = 1.e-12_rt;
+    constexpr amrex::Real riemann_integral_tol = 1.e-8_rt;
+    constexpr amrex::Real riemann_u_tol = 1.e-6_rt;
+    constexpr amrex::Real riemann_p_tol = 1.e-8_rt;
     constexpr int HISTORY_SIZE=40;
     constexpr int PSTAR_BISECT_FACTOR = 5;
 }

diff --git a/Source/mhd/Castro_mhd.H b/Source/mhd/Castro_mhd.H
@@ -99,7 +99,7 @@
               amrex::Array4<amrex::Real const> const& Ed2,
               const int d, const int d1, const int d2, const amrex::Real dt);
 
-    void
+    static void
     hlld(const amrex::Box& bx,
          amrex::Array4<amrex::Real const> const& qleft,
          amrex::Array4<amrex::Real const> const& qright,

diff --git a/Source/mhd/Castro_mhd.cpp b/Source/mhd/Castro_mhd.cpp
@@ -1,6 +1,7 @@
 #include <Castro.H>
 
 #include <advection_util.H>
+#include <flatten.H>
 
 using namespace amrex;
 
@@ -196,10 +197,6 @@ Castro::construct_ctu_mhd_source(Real time, Real dt)
           auto flatn_arr = flatn.array();
           auto elix_flatn = flatn.elixir();
 
-          flatg.resize(bxi, 1);
-          auto flatg_arr = flatg.array();
-          auto elix_flatg = flatg.elixir();
-
           if (use_flattening == 0) {
             amrex::ParallelFor(bxi,
             [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
@@ -209,13 +206,11 @@ Castro::construct_ctu_mhd_source(Real time, Real dt)
 
           } else {
 
-            uflatten(bxi, q_arr, flatn_arr, QPRES);
-            uflatten(bxi, q_arr, flatg_arr, QPTOT);
-
             amrex::ParallelFor(bxi,
             [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
             {
-              flatn_arr(i,j,k) = flatn_arr(i,j,k) * flatg_arr(i,j,k);
+                flatn_arr(i,j,k) = hydro::flatten(i, j, k, q_arr, QPRES) *
+                                   hydro::flatten(i, j, k, q_arr, QPTOT);
             });
 
           }