Merge branch 'main' of github.com:lanl/singularity-eos into apg/sap_i…

…ntegration

CHANGELOG.md

@@ -5,6 +5,7 @@
 
 ### Fixed (Repair bugs, etc)
 - [[PR278]](https://github.com/lanl/singularity-eos/pull/278) Fixed EOSPAC unit conversion errors for scalar lookups
+- [[PR316]](https://github.com/lanl/singularity-eos/pull/316) removed `fmax-errors=3` from `singularity-eos` compile flags
 - [[PR296]](https://github.com/lanl/singularity-eos/pull/296) changed `CMAKE_SOURCE_DIR` to `PROJECT_SOURCE_DIR` to fix downstream submodule build
 - [[PR291]](https://github.com/lanl/singularity-eos/pull/291) package.py updates to reflect new CMake options 
 - [[PR290]](https://github.com/lanl/singularity-eos/pull/290) Added target guards on export config
@@ -26,6 +27,7 @@
 - [[PR278]](https://github.com/lanl/singularity-eos/pull/278) Added EOSPAC option functionality in class constructor
 - [[PR278]](https://github.com/lanl/singularity-eos/pull/278) Added a new function for returning the minimum energy as a function of density for an EOS (only EOSPAC at the moment)
 - [[PR278]](https://github.com/lanl/singularity-eos/pull/278) Added a new Fortran API for simple pressure and bulk moduli lookups
+- [[PR306]](https://github.com/lanl/singularity-eos/pull/306) Added generic Evaluate method
 - [[PR304]](https://github.com/lanl/singularity-eos/pull/304) added a Newton-Raphson root find for use with the Helmholtz EoS
 - [[PR265]](https://github.com/lanl/singularity-eos/pull/265) Add missing UnitSystem modifier combinations to variant and EOSPAC
 - [[PR279]](https://github.com/lanl/singularity-eos/pull/279) added noble-abel EoS
@@ -43,6 +45,7 @@
 - [[PR177]](https://github.com/lanl/singularity-eos/pull/177) added EOSPAC vector functions
 
 ### Changed (changing behavior/API/variables/...)
+- [[PR310]](https://github.com/lanl/singularity-eos/pull/310) Speed up and clean up tests
 - [[PR295]](https://github.com/lanl/singularity-eos/pull/295) Add fast logs to singularity-eos
 - [[PR246]](https://github.com/lanl/singularity-eos/pull/246) Update CMake with upstream changes
 - [[PR223]](https://github.com/lanl/singularity-eos/pull/223) Update ports-of-call and add portable error handling

CMakeLists.txt

@@ -84,8 +84,9 @@ cmake_dependent_option(
   "SINGULARITY_BUILD_TESTS;SINGUALRITY_BUILD_STELLARCOLLAPSE2SPINER" OFF)
 cmake_dependent_option(SINGULARITY_TEST_PYTHON "Test the Python bindings" ON
                        "SINGULARITY_BUILD_TESTS;SINGULARITY_BUILD_PYTHON" OFF)
-cmake_dependent_option(SINGULARITY_TEST_HELMHOLTZ "Test the Helmholtz equation of state" ON
-                       "SINGULARITY_BUILD_TESTS;SINGULARITY_USE_HELMHOLTZ" OFF)
+cmake_dependent_option(
+  SINGULARITY_TEST_HELMHOLTZ "Test the Helmholtz equation of state" ON
+  "SINGULARITY_BUILD_TESTS;SINGULARITY_USE_HELMHOLTZ" OFF)
 
 # modify flags options
 option(SINGULARITY_BETTER_DEBUG_FLAGS "Better debug flags for singularity" ON)
@@ -96,11 +97,12 @@ option(SINGULARITY_USE_SINGLE_LOGS
        "Use single precision logs. Can harm accuracy." OFF)
 option(SINGULARITY_USE_TRUE_LOG_GRIDDING
        "Use grids that conform to log spacing." OFF)
-# TODO(JMM): Should this automatically be activated when true log
-# gridding is off?
-cmake_dependent_option(SINGULARITY_USE_HIGH_RISK_MATH
-        "Use integer aliased logs, may not be portable" OFF
-        "NOT SINGULARITY_USE_TRUE_LOG_GRIDDING" OFF)
+# TODO(JMM): Should this automatically be activated when true log gridding is
+# off?
+cmake_dependent_option(
+  SINGULARITY_USE_HIGH_RISK_MATH
+  "Use integer aliased logs, may not be portable" OFF
+  "NOT SINGULARITY_USE_TRUE_LOG_GRIDDING" OFF)
 
 # misc options
 option(SINGULARITY_FORCE_SUBMODULE_MODE "Submodule mode" OFF)
@@ -236,8 +238,9 @@ endif()
 if(SINGULARITY_USE_SINGLE_LOGS)
   target_compile_definitions(singularity-eos PUBLIC SINGULARITY_USE_SINGLE_LOGS)
 endif()
-if (SINGULARITY_USE_HIGH_RISK_MATH)
-  target_compile_definitions(singularity-eos PUBLIC SINGULARITY_USE_HIGH_RISK_MATH)
+if(SINGULARITY_USE_HIGH_RISK_MATH)
+  target_compile_definitions(singularity-eos
+                             PUBLIC SINGULARITY_USE_HIGH_RISK_MATH)
 endif()
 
 if(SINGULARITY_TEST_SESAME)
@@ -448,7 +451,7 @@ target_compile_options(
          -use_fast_math
          > # release
          > # cuda
-  PUBLIC -fmax-errors=3)
+)
 
 target_link_options(singularity-eos PRIVATE ${xlfix})
 

doc/sphinx/src/contributing.rst

@@ -233,10 +233,9 @@ test in the ``CMakeLists.txt`` file,
 
 .. code-block:: cmake
 
-    add_executable(eos_unit_tests
+    add_executable(eos_analytic_unit_tests
         catch2_define.cpp
         eos_unit_test_helpers.hpp
-        test_eos_unit.cpp
         test_eos_gruneisen.cpp
         test_eos_vinet.cpp
         test_my_new_eos.cpp
@@ -245,15 +244,23 @@ test in the ``CMakeLists.txt`` file,
 in order for the test to be compiled. If your EOS requires any special
 dependencies, be sure to block off the test using ``#IFDEF`` blocks.
 
-**Important:** this is a subtlety that highlights the importance of unit tests!
-Since our library is header only, the unit tests are often the only place where
-a specific EOS may be instantiated when ``singularity-eos`` is compiled. Unit
-tests _must_ make use of the ``EOS`` type, i.e.
+.. note::
+
+  Note that there are three executables, ``eos_analytic_unit_tests``,
+  ``eos_infrastructure_tests`` and ``eos_tabulated_unit_tests``. Pick
+  the executable that most closely matches what your model is.
+
+**Important:** Since our library is header only, the unit
+tests are often the only place where a specific EOS may be
+instantiated when ``singularity-eos`` is compiled.  Therefore to
+exercise all code paths, it is best to create an ``EOS`` type
+instantiated as
 
 .. code-block:: c++
 
     #include <singularity-eos/eos/eos.hpp>
-    EOS my_eos = my_new_eos(parameter1, parameter2, ...)
+    using EOS = singularity::Variant<MyNewEOS>;``.
+    EOS my_eos = MyNewEOS(parameter1, parameter2, ...)
 
 in order to properly test the functionality of a new EOS. Simply using the
 new class as the type such as

doc/sphinx/src/using-eos.rst

@@ -35,6 +35,10 @@ conditions. The type of parallelism used depends on how
 ``singularity-eos`` is compiled. If the ``Kokkos`` backend is used,
 any parallel dispatch supported by ``Kokkos`` is supported.
 
+A more generic version of the vector calls exists in the ``Evaluate``
+method, which allows the user to specify arbitrary parallel dispatch
+models by writing their own loops. See the relevant section below.
+
 .. _variant section:
 
 Variants
@@ -216,6 +220,121 @@ available member function.
    eos.TemperatureFromDensityInternalEnergy(density.data(), energy.data(), temperature.data(),
                                             scratch.data(), density.size());
 
+The Evaluate Method
+~~~~~~~~~~~~~~~~~~~
+
+A special call related to the vector calls is the ``Evaluate``
+method. The ``Evaluate`` method requests the EOS object to evaluate
+almost arbitrary code, but in a way where the type of the underlying
+EOS object is resolved *before* this arbitrary code is evaluated. This
+means the code required to resolve the type of the variant is only
+executed *once* per ``Evaluate`` call. This can enable composite EOS
+calls, non-standard vector calls, and vector calls with non-standard
+loop structure.
+
+The ``Evaluate`` call has the signature
+
+.. code-block:: cpp
+
+  template<typename Functor_t>
+  PORTABLE_INLINE_FUNCTION
+  void Evaluate(Functor_t f);
+
+where a ``Functor_t`` is a class that *must* provide a ``void
+operator() const`` method templated on EOS type. ``Evaluate`` is
+decorated so that it may be evaluated on either host or device,
+depending on desired use-case. Alternatively, you may use an anonymous
+function with an `auto` argument as the input, e.g.,
+
+.. code-block::
+
+   // equivalent to [=], but with device markings
+   eos.Evaluate(PORTABLE_LAMBDA(auto eos) { /* my code snippet */ });
+
+.. warning::
+
+  It can be dangerous to use functors with side-effects. Especially
+  with GPUs it can produce very unintuitive behaviour. We recommend
+  you only make the ``operator()`` non-const if you really know what
+  you're doing. And in the anonymous function case, we recommend you
+  capture by value, not reference.
+
+To see the utlity of the ``Evaluate`` function, it's probably just
+easiest to provide an example. The following code evaluates the EOS on
+device and compares against a tabulated pressure. The total difference
+is summed using the ``Kokkos::parallel_reduce`` functionality in the
+``Kokkos`` performance portability library.
+
+.. code-block:: cpp
+
+  // The functor we use is defined here.
+  // This class definition needs to be of appropriately global scope.
+  class CheckPofRE {
+   public:
+    CheckPofRE(Real *P, Real *rho, Real *sie, int N) : P_(P), rho_(rho), sie_(sie), N_(N) {}
+    template <typename T>
+    // this is a host-only call, but if you wanted to write
+    // a function that you wanted to evaluate on device
+    // you could add the
+    // PORTABLE_INLINE_FUNCTION
+    // decorator here.
+    void operator()(const T &eos) const {
+      // Capturing member functions of a class in a lambda typically causes problems
+      // when launching a GPU kernel.
+      // Better to pull out new variables to capture before launching a kernel.
+      Real *P = P_;
+      Real *rho = rho_;
+      Real *sie = sie_;
+      // reduction target
+      Real tot_diff;
+      // reduction op
+      Kokkos::parallel_reduce(
+          "MyCheckPofRE", N_,
+          KOKKOS_LAMBDA(const int i, Real &diff) {
+            diff += std::abs(P[i] - eos.PressureFromDensityInternalEnergy(rho[i], sie[i]));
+          },
+          tot_diff);
+      std::cout << "Total difference = " << tot_diff << std::endl;
+    }
+  
+   private:
+    int N_;
+    Real *P_;
+    Real *rho_;
+    Real *sie_;
+  };
+
+  // Here we construct our functor
+  // it is assumed the pointers to device memory P, rho, sie, are defined elsewhere.
+  CheckPofRE my_op(P, rho, sie, N);
+
+  // Here we call the evaluate function
+  eos.Evaluate(my_op);
+
+  // The above two lines could have been called "in-one" with:
+  // eos.Evaluate(CheckPofRE(P, rho, sie, N));
+
+Alternatively, you could eliminate the functor and use an anonymous
+function with:
+
+.. code-block:: cpp
+
+  eos.Evaluate([=](auto eos) {
+    Real tot_diff;
+    Kokkos::parallel_reduce(
+        "MyCheckPofRE", N_,
+        KOKKOS_LAMBDA(const int i, Real &diff) {
+          diff += std::abs(P[i] - eos.PressureFromDensityInternalEnergy(rho[i], sie[i]));
+        },
+        tot_diff);
+    std::cout << "Total difference = " << tot_diff << std::endl;
+  });
+
+This is not functionality that would be available with the standard
+vector calls provided by ``singularity-eos``, at least not without
+chaining multiple parallel dispatch calls. Here we can do it in a
+single call.
+
 Lambdas and Optional Parameters
 --------------------------------
 

singularity-eos/CMakeLists.txt

@@ -33,6 +33,7 @@ set(EOS_HEADERS
     eos/eos_builder.hpp
     eos/eos_jwl.hpp
     eos/eos_helmholtz.hpp
+    eos/eos_sap_polynomial.hpp
     eos/modifiers/relativistic_eos.hpp
     eos/modifiers/scaled_eos.hpp
     eos/modifiers/ramps_eos.hpp

singularity-eos/closure/mixed_cell_models.hpp

@@ -272,9 +272,10 @@ class PTESolverBase {
     return Tguess;
   }
 
-  PORTABLE_FORCEINLINE_FUNCTION
-  static Real GetPressureFromPreferred(const EOS &eos, const Real rho, const Real T,
-                                       Real sie, Real *lambda, const bool do_e_lookup) {
+  template <typename EOS_t>
+  PORTABLE_FORCEINLINE_FUNCTION static Real
+  GetPressureFromPreferred(const EOS_t &eos, const Real rho, const Real T, Real sie,
+                           Real *lambda, const bool do_e_lookup) {
     Real P{};
     if (eos.PreferredInput() ==
         (thermalqs::density | thermalqs::specific_internal_energy)) {

singularity-eos/eos/eos_base.hpp

@@ -130,6 +130,13 @@ class EosBase {
   struct is_raw_pointer
       : std::is_same<std::remove_reference_t<std::remove_cv_t<T>>, R *> {};
 
+  // Generic evaluator
+  template <typename Functor_t>
+  constexpr void Evaluate(Functor_t &f) const {
+    CRTP copy = *(static_cast<CRTP const *>(this));
+    f(copy);
+  }
+
   // Vector member functions
   template <typename RealIndexer, typename ConstRealIndexer, typename LambdaIndexer>
   inline void

singularity-eos/eos/eos_variant.hpp

@@ -62,7 +62,7 @@ class Variant {
             typename std::enable_if<
                 !std::is_same<Variant, typename std::decay<EOSChoice>::type>::value,
                 bool>::type = true>
-  EOSChoice get() {
+  PORTABLE_INLINE_FUNCTION EOSChoice get() {
     return mpark::get<EOSChoice>(eos_);
   }
 
@@ -72,6 +72,11 @@ class Variant {
   }
 
   // Place member functions here
+  template <typename Functor_t>
+  constexpr void Evaluate(Functor_t &f) const {
+    return mpark::visit([&f](const auto &eos) { return eos.Evaluate(f); }, eos_);
+  }
+
   PORTABLE_INLINE_FUNCTION
   Real TemperatureFromDensityInternalEnergy(const Real rho, const Real sie,
                                             Real *lambda = nullptr) const {

test/CMakeLists.txt

@@ -12,35 +12,58 @@
 # ------------------------------------------------------------------------------
 
 add_executable(
-  eos_unit_tests
+  eos_analytic_unit_tests
   catch2_define.cpp
   eos_unit_test_helpers.hpp
-  test_eos_unit.cpp
+  test_eos_ideal.cpp
   test_eos_gruneisen.cpp
   test_eos_sap_polynomial.cpp
-  test_eos_vinet.cpp
   test_eos_noble_abel.cpp
   test_eos_stiff.cpp
+  )
+
+add_executable(
+  eos_infrastructure_tests
+  catch2_define.cpp
+  eos_unit_test_helpers.hpp
+  test_eos_modifiers.cpp
+  test_eos_vector.cpp
+  test_math_utils.cpp
+  )
+
+add_executable(
+  eos_tabulated_unit_tests
+  catch2_define.cpp
+  eos_unit_test_helpers.hpp
+  test_eos_vinet.cpp
   test_eos_helmholtz.cpp
+  test_eos_tabulated.cpp
+  test_eos_stellar_collapse.cpp
   )
 
 if(SINGULARITY_TEST_HELMHOLTZ)
     configure_file(${PROJECT_SOURCE_DIR}/data/helmholtz/helm_table.dat ${CMAKE_BINARY_DIR}/data/helmholtz/helm_table.dat COPYONLY)
-  target_compile_definitions(eos_unit_tests PRIVATE SINGULARITY_TEST_HELMHOLTZ SINGULARITY_USE_HELMHOLTZ)
+  target_compile_definitions(eos_tabulated_unit_tests PRIVATE SINGULARITY_TEST_HELMHOLTZ SINGULARITY_USE_HELMHOLTZ)
 endif()
 
 if(SINGULARITY_TEST_SESAME)
-  target_compile_definitions(eos_unit_tests PRIVATE SINGULARITY_TEST_SESAME)
+  target_compile_definitions(eos_tabulated_unit_tests PRIVATE SINGULARITY_TEST_SESAME)
 endif()
 if(SINGULARITY_TEST_STELLAR_COLLAPSE)
-  target_compile_definitions(eos_unit_tests
+  target_compile_definitions(eos_tabulated_unit_tests
                              PRIVATE SINGULARITY_TEST_STELLAR_COLLAPSE)
 endif()
 
-target_link_libraries(eos_unit_tests PRIVATE Catch2::Catch2
-                                             singularity-eos::singularity-eos)
+target_link_libraries(eos_analytic_unit_tests PRIVATE Catch2::Catch2
+  singularity-eos::singularity-eos)
+target_link_libraries(eos_infrastructure_tests PRIVATE Catch2::Catch2
+  singularity-eos::singularity-eos)
+target_link_libraries(eos_tabulated_unit_tests PRIVATE Catch2::Catch2
+  singularity-eos::singularity-eos)
 include(Catch)
-catch_discover_tests(eos_unit_tests PROPERTIES TIMEOUT 60)
+catch_discover_tests(eos_analytic_unit_tests PROPERTIES TIMEOUT 60)
+catch_discover_tests(eos_infrastructure_tests PROPERTIES TIMEOUT 60)
+catch_discover_tests(eos_tabulated_unit_tests PROPERTIES TIMEOUT 60)
 
 if(SINGULARITY_USE_EOSPAC AND SINGULARITY_TEST_SESAME)
   add_executable(compare_to_eospac compare_to_eospac.cpp)

test/catch2_define.cpp

@@ -16,6 +16,7 @@
 
 #ifndef CATCH_CONFIG_RUNNER
 #define CATCH_CONFIG_RUNNER
+#define CATCH_CONFIG_FAST_COMPILE
 #include "catch2/catch.hpp"
 #endif