Analytic electron EOS and zsplit modifier

.github/workflows/formatting.yml

@@ -9,7 +9,7 @@ on:
 jobs:
     formatting:
       name: Check Formatting
-      runs-on: ubuntu-latest
+      runs-on: ubuntu-noble
 
       steps:      
         - name: Checkout code

CHANGELOG.md

@@ -3,6 +3,7 @@
 ## Current develop
 
 ### Added (new features/APIs/variables/...)
+- [[PR444]](https://github.com/lanl/singularity-eos/pull/444) Add Z split modifier and electron ideal gas EOS
 
 ### Fixed (Repair bugs, etc)
 - [[PR439]](https://github.com/lanl/singularity-eos/pull/439) Add mean atomic mass and number to EOS API

CMakeLists.txt

@@ -440,8 +440,8 @@ if(SINGULARITY_BUILD_TESTS)
     FetchContent_Declare(
       Catch2
       GIT_REPOSITORY https://github.com/catchorg/Catch2.git
-      # or later is fine too
-      GIT_TAG v3.0.1)
+      # jmm: updated Dec 17, 2024 to avoid build errors on modern gcc
+      GIT_TAG v3.7.1)
     FetchContent_MakeAvailable(Catch2)
     list(APPEND CMAKE_MODULE_PATH ${Catch2_SOURCE_DIR}/contrib)
   endif()
@@ -560,6 +560,10 @@ target_compile_options(
          -use_fast_math
          > # release
          > # cuda
+
+         # Suppresses annoying ABI notes. See:
+         # https://stackoverflow.com/questions/52020305/what-exactly-does-gccs-wpsabi-option-do-what-are-the-implications-of-supressi
+         $<$<CXX_COMPILER_ID:GNU>:-Wno-psabi>
 )
 if (SINGULARITY_STRICT_WARNINGS)
   target_compile_options(singularity-eos_Interface INTERFACE

doc/sphinx/src/models.rst

@@ -164,6 +164,61 @@ functional forms for :math:`\Gamma` and the reference curves, the task of
 calculating a **thermodynamically consistent** temperature becomes more
 complicated.
 
+The 3T Model
+``````````````
+
+In high-energy denisty physics, the so-called three-temperature (or
+3T) model is often employed. In this model, a material is fully or
+partially ionized. Due to the separation in masses between electrons
+and atomic nuclei, electrons come into thermal equilibrium and ions
+come into thermal equilibrium, but the two populations are not in
+thermal equilibrium with each other. This introduces separate electron
+and ion temperatures. The third temperature is potentially a radiation
+temperature.
+
+By convention in 3T physics, density :math:`\rho` is always the mass
+density of the ions. The electron number density may then be computed
+as
+
+.. math::
+
+  n = \left\langle Z \right\rangle \frac{\rho}{\left\langle A\right\rangle m_p}
+
+where here :math:`\left\langle Z\right\rangle` is the average number
+of electrons contributed per atom, also called mean ionization state,
+:math:`\rho` is the ion mass density, :math:`\bar{A}` is the mean
+atomic mass (in grams per mole) of a given material and :math:`m_p` is
+the proton mass.
+
+.. note::
+
+  Note that there is notational ambiguity between the average
+  ionization state and the average atomic number, as the symbol for
+  both is :math:`Z`. To disambiguate in ``singularity-eos``, we use
+  overbars to reference mean atomic properties such as mean atomic
+  mass :math:`\bar{A}` and mean atomic number :math:`\bar{Z}` while we
+  use :math:`\left\langle Z\right\rangle` to denote mean ionizaiton
+  state.
+
+Also by convention, the specific internal energy carried by electrons
+:math:`\varepsilon_e` is specific with respect to the *ion* mass. In
+particular:
+
+.. math::
+
+  u_e = \varepsilon_e \rho
+
+where here :math:`u_e` is the internal energy per unit volume carried
+by electrons and, as discussed above, :math:`\rho` is ion mass
+density.
+
+``singularity-eos`` assumes that, when 3T physics is active, electrons
+and ions are each described by a separate equation of state
+object. Several models are specifically designed to represent, e.g.,
+the electron equation of state or ion equation of state. The tabulated
+models may also support loading tables specifically for electron or
+ion equations of state.
+
 Available EOS Information and Nomenclature
 ------------------------------------------
 
@@ -469,6 +524,71 @@ the atomic mass and number as a final optional parameter, e.g.,
    IdealGas(Real gm1, Real Cv, MeanAtomicProperties(Abar, Zbar));
    IdealGas(Real gm1, Real Cv, Real EntropyT0, Real EntropyRho0, MeanAtomicProperties(Abar, Zbar));
 
+Ideal Electron Gas
+```````````````````
+
+The ideal electron gas equation of state is designed for use with 3T
+physics. It is an ideal Boltzmann gas of electrons. As such, each
+electron is assumed to have three translational degrees of freedom
+:math:`f`, such that
+
+.. math::
+
+  \Gamma = \frac{2}{f} = \frac{2}{3}
+
+and
+
+.. math::
+
+  \gamma = \Gamma + 1 = \frac{5}{3}.
+
+The pressure is given by the number density of electrons times
+:math:`k_b T` for Boltzmann constant :math:`k_b` and temperature
+:math:`T`:
+
+.. math::
+
+  P = \frac{\rho}{m_p \bar{A}}\left\langle Z \right\rangle k_b T
+
+The specific heat is then
+
+.. math::
+
+  C_V = \frac{\left\langle Z \right\rangle k_b}{\Gamma m_p \bar{A}}
+
+so that
+
+.. math::
+
+  P = \Gamma \rho C_V T
+
+as expected.
+
+The constructor takes only the ``MeanAtomicProperties`` struct, which
+is a required input:
+
+.. code-block::cpp
+
+  IdealElectrons(const MeanAtomicProperties &AZbar);
+
+Optionally reference values may be provided for the entropy
+calculation, which is computed in the same way as the standard ideal
+gas.
+
+Calls to compute state variables require the mean ionization state,
+which must be passed in the ``lambda`` parameter, e.g.,
+
+.. code-block:: cpp
+
+  Real lambda[1] = {Z};
+  Real P = eos.PressureFromDensityTemperature(rho, T, lambda);
+
+.. note::
+
+  For now, the ideal electron gas is not in the default variant
+  provided by singularity-eos. If you would like to use it, you must
+  implement your own custom variant.
+
 Stiffened Gas
 `````````````
 

doc/sphinx/src/modifiers.rst

@@ -8,14 +8,16 @@ Equation of State Modifiers
 provides some internal transformation on inputs and outputs. For
 example the ``ShiftedEOS`` modifier changes the zero point energy of a
 given EOS model by shifting all energies up or down. Modifiers can be
-used to, for example, production-harden a model. Only certain
-combinations of ``EOS`` and ``modifier`` are permitted by the defualt
-``Variant``. For example, only ``IdealGas``, ``SpinerEOS``, and
-``StellarCollapse`` support the ``RelativisticEOS`` and ``UnitSystem``
-modifiers. All models support the ``ShiftedEOS`` and ``ScaledEOS``
-modifiers. However, note that modifiers do not commute, and only one
-order is supported. The ordering, inside-out, is ``UnitSystem`` or
-``RelativisticEOS``, then ``ScaledEOS``, then ``ShiftedEOS``.
+used to, for example, production-harden a model.
+
+Only certain combinations of ``EOS`` and ``modifier`` are permitted by
+the defualt ``Variant``. For example, only ``IdealGas``,
+``SpinerEOS``, and ``StellarCollapse`` support the ``RelativisticEOS``
+and ``UnitSystem`` modifiers. All models support the ``ShiftedEOS``
+and ``ScaledEOS`` modifiers. However, note that modifiers do not
+commute, and only one order is supported. The ordering, inside-out, is
+``UnitSystem`` or ``RelativisticEOS``, then ``ScaledEOS``, then
+``ShiftedEOS``.
 
 We list below the available modifiers and their constructors.
 
@@ -187,6 +189,99 @@ internal energy, and temperature. On the other hand,
 specifies the unit system by specifying units for time, mass, length,
 and temperature.
 
+Z-Split EOS
+-------------
+
+For 3T physics (as described in the models section) it is often
+desirable to have a separate equation of state for electrons and a
+separate equation of state for ions. The Z-split model takes a total
+equation of state and splits it into electron and ion components. The
+"Z" here signifies mean ionization state, or the average number of
+free electrons contributed per atomic nucleus, which is the mean
+atomic number in the case of full ionization, but could be smaller in
+the case of partial ionization. (It is zero for an unionized gas.) The
+physical model of Z-split can be derived from an ideal gas equation of
+state. For an ideal gas made up of electrons and ions, where all
+molecular bonds have been broken, the total pressure is given by
+
+.. math::
+
+  P_t = (\left\langle Z\right\rangle + 1) \frac{\rho}{m_p \bar{A}} k_b T
+
+where :math:`\left\langle Z\right\rangle` is the mean ionization
+state, :math:`rho` is the ion mass density (the electron ion mass
+density is negligible), :math:`m_p` is the proton mass,
+:math:`\bar{A}` is the mean atomic mass, :math:`k_b` is Boltzmann's
+constant, and :math:`T` is temperature. The contribution from
+electrons is proportional to :math:`\left\langle Z\right\rangle`.
+
+The split simply splits the total pressure into normalized
+contributions for the electrons:
+
+.. math::
+
+  P_e = \frac{\left\langle Z\right\rangle}{1 + \left\langle Z\right\rangle} P_t(\rho, T_e)
+
+and ions:
+
+.. math::
+
+  P_i = \frac{1}{1 + \left\langle Z\right\rangle} P_t(\rho, T_i)
+
+where here :math:`T_e` is the electron temperature and :math:`T_i` the
+ion temperature such that, when the temperatures are equal,
+
+.. math::
+
+  P_e + P_i = P_t
+
+The same split is applied to the specific internal energy:
+
+.. math::
+
+  \varepsilon_e = \frac{\left\langle Z\right\rangle}{1 + \left\langle Z\right\rangle} \varepsilon_t(\rho, T_e)\\
+  \varepsilon_i = \frac{1}{1 + \left\langle Z\right\rangle} \varepsilon_t(\rho, T_i)
+
+and the remaining state variables and thermodynamic derivatives can be
+derived from these relations.
+
+In ``singularity-eos``, the z-split is implemented as the templated class
+
+.. code-block:: cpp
+
+  template<ZSplitComponent ztype, typename T>
+  class ZSplit
+
+where ``ZSPlitComponent`` may either be ``ZSplitComponent::Electrons``
+or ``ZSplitComponent::Ions``. As syntactic sugar, ``ZSplitE<T>`` and
+``ZSplitI<T>`` are available. The Z-split constructor does not require
+any additional parameters, so you may construct one as, e.g.,
+
+.. code-block:: cpp
+
+  using namespace singularity
+  auto ion_eos = ZSPlitI<IdealGas>(IdealGas(gm1, Cv);
+
+and similarly for electrons,
+
+.. code-block:: cpp
+
+  auto electron_eos = ZSPlitE<IdealGas>(IdealGas(gm1, Cv);
+
+The Z-split modifier takes the ionization state as an additional
+parameter via the lambda. For example:
+
+.. code-block:: cpp
+
+  Real lambda[1] = {Z};
+  Real Pe = electron_eos.PressureFromDensityTemperature(rho, temperature, lambda);
+
+.. note::
+
+  For now, the Z-split EOS is not in the default variant provided by
+  singularity-eos. If you would like to use it, you must implement
+  your own custom variant.
+
 Composing Modifiers
 --------------------
 

doc/sphinx/src/using-eos.rst

@@ -737,6 +737,18 @@ For more details on modifiers, see the :ref:`modifiers<modifiers>`
 section. If you need a combination of modifiers not supported by
 default, we recommend building a custom variant as described above.
 
+Modifiers and Lambdas
+-----------------------
+
+Modifiers may require lambdas. When this is the case, the lambda
+required by the modifier is appended to the end of the lambda
+indexer. For example, the ``StellarCollapse`` EOS model requires
+``nlambda=2``. The ``ZSplitI`` modifier rquires
+``nlambda=1``. Together, ``ZSplitI<StellarCollapse>`` requires a
+lambda indexer of length 3, an the ordering is two parameters for
+``StellarCollapse`` first, and then the parameter required by
+``ZSplitI``.
+
 Preferred Inputs
 -----------------
 

singularity-eos/CMakeLists.txt

@@ -57,13 +57,15 @@ register_headers(
     eos/modifiers/ramps_eos.hpp
     eos/modifiers/shifted_eos.hpp
     eos/modifiers/eos_unitsystem.hpp
+    eos/modifiers/zsplit_eos.hpp
     eos/eos_base.hpp
     eos/eos_eospac.hpp
     eos/eos_noble_abel.hpp
     eos/eos_stiff.hpp
     eos/singularity_eos_init_utils.hpp
     eos/variant_utils.hpp
     eos/eos_carnahan_starling.hpp
+    eos/eos_electrons.hpp
 )
 
 if (SINGULARITY_BUILD_CLOSURE)

singularity-eos/eos/eos_base.hpp

@@ -66,23 +66,24 @@ char *StrCat(char *destination, const char *source) {
 // VECTOR functionality to overload the scalar implementations in the derived
 // classes. Do not add functions here that are not overloads of derived class features.
 // TODO(JMM): Should we have more macros that capture just some of these?
-#define SG_ADD_BASE_CLASS_USINGS(EOSDERIVED)                                             \
-  using EosBase<EOSDERIVED>::TemperatureFromDensityInternalEnergy;                       \
-  using EosBase<EOSDERIVED>::InternalEnergyFromDensityTemperature;                       \
-  using EosBase<EOSDERIVED>::PressureFromDensityTemperature;                             \
-  using EosBase<EOSDERIVED>::PressureFromDensityInternalEnergy;                          \
-  using EosBase<EOSDERIVED>::MinInternalEnergyFromDensity;                               \
-  using EosBase<EOSDERIVED>::SpecificHeatFromDensityTemperature;                         \
-  using EosBase<EOSDERIVED>::SpecificHeatFromDensityInternalEnergy;                      \
-  using EosBase<EOSDERIVED>::BulkModulusFromDensityTemperature;                          \
-  using EosBase<EOSDERIVED>::BulkModulusFromDensityInternalEnergy;                       \
-  using EosBase<EOSDERIVED>::GruneisenParamFromDensityTemperature;                       \
-  using EosBase<EOSDERIVED>::GruneisenParamFromDensityInternalEnergy;                    \
-  using EosBase<EOSDERIVED>::FillEos;                                                    \
-  using EosBase<EOSDERIVED>::EntropyFromDensityTemperature;                              \
-  using EosBase<EOSDERIVED>::EntropyFromDensityInternalEnergy;                           \
-  using EosBase<EOSDERIVED>::GibbsFreeEnergyFromDensityTemperature;                      \
-  using EosBase<EOSDERIVED>::GibbsFreeEnergyFromDensityInternalEnergy;
+// JMM: Use VA_ARGS to capture more complex template types
+#define SG_ADD_BASE_CLASS_USINGS(...)                                                    \
+  using EosBase<__VA_ARGS__>::TemperatureFromDensityInternalEnergy;                      \
+  using EosBase<__VA_ARGS__>::InternalEnergyFromDensityTemperature;                      \
+  using EosBase<__VA_ARGS__>::PressureFromDensityTemperature;                            \
+  using EosBase<__VA_ARGS__>::PressureFromDensityInternalEnergy;                         \
+  using EosBase<__VA_ARGS__>::MinInternalEnergyFromDensity;                              \
+  using EosBase<__VA_ARGS__>::SpecificHeatFromDensityTemperature;                        \
+  using EosBase<__VA_ARGS__>::SpecificHeatFromDensityInternalEnergy;                     \
+  using EosBase<__VA_ARGS__>::BulkModulusFromDensityTemperature;                         \
+  using EosBase<__VA_ARGS__>::BulkModulusFromDensityInternalEnergy;                      \
+  using EosBase<__VA_ARGS__>::GruneisenParamFromDensityTemperature;                      \
+  using EosBase<__VA_ARGS__>::GruneisenParamFromDensityInternalEnergy;                   \
+  using EosBase<__VA_ARGS__>::FillEos;                                                   \
+  using EosBase<__VA_ARGS__>::EntropyFromDensityTemperature;                             \
+  using EosBase<__VA_ARGS__>::EntropyFromDensityInternalEnergy;                          \
+  using EosBase<__VA_ARGS__>::GibbsFreeEnergyFromDensityTemperature;                     \
+  using EosBase<__VA_ARGS__>::GibbsFreeEnergyFromDensityInternalEnergy;
 
 // This macro adds these methods to a derived class. Due to scope,
 // these can't be implemented in the base class, unless we make
@@ -111,7 +112,9 @@ char *StrCat(char *destination, const char *source) {
   }                                                                                      \
   constexpr bool AllDynamicMemoryIsShareable() const {                                   \
     return t_.AllDynamicMemoryIsShareable();                                             \
-  }                                                                                      \
+  }
+
+#define SG_ADD_MODIFIER_MEAN_METHODS(t_)                                                 \
   PORTABLE_INLINE_FUNCTION                                                               \
   Real MeanAtomicMass() const { return t_.MeanAtomicMass(); }                            \
   PORTABLE_INLINE_FUNCTION                                                               \