Run on Frontier: update singularity, relocatable device code, misc fixes

src/geometry/mckinney_gammie_ryan.hpp

@@ -54,7 +54,7 @@ class McKinneyGammieRyan {
                      Real smooth, Real hexp_br, Real hexp_nsq, Real hexp_csq)
       : derefine_poles_(derefine_poles), h_(h), xt_(xt), alpha_(alpha), x0_(x0),
         smooth_(smooth), norm_(GetNorm_(alpha_, xt_)), hexp_br_(hexp_br),
-        hexp_nsq_(hexp_nsq), hexp_csq_(hexp_csq_) {}
+        hexp_nsq_(hexp_nsq), hexp_csq_(hexp_csq) {}
   McKinneyGammieRyan(bool derefine_poles, Real h, Real xt, Real alpha, Real x0,
                      Real smooth)
       : derefine_poles_(derefine_poles), h_(h), xt_(xt), alpha_(alpha), x0_(x0),

src/microphysics/eos_phoebus/eos_phoebus.cpp

@@ -87,7 +87,6 @@ std::shared_ptr<StateDescriptor> Initialize(ParameterInput *pin) {
     // Can specify rho_min, etc, in <eos>
     rho_min = pin->GetOrAddReal("eos", "rho_min", 0.0);
     sie_min = pin->GetOrAddReal("eos", "sie_min", 0.0);
-    lambda[2] = {0.};
     T_min = eos_host.TemperatureFromDensityInternalEnergy(rho_min, sie_min, lambda);
     rho_max = pin->GetOrAddReal("eos", "rho_max", 1e18);
     sie_max = pin->GetOrAddReal("eos", "sie_max", 1e35);

src/pgen/pgen.cpp

@@ -14,7 +14,6 @@
 #include "geometry/geometry.hpp"
 
 #include "pgen.hpp"
-#include "phoebus_utils/root_find.hpp"
 
 namespace phoebus {
 
@@ -56,32 +55,4 @@ void PostInitializationModifier(ParameterInput *pin, Mesh *pmesh) {
   f(pin, pmesh);
 }
 
-class PressResidual {
- public:
-  KOKKOS_INLINE_FUNCTION
-  PressResidual(const EOS &eos, const Real rho, const Real P, const Real Ye)
-      : eos_(eos), rho_(rho), P_(P) {
-    lambda_[0] = Ye;
-  }
-  KOKKOS_INLINE_FUNCTION
-  Real operator()(const Real e) {
-    return eos_.PressureFromDensityInternalEnergy(rho_, e, lambda_) - P_;
-  }
-
- private:
-  const EOS &eos_;
-  Real rho_, P_;
-  Real lambda_[2];
-};
-
-KOKKOS_FUNCTION
-Real energy_from_rho_P(const EOS &eos, const Real rho, const Real P, const Real emin,
-                       const Real emax, const Real Ye) {
-  PARTHENON_REQUIRE(P >= 0, "Pressure is negative!");
-  PressResidual res(eos, rho, P, Ye);
-  root_find::RootFind root;
-  Real eroot = root.regula_falsi(res, emin, emax, 1.e-10 * P, emin - 1.e10);
-  return rho * eroot;
-}
-
 } // namespace phoebus

src/pgen/pgen.hpp

@@ -22,6 +22,7 @@ using namespace parthenon::package::prelude;
 // internal includes
 #include "fluid/fluid.hpp"
 #include "geometry/geometry.hpp"
+#include "phoebus_utils/root_find.hpp"
 #include "phoebus_utils/variables.hpp"
 #include "radiation/radiation.hpp"
 #include "tracers/tracers.hpp"
@@ -117,9 +118,34 @@ static std::map<std::string, std::function<void(ParameterInput *pin, Mesh *pmesh
 // END OF UNTOUCHABLE MACRO SECTION
 */
 
-KOKKOS_FUNCTION
-Real energy_from_rho_P(const Microphysics::EOS::EOS &eos, const Real rho, const Real P,
-                       const Real emin, const Real emax, const Real Ye = 0.0);
+class PressResidual {
+ public:
+  KOKKOS_INLINE_FUNCTION
+  PressResidual(const EOS &eos, const Real rho, const Real P, const Real Ye)
+      : eos_(eos), rho_(rho), P_(P) {
+    lambda_[0] = Ye;
+  }
+  KOKKOS_INLINE_FUNCTION
+  Real operator()(const Real e) {
+    return eos_.PressureFromDensityInternalEnergy(rho_, e, lambda_) - P_;
+  }
+
+ private:
+  const EOS &eos_;
+  Real rho_, P_;
+  Real lambda_[2];
+};
+
+template <typename T>
+KOKKOS_FUNCTION Real energy_from_rho_P(T &eos, const Real rho, const Real P,
+                                       const Real emin, const Real emax,
+                                       const Real Ye = 0.0) {
+  PARTHENON_REQUIRE(P >= 0, "Pressure is negative!");
+  PressResidual res(eos, rho, P, Ye);
+  root_find::RootFind root;
+  Real eroot = root.regula_falsi(res, emin, emax, 1.e-10 * P, emin - 1.e10);
+  return rho * eroot;
+}
 
 } // namespace phoebus
 

src/pgen/sedov.cpp

@@ -51,7 +51,7 @@ void ProblemGenerator(MeshBlock *pmb, ParameterInput *pin) {
 
   const Real rhoa = pin->GetOrAddReal("sedov", "rho_ambient", 1.0);
   const Real rinner = pin->GetOrAddReal("sedov", "rinner", 0.01);
-  const bool spherical = pin->GetOrAddReal("sedov", "spherical_coords", true);
+  const bool spherical = pin->GetOrAddBoolean("sedov", "spherical_coords", true);
 
   auto &coords = pmb->coords;
   auto pmesh = pmb->pmy_mesh;