Create minimum energy modifier

CHANGELOG.md

@@ -14,6 +14,7 @@
 - [[PR449]](https://github.com/lanl/singularity-eos/pull/449) Ensure that DensityEnergyFromPressureTemperature works for all equations of state and is properly tested
 - [[PR439]](https://github.com/lanl/singularity-eos/pull/439) Add mean atomic mass and number to EOS API
 - [[PR437]](https://github.com/lanl/singularity-eos/pull/437) Fix segfault on HIP, clean up warnings, add strict sanitizer test
+- [[PR375]](https://github.com/lanl/singularity-eos/pull/375) Added energy floor modifier
 
 ### Changed (changing behavior/API/variables/...)
 

doc/sphinx/src/modifiers.rst

@@ -301,6 +301,42 @@ parameter via the lambda. For example:
   singularity-eos. If you would like to use it, you must implement
   your own custom variant.
 
+Energy Floor Modifier
+---------------------
+
+The models in ``singularity-eos`` do not guarantee that a lookup will lie on the
+EOS surface. The laws of thermodynamics require that volume (or density) and
+temperature be positive, but there are no natural lower bounds on other
+variables, particularly energy.
+
+In general, this modifer imposes a lower bound on the energy by finding the
+energy along the zero kelvin isotherm, i.e. :math:`e(\rho, T=0)`, and then
+flooring the energy to that value when performing a lookup. As a result, any
+lookups in energy space that would otherwise be below the zero K isotherm return
+whatever value would lie upon the isotherm.
+
+Notably though, tabular SESAME tables can result in unexpected behavior. The
+``EOSPAC`` and ``SpinerEOSDependsRhoT`` both rely on the result when EOSPAC
+querries the cold energy table (306). This energy *can* differ from the zero K
+isotherm for two reasons:
+
+#. The zero K isotherm includes the zero point energy. In other words, it is the
+   energy of the lowest energy level whereas the cold curve is the energy of the
+   potential well itself. This effect may be insignificant.
+
+#. In the absence of a 306 table for the material, EOSPAC will instead report
+   the lowest isotherm's energy in place of a cold curve. If a given table does
+   not extend to zero K, then an extrapolaiton is *not* performed. This results
+   in an energy floor from the lowest temperature on the table, **not** from
+   an extrapolated energy at zero K.
+
+.. note::
+
+  This modifier is only compatible with EOS where the ``MinInternalEnergyFromDensity``
+  function is properly defined. Where it isn't defined, an error will occur at
+  runtime.
+
+
 Composing Modifiers
 --------------------
 

singularity-eos/CMakeLists.txt

@@ -1,5 +1,5 @@
 #------------------------------------------------------------------------------#
-# © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+# © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 # program was produced under U.S. Government contract 89233218CNA000001
 # for Los Alamos National Laboratory (LANL), which is operated by Triad
 # National Security, LLC for the U.S.  Department of Energy/National
@@ -58,6 +58,7 @@ register_headers(
     eos/modifiers/ramps_eos.hpp
     eos/modifiers/shifted_eos.hpp
     eos/modifiers/eos_unitsystem.hpp
+    eos/modifiers/floored_energy.hpp
     eos/modifiers/zsplit_eos.hpp
     eos/eos_base.hpp
     eos/eos_eospac.hpp

singularity-eos/eos/eos_base.hpp

@@ -841,7 +841,18 @@ class EosBase {
     PORTABLE_ALWAYS_THROW_OR_ABORT(msg);
   }
 
-  // Default MinInternalEnergyFromDensity behavior is to cause an error
+  // Default MinInternalEnergyFromDensity behavior is to just return the zero-K isotherm.
+  // This should be fine for all thermodynamically consistent EOS, but could cause issues
+  // with EOS that aren't thermodynamically consistent.
+  template <typename Indexer_t = Real *>
+  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
+      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const {
+    CRTP copy = *(static_cast<CRTP const *>(this));
+    return copy.InternalEnergyFromDensityTemperature(rho, 0.);
+  }
+
+  // This error is useful for EOS where the zero-K approximation is invalid for whatever
+  // reason
   PORTABLE_FORCEINLINE_FUNCTION
   void MinInternalEnergyIsNotEnabled(const char *eosname) const {
     // Construct the error message using char* so it works on device

singularity-eos/eos/eos_gruneisen.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -377,9 +377,8 @@ PORTABLE_INLINE_FUNCTION Real Gruneisen::PressureFromDensityInternalEnergy(
 template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real
 Gruneisen::MinInternalEnergyFromDensity(const Real rho_in, Indexer_t &&lambda) const {
-  // const Real rho = std::min(rho_in, _rho_max);
-  MinInternalEnergyIsNotEnabled("Gruneisen");
-  return 0.0;
+  const Real rho = std::min(rho_in, _rho_max);
+  return EosBase<Gruneisen>::MinInternalEnergyFromDensity(rho);
 }
 template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real Gruneisen::EntropyFromDensityInternalEnergy(

singularity-eos/eos/eos_helmholtz.hpp

@@ -558,12 +558,6 @@ class Helmholtz : public EosBase<Helmholtz> {
             thermalqs::pressure | thermalqs::temperature, lambda);
     return p;
   }
-  template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const {
-    MinInternalEnergyIsNotEnabled("Helmholtz");
-    return 0.0;
-  }
 
   template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real

singularity-eos/eos/eos_jwl.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -73,9 +73,6 @@ class JWL : public EosBase<JWL> {
       const Real rho, const Real sie,
       Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
   template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
-  template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
   EntropyFromDensityTemperature(const Real rho, const Real temperature,
                                 Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
@@ -171,12 +168,6 @@ PORTABLE_INLINE_FUNCTION Real JWL::PressureFromDensityInternalEnergy(
   return ReferencePressure(rho) + _w * rho * (sie - ReferenceEnergy(rho));
 }
 template <typename Indexer_t>
-PORTABLE_INLINE_FUNCTION Real
-JWL::MinInternalEnergyFromDensity(const Real rho, Indexer_t &&lambda) const {
-  MinInternalEnergyIsNotEnabled("JWL");
-  return 0.0;
-}
-template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real JWL::EntropyFromDensityInternalEnergy(
     const Real rho, const Real sie, Indexer_t &&lambda) const {
   EntropyIsNotEnabled("JWL");

singularity-eos/eos/eos_mgusup.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -64,9 +64,6 @@ class MGUsup : public EosBase<MGUsup> {
       const Real rho, const Real sie,
       Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
   template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
-  template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
   EntropyFromDensityTemperature(const Real rho, const Real temp,
                                 Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
@@ -336,12 +333,6 @@ PORTABLE_INLINE_FUNCTION Real MGUsup::PressureFromDensityInternalEnergy(
   return value;
 }
 template <typename Indexer_t>
-PORTABLE_INLINE_FUNCTION Real
-MGUsup::MinInternalEnergyFromDensity(const Real rho, Indexer_t &&lambda) const {
-  MinInternalEnergyIsNotEnabled("MGUsup");
-  return 0.0;
-}
-template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real MGUsup::EntropyFromDensityInternalEnergy(
     const Real rho, const Real sie, Indexer_t &&lambda) const {
   Real eta = 1.0 - robust::ratio(_rho0, rho);

singularity-eos/eos/eos_noble_abel.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -97,7 +97,6 @@ class NobleAbel : public EosBase<NobleAbel> {
       const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const {
     return _qq;
   }
-
   template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
   EntropyFromDensityTemperature(const Real rho, const Real temperature,

singularity-eos/eos/eos_powermg.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -70,9 +70,6 @@ class PowerMG : public EosBase<PowerMG> {
       const Real rho, const Real sie,
       Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
   template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
-  template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
   EntropyFromDensityTemperature(const Real rho, const Real temp,
                                 Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
@@ -428,12 +425,6 @@ PORTABLE_INLINE_FUNCTION Real PowerMG::PressureFromDensityInternalEnergy(
   return value;
 }
 template <typename Indexer_t>
-PORTABLE_INLINE_FUNCTION Real
-PowerMG::MinInternalEnergyFromDensity(const Real rho, Indexer_t &&lambda) const {
-  MinInternalEnergyIsNotEnabled("PowerMG");
-  return 0.0;
-}
-template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real PowerMG::EntropyFromDensityInternalEnergy(
     const Real rho, const Real sie, Indexer_t &&lambda) const {
   const Real eta = 1.0 - robust::ratio(_rho0, rho);

singularity-eos/eos/eos_sap_polynomial.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -85,13 +85,6 @@ class SAP_Polynomial : public EosBase<SAP_Polynomial> {
              sie * (_b0 + _b1 * mu + _b2e * mu * mu + _b3 * mu * mu * mu);
   }
 
-  template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const {
-    MinInternalEnergyIsNotEnabled("SAP Polynomial");
-    return 0.0;
-  };
-
   template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
   EntropyFromDensityTemperature(const Real rho, const Real temperature,

singularity-eos/eos/eos_spiner.hpp

@@ -1047,8 +1047,8 @@ PORTABLE_INLINE_FUNCTION Real SpinerEOSDependsRhoT::PressureFromDensityInternalE
 template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real SpinerEOSDependsRhoT::MinInternalEnergyFromDensity(
     const Real rho, Indexer_t &&lambda) const {
-  MinInternalEnergyIsNotEnabled("SpinerEOSDependsRhoT");
-  return 0.0;
+  const Real lRho = lRho_(rho);
+  return sieCold_.interpToReal(lRho);
 }
 template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real SpinerEOSDependsRhoT::EntropyFromDensityTemperature(

singularity-eos/eos/eos_stiff.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -90,13 +90,6 @@ class StiffGas : public EosBase<StiffGas> {
     return std::max(-_Pinf, _gm1 * rho * (sie - _qq) - (_gm1 + 1.0) * _Pinf);
   }
 
-  template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const {
-    MinInternalEnergyIsNotEnabled("StiffGas");
-    return 0.0;
-  };
-
   template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
   EntropyFromDensityTemperature(const Real rho, const Real temperature,

singularity-eos/eos/eos_vinet.hpp

@@ -1,5 +1,5 @@
 //------------------------------------------------------------------------------
-// © 2021-2024. Triad National Security, LLC. All rights reserved.  This
+// © 2021-2025. Triad National Security, LLC. All rights reserved.  This
 // program was produced under U.S. Government contract 89233218CNA000001
 // for Los Alamos National Laboratory (LANL), which is operated by Triad
 // National Security, LLC for the U.S.  Department of Energy/National
@@ -67,9 +67,6 @@ class Vinet : public EosBase<Vinet> {
   PORTABLE_INLINE_FUNCTION Real PressureFromDensityInternalEnergy(
       const Real rho, const Real sie,
       Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
-  template <typename Indexer_t = Real *>
-  PORTABLE_INLINE_FUNCTION Real MinInternalEnergyFromDensity(
-      const Real rho, Indexer_t &&lambda = static_cast<Real *>(nullptr)) const;
   // Entropy added AEM Dec. 2022
   template <typename Indexer_t = Real *>
   PORTABLE_INLINE_FUNCTION Real
@@ -359,12 +356,6 @@ PORTABLE_INLINE_FUNCTION Real Vinet::PressureFromDensityInternalEnergy(
   return output[1];
 }
 template <typename Indexer_t>
-PORTABLE_INLINE_FUNCTION Real
-Vinet::MinInternalEnergyFromDensity(const Real rho, Indexer_t &&lambda) const {
-  MinInternalEnergyIsNotEnabled("Vinet");
-  return 0.0;
-}
-template <typename Indexer_t>
 PORTABLE_INLINE_FUNCTION Real Vinet::EntropyFromDensityInternalEnergy(
     const Real rho, const Real sie, Indexer_t &&lambda) const {
   Real temp;