1) 3 Tests for testing PTESolverRhoT versus Ann's Flag PTsolv 2) Implementing Carl Grief's mass fraction update model.

test/CMakeLists.txt

@@ -99,6 +99,20 @@ if(SINGULARITY_BUILD_CLOSURE)
     target_include_directories(test_pte PRIVATE ${CMAKE_CURRENT_SOURCE_DIR})
     add_test(pte test_pte)
   endif()
+  if(SINGULARITY_USE_SPINER)
+    add_executable(test_pte_2phase test_pte_2phase.cpp)
+    target_link_libraries(test_pte_2phase PRIVATE Catch2::Catch2
+                                         singularity-eos::singularity-eos)
+    target_include_directories(test_pte_2phase PRIVATE ${CMAKE_CURRENT_SOURCE_DIR})
+    add_test(pte_2phase test_pte_2phase)
+  endif()
+   if(SINGULARITY_USE_SPINER)
+    add_executable(test_pte_3phase test_pte_3phase.cpp)
+    target_link_libraries(test_pte_3phase PRIVATE Catch2::Catch2
+                                         singularity-eos::singularity-eos)
+    target_include_directories(test_pte_3phase PRIVATE ${CMAKE_CURRENT_SOURCE_DIR})
+    add_test(pte_3phase test_pte_3phase)
+  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

test/pte_longtest_2phaseVinetSn.hpp

@@ -0,0 +1,135 @@
+//------------------------------------------------------------------------------
+// © 2021-2024. 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_2PHASEVINETSN_
+#define _SINGULARITY_EOS_TEST_PTE_TEST_2PHASEVINETSN_
+
+#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 = 2;
+constexpr int NTRIAL = 20;
+constexpr int NPTS = NTRIAL * NMAT;
+constexpr int HIST_SIZE = 10;
+
+constexpr Real in_rho_tot[NTRIAL] = {
+    7.278,      7.46221983, 7.65110201, 7.84476448, 8.04332821, 8.24691719, 8.45565857,
+    8.66968272, 8.88912327, 9.11411728, 9.34480523, 9.58133117, 9.8238428,  10.0724915,
+    10.3274327, 10.5888253, 10.8568327, 11.1316222, 11.4133653, 11.702238};
+constexpr Real in_sie_tot[NTRIAL] = {
+    8.41323509e8, 8.62513154e8,         9.28665206e8,  1.04377144e9,  1.21199344e9,
+    1.43767572e9, 1.72535639e9,         2.07977629e9,  2.50588681e9,  3.00885704e9,
+    3.59408011e9, 4.2671792100000005e9, 5.03401308e9,  5.90068122e9,  6.87352878e9,
+    7.95915117e9, 9.16439846e9,         1.04963795e10, 1.19624656e10, 1.35702945e10};
+constexpr Real in_lambda0 = 0.500480901;
+constexpr Real in_lambda1 = 0.499519099;
+
+constexpr Real trial_vfrac0 = 47.0 / 100.0;
+constexpr Real trial_vfrac1 = 53.0 / 100.0;
+
+constexpr Real out_press[NTRIAL] = {
+    -3.29164604e-6,        1.284232715e10,        2.753234765e10,
+    4.423652945000001e10,  6.313670939999999e10,  8.443021825e10,
+    1.083303075e11,        1.3506674800000002e11, 1.64886557e11,
+    1.9805482549999997e11, 2.3485562650000003e11, 2.755929975e11,
+    3.20591986e11,         3.70199756e11,         4.247867485e11,
+    4.84747905e11,         5.505039405000001e11,  6.225026735e11,
+    7.012204135e11,        7.871634035e11};
+constexpr Real out_temp[NTRIAL] = {298.,        317.657834,         338.05429449999997,
+                                   359.18724,   381.0524905,        403.64396550000004,
+                                   426.953795,  450.97240750000003, 475.6886215,
+                                   501.0897275, 527.1615735,        553.8886454999999,
+                                   581.254151,  609.2401025,        637.8273995,
+                                   666.995915,  696.7245780000001,  726.9914615,
+                                   757.7738645, 789.048397};
+constexpr Real out_rho0[NTRIAL] = {
+    7.285,      7.44383086, 7.61073863, 7.78523193, 7.9669718,  8.15572753, 8.35135008,
+    8.55375288, 8.76289814, 8.97878664, 9.20145031, 9.43094663, 9.6673544,  9.91077059,
+    10.1613078, 10.4190924, 10.6842631, 10.9569698, 11.2373726, 11.5256413};
+constexpr Real out_rho1[NTRIAL] = {
+    7.271,      7.48073556, 7.69197476, 7.90533181, 8.12131372, 8.34035068, 8.56281419,
+    8.78903065, 9.01929178, 9.25386249, 9.49298694, 9.73689325, 9.98579716, 10.2399049,
+    10.4994157, 10.7645232, 11.0354173, 11.3122855, 11.5953136, 11.8846866};
+constexpr Real out_vfrac0[NTRIAL] = {0.5,         0.501717271, 0.50313519,  0.504308007,
+                                     0.50527757,  0.506076807, 0.506731916, 0.507263967,
+                                     0.507690077, 0.508024281, 0.508278194, 0.508461499,
+                                     0.508582337, 0.508647597, 0.508663148, 0.50863402,
+                                     0.508564548, 0.508458492, 0.508319122, 0.508149303};
+constexpr Real out_vfrac1[NTRIAL] = {0.5,         0.498282729, 0.49686481,  0.495691993,
+                                     0.49472243,  0.493923193, 0.493268084, 0.492736033,
+                                     0.492309923, 0.491975719, 0.491721806, 0.491538501,
+                                     0.491417663, 0.491352403, 0.491336852, 0.49136598,
+                                     0.491435452, 0.491541508, 0.491680878, 0.491850697};
+
+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) {
+  constexpr Real d2to40[39] = {0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+                               0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+                               0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.};
+
+  singularity::EOS Snbeta =
+      singularity::Vinet(7.285, 298.0, 0.529e12, 5.3345, 0.000072977, 0.2149e07, 0.658e09,
+                         0.4419e07, d2to40);
+  singularity::EOS Sngamma =
+      singularity::Vinet(7.271, 298.0, 0.3878e12, 6.0532, 0.0001085405, 0.2161e07,
+                         1.025e09, 0.5051e07, d2to40);
+  eos[0] = Snbeta.GetOnDevice();
+  eos[1] = Sngamma.GetOnDevice();
+  return;
+}
+
+template <typename RealIndexer, typename EOSIndexer>
+inline void set_state(RealIndexer &&rho, RealIndexer &&vfrac, RealIndexer &&sie,
+                      RealIndexer &&temp, EOSIndexer &&eos) {
+
+  vfrac[0] = 47.0 / 100.0;
+  vfrac[1] = 53.0 / 100.0;
+  rho[0] = 7.278 * 0.500480901 / vfrac[0];
+  rho[1] = 7.278 * 0.499519099 / vfrac[1];
+  sie[0] = 8.41323509e08;
+  sie[1] = 8.41323509e08;
+
+  return;
+}
+
+#endif // _SINGULARITY_EOS_TEST_PTE_TEST_2PHASEVINETSN_

