Test refactor. The getsgeos tests are now in their own file that resp…

…ect build restrictions that reduce compile times. PTE tests still built under normal conditions.
lanl · Aug 14, 2023 · 9774621 · 9774621
1 parent ed4e55a
commit 9774621
Show file tree

Hide file tree

Showing 4 changed files with 262 additions and 193 deletions.
diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt
@@ -40,11 +40,21 @@ if(SINGULARITY_USE_EOSPAC AND SINGULARITY_TEST_SESAME)
                         PRIVATE Catch2::Catch2 singularity-eos::singularity-eos)
 endif()
 
-if(SINGULARITY_BUILD_CLOSURE AND SINGULARITY_USE_SPINER)
-  add_executable(test_pte test_pte.cpp)
-  target_link_libraries(test_pte PRIVATE Catch2::Catch2
+if(SINGULARITY_BUILD_CLOSURE)
+  if(SINGULARITY_USE_SPINER)
+    add_executable(test_pte test_pte.cpp)
+    target_link_libraries(test_pte PRIVATE Catch2::Catch2
                                          singularity-eos::singularity-eos)
-  add_test(pte test_pte)
+    target_include_directories(test_pte PRIVATE ${CMAKE_CURRENT_SOURCE_DIR})
+    add_test(pte test_pte)
+  endif()
+  if(SINGULARITY_USE_KOKKOS AND SINGULARITY_USE_FORTRAN)
+    add_executable(test_get_sg_eos test_get_sg_eos.cpp)
+    target_link_libraries(test_get_sg_eos PRIVATE Catch2::Catch2
+                                         singularity-eos::singularity-eos)
+    target_include_directories(test_get_sg_eos PRIVATE ${CMAKE_CURRENT_SOURCE_DIR})
+    add_test(get_sg_eos test_get_sg_eos)
+  endif()
 endif()
 
 if(SINGULARITY_TEST_STELLAR_COLLAPSE)

diff --git a/test/pte_test_utils.hpp b/test/pte_test_utils.hpp
@@ -0,0 +1,92 @@
+//------------------------------------------------------------------------------
+// © 2021-2023. 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
+// Nuclear Security Administration. All rights in the program are
+// reserved by Triad National Security, LLC, and the U.S. Department of
+// Energy/National Nuclear Security Administration. The Government is
+// granted for itself and others acting on its behalf a nonexclusive,
+// paid-up, irrevocable worldwide license in this material to reproduce,
+// prepare derivative works, distribute copies to the public, perform
+// publicly and display publicly, and to permit others to do so.
+//------------------------------------------------------------------------------
+
+#ifndef _SINGULARITY_EOS_TEST_PTE_TEST_UTILS_
+#define _SINGULARITY_EOS_TEST_PTE_TEST_UTILS_
+
+#include <stdlib.h>
+
+#include <ports-of-call/portability.hpp>
+#include <ports-of-call/portable_arrays.hpp>
+#include <singularity-eos/eos/eos.hpp>
+
+constexpr int NMAT = 3;
+constexpr int NTRIAL = 100;
+constexpr int NPTS = NTRIAL * NMAT;
+constexpr int HIST_SIZE = 10;
+
+template <typename T>
+class LinearIndexer {
+ public:
+  PORTABLE_FUNCTION LinearIndexer() = default;
+  LinearIndexer(const T &t) : data_(t) {}
+  PORTABLE_INLINE_FUNCTION
+  auto &operator[](const int i) const { return data_(i); }
+
+ private:
+  T data_;
+};
+
+template <typename T>
+class Indexer2D {
+ public:
+  PORTABLE_FUNCTION Indexer2D() = default;
+  PORTABLE_FUNCTION Indexer2D(const int j, const T &t) : j_(j), data_(t) {}
+  Indexer2D(const int j, const T &&t) = delete; // prevents r-value binding
+  PORTABLE_INLINE_FUNCTION
+  auto &operator[](const int i) const { return data_(j_, i); }
+
+ private:
+  const int j_;
+  const T &data_;
+};
+
+template <typename T>
+inline void set_eos(T *eos) {
+  singularity::EOS gr = 
+      singularity::Gruneisen(394000.0, 1.489, 0.0, 0.0, 2.02, 0.47, 8.93, 297.0, 1.0e6, 0.383e7);
+  singularity::EOS dr =
+      singularity::DavisReactants(1.890, 4.115e10, 1.0e6, 297.0, 1.8e5, 4.6, 0.34, 0.56, 0.0,
+                                  0.4265, 0.001074e10);
+  singularity::EOS dp =
+      singularity::DavisProducts(0.798311, 0.58, 1.35, 2.66182, 0.75419, 3.2e10, 0.001072e10, 0.0);
+  eos[0] = gr.GetOnDevice();
+  eos[1] = dr.GetOnDevice();
+  eos[2] = dp.GetOnDevice();
+  return;
+}
+
+template <typename RealIndexer, typename EOSIndexer>
+inline void set_state(RealIndexer &&rho, RealIndexer &&vfrac, RealIndexer &&sie,
+                      RealIndexer &&temp, EOSIndexer &&eos) {
+
+  rho[0] = 8.93;
+  rho[1] = 1.89;
+  rho[2] = 2.5;
+
+  Real vsum = 0.;
+  for (int i = 0; i < NMAT; i++) {
+    temp[i] = 600.0;
+    sie[i] = eos[i].InternalEnergyFromDensityTemperature(rho[i], temp[i]);
+    vfrac[i] = rand() / (1.0 * RAND_MAX);
+    vsum += vfrac[i];
+  }
+
+  for (int i = 0; i < NMAT; i++)
+    vfrac[i] *= 1.0 / vsum;
+
+  return;
+}
+
+#endif // _SINGULARITY_EOS_TEST_PTE_TEST_UTILS_
diff --git a/test/test_get_sg_eos.cpp b/test/test_get_sg_eos.cpp
@@ -0,0 +1,154 @@
+//------------------------------------------------------------------------------
+// © 2021-2023. 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
+// Nuclear Security Administration. All rights in the program are
+// reserved by Triad National Security, LLC, and the U.S. Department of
+// Energy/National Nuclear Security Administration. The Government is
+// granted for itself and others acting on its behalf a nonexclusive,
+// paid-up, irrevocable worldwide license in this material to reproduce,
+// prepare derivative works, distribute copies to the public, perform
+// publicly and display publicly, and to permit others to do so.
+//------------------------------------------------------------------------------
+
+#include <iostream>
+#include <memory>
+#include <stdlib.h>
+
+#include <pte_test_utils.hpp>
+#include <ports-of-call/portability.hpp>
+#include <singularity-eos/closure/mixed_cell_models.hpp>
+#include <singularity-eos/eos/eos.hpp>
+#include <singularity-eos/eos/singularity_eos.hpp>
+
+using namespace singularity;
+
+#ifdef PORTABILITY_STRATEGY_KOKKOS
+// TODO DAH: when the get_sg_eos function is moved out of sg,
+// this function will have to change how the PTE solutions
+// are obtained.
+int run_sg_get_eos_tests() {
+  int nfails = 0;
+  // initialize inputs outputs
+  static constexpr const double ev2k = 1.160451930280894026e4;
+  EOS eoss[NMAT];
+  set_eos(eoss);
+  // set volume fractions to equipartition
+  // all variables needed to check the pte solve
+  Real mfrac[NMAT];
+  Real P_true = 5.e10, T_true = 800.0;
+  Real T_true_ev = T_true / ev2k;
+  Real mass_tot = 0.0, sie_tot_true = 0.0, v_true = 1.0, rho_tot, spvol;
+  Real pmax, bmod, dpde, cv;
+  mfrac[0] = 0.67;
+  mfrac[1] = 0.14;
+  mfrac[2] = 0.19;
+  for (int i = 0; i < NMAT; ++i) {
+    mass_tot += mfrac[i];
+  }
+  rho_tot = mass_tot / v_true;
+  spvol = 1.0 / rho_tot;
+  int cell_offset = 1;
+  int eos_offset[NMAT];
+  for (int m = 0; m < NMAT; ++m) {
+    eos_offset[m] = m + 1;
+  }
+  // do a rho/e(P,T) solve
+  Real vfrac_true[NMAT], ie_true[NMAT];
+  get_sg_eos(NMAT, 1, 1, -1, eos_offset, eoss, &cell_offset, &P_true, &pmax, &v_true,
+             &spvol, &sie_tot_true, &T_true_ev, &bmod, &dpde, &cv, mfrac, vfrac_true,
+             ie_true, nullptr, nullptr, nullptr);
+  Real sie_tot_check = 0.0;
+  for (int m = 0; m < NMAT; ++m) {
+    const Real r_m = mfrac[m] / vfrac_true[m];
+    const Real sie_m = ie_true[m] / mfrac[m];
+    const Real p_eos = eoss[m].PressureFromDensityInternalEnergy(r_m, sie_m);
+    const Real p_resid = std::abs(P_true - p_eos) / P_true;
+    const Real t_eos = eoss[m].TemperatureFromDensityInternalEnergy(r_m, sie_m);
+    const Real t_resid = std::abs(T_true - t_eos) / T_true;
+    // check internal consitency of P-T input function
+    if (t_resid > 1.e-5 || p_resid > 1.e-5) {
+      printf("P-T: t_resid: %e | p_resid: %e\n", t_resid, p_resid);
+      nfails += 1;
+    }
+    sie_tot_check += ie_true[m];
+  }
+  sie_tot_check /= mass_tot;
+  // further check of internal consintency of P-T input function
+  if (std::abs(sie_tot_check - sie_tot_true) / std::abs(sie_tot_true) > 1.e-5) {
+    printf("P-T: sie_tot_true: %e | sie_tot_check: %e\n", sie_tot_true, sie_tot_check);
+    nfails += 1;
+  }
+  // obtain converged and consistent PTE solution
+  // do rho-T input solve
+  Real p_check, vfrac_check[NMAT], ie_check[NMAT];
+  get_sg_eos(NMAT, 1, 1, -3, eos_offset, eoss, &cell_offset, &p_check, &pmax, &v_true,
+             &spvol, &sie_tot_check, &T_true_ev, &bmod, &dpde, &cv, mfrac, vfrac_check,
+             ie_check, nullptr, nullptr, nullptr);
+  // check output pressure and sie, indicate failure if relative err is too large
+  if (std::abs(P_true - p_check) / std::abs(P_true) > 1.e-5 ||
+      std::abs(sie_tot_true - sie_tot_check) / std::abs(sie_tot_true) > 1.e-5) {
+    printf("r-T: p_true: %e | p_check: %e\n", P_true, p_check);
+    printf("r-T: sie_tot_true: %e | sie_tot_check: %e\n", sie_tot_true, sie_tot_check);
+    nfails += 1;
+  }
+  Real max_vfrac_resid = 0.0;
+  Real max_sie_resid = 0.0;
+  for (int m = 0; m < NMAT; ++m) {
+    max_vfrac_resid = std::max(max_vfrac_resid, std::abs(vfrac_true[m] - vfrac_check[m]) /
+                                                    std::abs(vfrac_true[m]));
+    max_sie_resid = std::max(max_sie_resid,
+                             std::abs(ie_true[m] - ie_check[m]) / std::abs(ie_true[m]));
+  }
+  if (max_vfrac_resid > 1.e-5 || max_sie_resid > 1.e-5) {
+    printf("r-T: vr: %e | sr: %e\n", max_vfrac_resid, max_sie_resid);
+    nfails += 1;
+  }
+  // do rho-P input solve
+  Real t_check;
+  get_sg_eos(NMAT, 1, 1, -2, eos_offset, eoss, &cell_offset, &P_true, &pmax, &v_true,
+             &spvol, &sie_tot_check, &t_check, &bmod, &dpde, &cv, mfrac, vfrac_check,
+             ie_check, nullptr, nullptr, nullptr);
+  // check output temperature and sie, indicate failure if relative err is too large
+  if (std::abs(T_true_ev - t_check) / std::abs(T_true_ev) > 1.e-5 ||
+      std::abs(sie_tot_true - sie_tot_check) / std::abs(sie_tot_true) > 1.e-5) {
+    printf("p-T: t_true: %e | t_check: %e\n", T_true_ev, t_check);
+    printf("p-T: sie_tot_true: %e | sie_tot_check: %e\n", sie_tot_true, sie_tot_check);
+    nfails += 1;
+  }
+  max_vfrac_resid = 0.0;
+  max_sie_resid = 0.0;
+  for (int m = 0; m < NMAT; ++m) {
+    max_vfrac_resid = std::max(max_vfrac_resid, std::abs(vfrac_true[m] - vfrac_check[m]) /
+                                                    std::abs(vfrac_true[m]));
+    max_sie_resid = std::max(max_sie_resid,
+                             std::abs(ie_true[m] - ie_check[m]) / std::abs(ie_true[m]));
+  }
+  if (max_vfrac_resid > 1.e-5 || max_sie_resid > 1.e-5) {
+    printf("p-T: vr: %e | sr: %e\n", max_vfrac_resid, max_sie_resid);
+    nfails += 1;
+  }
+  return nfails;
+}
+#endif
+
+int main(int argc, char *argv[]) {
+
+  int nfails_get_sg_eos = 1;
+#ifdef PORTABILITY_STRATEGY_KOKKOS
+  Kokkos::initialize();
+  {
+    // run the various solver tests:
+    // since this uses the sg_get_eos, only run
+    // if kokkos enable since that function requires it
+    // to run the solvers
+    nfails_get_sg_eos = run_sg_get_eos_tests();
+    if (nfails_get_sg_eos > 0) {
+      printf("nfails of fixed T/P solvers = %i\n", nfails_get_sg_eos);
+    }
+  }
+  Kokkos::finalize();
+#endif // PORTABILITY_STRATEGY_KOKKOS
+  return nfails_get_sg_eos;
+}