test/pte_test_2phaseVinetSn.hpp

@@ -0,0 +1,91 @@
+//------------------------------------------------------------------------------
+// © 2021-2024. 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_2PHASEVINETSN_
+#define _SINGULARITY_EOS_TEST_PTE_TEST_2PHASEVINETSN_
+
+#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 = 2;
+constexpr int NTRIAL = 5;
+constexpr int NPTS = NTRIAL * NMAT;
+constexpr int HIST_SIZE = 10;
+
+constexpr Real in_rho_tot[5] = {7.27800000, 7.27981950, 7.28163945, 7.28345986,
+                                7.28528073};
+constexpr Real in_sie_tot[5] = {8.41323509e08, 8.41325565e08, 8.41331734e08,
+                                8.41342022e08, 8.41356432e08};
+constexpr Real in_lambda0 = 0.500480901;
+constexpr Real in_lambda1 = 0.499519099;
+
+constexpr Real trial_vfrac0 = 47.0 / 100.0;
+constexpr Real trial_vfrac1 = 53.0 / 100.0;
+
+constexpr Real out_press[5] = {-3.29164604e-6, 1.19722694e8, 2.3968114450000003e8,
+                               3.598077865e8, 4.80104422e8};
+constexpr Real out_temp[5] = {298., 298.192952, 298.38594450000005, 298.5790125,
+                              298.7721535};
+constexpr Real out_rho0[5] = {7.28500000, 7.28653963, 7.28808705, 7.28963516, 7.29118416};
+constexpr Real out_rho1[5] = {7.27100000, 7.27309885, 7.27519088, 7.27728316, 7.27937551};
+constexpr Real out_vfrac0[5] = {0.5, 0.500019325, 0.500038138, 0.500056927, 0.500075679};
+constexpr Real out_vfrac1[5] = {0.5, 0.499980675, 0.499961862, 0.499943073, 0.499924321};
+
+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) {
+  constexpr Real d2to40[39] = {0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+                               0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.,
+                               0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.};
+
+  singularity::EOS Snbeta =
+      singularity::Vinet(7.285, 298.0, 0.529e12, 5.3345, 0.000072977, 0.2149e07, 0.658e09,
+                         0.4419e07, d2to40);
+  singularity::EOS Sngamma =
+      singularity::Vinet(7.271, 298.0, 0.3878e12, 6.0532, 0.0001085405, 0.2161e07,
+                         1.025e09, 0.5051e07, d2to40);
+  eos[0] = Snbeta.GetOnDevice();
+  eos[1] = Sngamma.GetOnDevice();
+  return;
+}
+
+#endif // _SINGULARITY_EOS_TEST_PTE_TEST_2PHASEVINETSN_

test/pte_test_5phaseSesameSn.hpp

@@ -0,0 +1,113 @@
+//------------------------------------------------------------------------------
+// © 2021-2024. 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_5PHASESESAMESN_
+#define _SINGULARITY_EOS_TEST_PTE_TEST_5PHASESESAMESN_
+
+#include <stdlib.h>
+
+#include <unistd.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 = 5;
+constexpr int NPTS = NTRIAL * NMAT;
+constexpr int HIST_SIZE = 30;
+
+constexpr Real in_rho_tot[5] = {8.41382588, 8.41550873, 8.41719191, 8.41887543,
+                                8.42055928};
+constexpr Real in_sie_tot[5] = {1.67602914e09, 1.67876139e09, 1.68149578e09,
+                                1.68423198e09, 1.68696932e09};
+constexpr Real in_lambda0[5] = {0.663852654, 0.660579238, 0.656880511, 0.652218675,
+                                0.646232106};
+constexpr Real in_lambda1[5] = {0.336093276, 0.339262144, 0.342703726, 0.346800126,
+                                0.352135167};
+constexpr Real in_lambda2[5] = {0.000054070, 0.000158618, 0.000415763, 0.000982199,
+                                0.001632727};
+
+constexpr Real trial_vfrac0 = 650.0 / 1000.0;
+constexpr Real trial_vfrac1 = 349.5 / 1000.0;
+constexpr Real trial_vfrac2 = 0.5 / 1000.0;
+
+constexpr Real out_press[5] = {1.0946046e11, 1.0959183e11, 1.0971241e11, 1.0980857e11,
+                               1.0987166e11};
+constexpr Real out_temp[5] = {438.96969, 438.98046, 438.95507, 438.84886, 438.64322};
+constexpr Real out_rho0[5] = {8.35319483, 8.35425038, 8.35523024, 8.35603912, 8.35661359};
+constexpr Real out_rho1[5] = {8.53618806, 8.53734120, 8.53841145, 8.53929455, 8.53992118};
+constexpr Real out_rho2[5] = {8.53854993, 8.53970495, 8.54077710, 8.54166211, 8.54229063};
+constexpr Real out_sie0[5] = {1.54379257e09, 1.54520613e09, 1.54644403e09, 1.54728448e09,
+                              1.54760327e09};
+constexpr Real out_sie1[5] = {1.93716882e09, 1.93864981e09, 1.93994981e09, 1.94084038e09,
+                              1.94119298e09};
+constexpr Real out_sie2[5] = {2.01473728e09, 2.01621655e09, 2.01751522e09, 2.01840528e09,
+                              2.01875845e09};
+
+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) {
+
+  symlink("/projects/shavano/dev/aematts/SEOS/PTEsolver_test_cases/sn2162-v01.bin",
+          "sesameu");
+
+  int SnbetaID = 2102;
+  int SngammaID = 2103;
+  int SndeltaID = 2104;
+  int SnhcpID = 2105;
+  int SnliquidID = 2106;
+
+  bool invert_at_setup = true;
+
+  singularity::EOS Snbeta = singularity::EOSPAC(SnbetaID);
+  singularity::EOS Sngamma = singularity::EOSPAC(SngammaID);
+  //  singularity::EOS Sndelta = singularity::EOSPAC(SndeltaID);
+  singularity::EOS Snhcp = singularity::EOSPAC(SnhcpID);
+  //  singularity::EOS Snliquid = singularity::EOSPAC(SnliquidID);
+
+  eos[0] = Snbeta.GetOnDevice();
+  eos[1] = Sngamma.GetOnDevice();
+  //  eos[x] = Sndelta.GetOnDevice();
+  eos[2] = Snhcp.GetOnDevice();
+  //  eos[x] = Snliquid.GetOnDevice();
+  return;
+}
+
+#endif // _SINGULARITY_EOS_TEST_PTE_TEST_5PHASESESAMESN_