Skip to content

Commit

Permalink
Merge pull request #466 from lanl/jmm/thread-lambda-through-pte
Browse files Browse the repository at this point in the history 
thread lambda through pte
  • Loading branch information
@Yurlungur
Yurlungur authored Feb 13, 2025
2 parents 129b49f + 59b0034 commit b2c6539
Show file tree
Hide file tree
Showing 5 changed files with 206 additions and 153 deletions.
1 change: 1 addition & 0 deletions CHANGELOG.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@
- [[PR444]](https://github.com/lanl/singularity-eos/pull/444) Add Z split modifier and electron ideal gas EOS

### Fixed (Repair bugs, etc)
- [[PR466]](https://github.com/lanl/singularity-eos/pull/466) Fully thread lambda inputs through the PTE solver
- [[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
Expand Down
85 changes: 56 additions & 29 deletions singularity-eos/closure/mixed_cell_models.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -478,10 +478,11 @@ class PTESolverBase {

} // namespace mix_impl

template <typename EOSIndexer, typename RealIndexer>
template <typename EOSIndexer, typename RealIndexer, typename LambdaIndexer = NullIndexer>
PORTABLE_INLINE_FUNCTION Real ApproxTemperatureFromRhoMatU(
const std::size_t nmat, EOSIndexer &&eos, const Real u_tot, RealIndexer &&rho,
RealIndexer &&vfrac, const Real Tguess = 0.0) {
RealIndexer &&vfrac, const Real Tguess = 0.0,
LambdaIndexer &&lambda = NullIndexer()) {
// should these be passed in?
constexpr Real minimum_temperature = 1.e-9;
constexpr Real maximum_temperature = 1.e9;
Expand All @@ -493,7 +494,8 @@ PORTABLE_INLINE_FUNCTION Real ApproxTemperatureFromRhoMatU(
auto ufunc = [&](const Real T) {
Real usum = 0.0;
for (std::size_t m = 0; m < nmat; m++) {
usum += rho[m] * vfrac[m] * eos[m].InternalEnergyFromDensityTemperature(rho[m], T);
usum += rho[m] * vfrac[m] *
eos[m].InternalEnergyFromDensityTemperature(rho[m], T, lambda[m]);
}
return usum;
};
Expand Down Expand Up @@ -593,10 +595,15 @@ class PTESolverRhoT : public mix_impl::PTESolverBase<EOSIndexer, RealIndexer> {
dedv = AssignIncrement(scratch, nmat);
dpdT = AssignIncrement(scratch, nmat);
vtemp = AssignIncrement(scratch, nmat);
// TODO(JCD): use whatever lambdas are passed in
/*for (int m = 0; m < nmat; m++) {
if (!variadic_utils::is_nullptr(lambda[m])) Cache[m] = lambda[m];
}*/
// JMM: Copy lambda data into cache. Since the cache uses
// pointers, it must be done this way.
for (std::size_t m = 0; m < nmat; ++m) {
if (!variadic_utils::is_nullptr(lambda[m])) {
for (int l = 0; l < eos[m].nlambda(); ++l) {
Cache[m][l] = lambda[m][l];
}
}
}
}

PORTABLE_INLINE_FUNCTION
Expand Down Expand Up @@ -767,14 +774,14 @@ class PTESolverRhoT : public mix_impl::PTESolverBase<EOSIndexer, RealIndexer> {
// ======================================================================
// PT space solver
// ======================================================================
inline int PTESolverPTRequiredScratch(const int nmat) {
inline int PTESolverPTRequiredScratch(const std::size_t nmat) {
constexpr int neq = 2;
return neq * neq // jacobian
+ 4 * neq // dx, residual, and sol_scratch
+ 2 * nmat // all the nmat sized arrays
+ MAX_NUM_LAMBDAS * nmat; // the cache
}
inline size_t PTESolverPTRequiredScratchInBytes(const int nmat) {
inline size_t PTESolverPTRequiredScratchInBytes(const std::size_t nmat) {
return PTESolverPTRequiredScratch(nmat) * sizeof(Real);
}
template <typename EOSIndexer, typename RealIndexer, typename LambdaIndexer>
Expand Down Expand Up @@ -812,17 +819,22 @@ class PTESolverPT : public mix_impl::PTESolverBase<EOSIndexer, RealIndexer> {
// template the ctor to get type deduction/universal references prior to c++17
template <typename EOS_t, typename Real_t, typename Lambda_t>
PORTABLE_INLINE_FUNCTION
PTESolverPT(const int nmat, EOS_t &&eos, const Real vfrac_tot, const Real sie_tot,
Real_t &&rho, Real_t &&vfrac, Real_t &&sie, Real_t &&temp, Real_t &&press,
Lambda_t &&lambda, Real *scratch, const Real Tnorm = 0.0,
const MixParams &params = MixParams())
PTESolverPT(const std::size_t nmat, EOS_t &&eos, const Real vfrac_tot,
const Real sie_tot, Real_t &&rho, Real_t &&vfrac, Real_t &&sie,
Real_t &&temp, Real_t &&press, Lambda_t &&lambda, Real *scratch,
const Real Tnorm = 0.0, const MixParams &params = MixParams())
: mix_impl::PTESolverBase<EOSIndexer, RealIndexer>(nmat, 2, eos, vfrac_tot, sie_tot,
rho, vfrac, sie, temp, press,
scratch, Tnorm, params) {
// TODO(JCD): use whatever lambdas are passed in
/*for (int m = 0; m < nmat; m++) {
if (!variadic_utils::is_nullptr(lambda[m])) Cache[m] = lambda[m];
}*/
// JMM: Copy lambda data into cache. Since the cache uses
// pointers, it must be done this way.
for (std::size_t m = 0; m < nmat; ++m) {
if (!variadic_utils::is_nullptr(lambda[m])) {
for (int l = 0; l < eos[m].nlambda(); ++l) {
Cache[m][l] = lambda[m][l];
}
}
}
}

PORTABLE_INLINE_FUNCTION
Expand Down Expand Up @@ -1066,10 +1078,15 @@ class PTESolverFixedT : public mix_impl::PTESolverBase<EOSIndexer, RealIndexer>
Tequil = T_true;
Ttemp = T_true;
Tnorm = 1.0;
// TODO(JCD): use whatever lambdas are passed in
/*for (int m = 0; m < nmat; m++) {
if (variadic_utils::is_nullptr(lambda[m])) Cache[m] = lambda[m];
}*/
// JMM: Copy lambda data into cache. Since the cache uses
// pointers, it must be done this way.
for (std::size_t m = 0; m < nmat; ++m) {
if (!variadic_utils::is_nullptr(lambda[m])) {
for (int l = 0; l < eos[m].nlambda(); ++l) {
Cache[m][l] = lambda[m][l];
}
}
}
}

PORTABLE_INLINE_FUNCTION
Expand Down Expand Up @@ -1275,10 +1292,15 @@ class PTESolverFixedP : public mix_impl::PTESolverBase<EOSIndexer, RealIndexer>
dpdT = AssignIncrement(scratch, nmat);
vtemp = AssignIncrement(scratch, nmat);
Pequil = P;
// TODO(JCD): use whatever lambdas are passed in
/*for (int m = 0; m < nmat; m++) {
if (variadic_utils::is_nullptr(lambda[m])) Cache[m] = lambda[m];
}*/
// JMM: Copy lambda data into cache. Since the cache uses
// pointers, it must be done this way.
for (std::size_t m = 0; m < nmat; ++m) {
if (!variadic_utils::is_nullptr(lambda[m])) {
for (int l = 0; l < eos[m].nlambda(); ++l) {
Cache[m][l] = lambda[m][l];
}
}
}
}

PORTABLE_INLINE_FUNCTION
Expand Down Expand Up @@ -1514,10 +1536,15 @@ class PTESolverRhoU : public mix_impl::PTESolverBase<EOSIndexer, RealIndexer> {
dtde = AssignIncrement(scratch, nmat);
vtemp = AssignIncrement(scratch, nmat);
utemp = AssignIncrement(scratch, nmat);
// TODO(JCD): use whatever lambdas are passed in
/*for (int m = 0; m < nmat; m++) {
if (variadic_utils::is_nullptr(lambda[m])) Cache[m] = lambda[m];
}*/
// JMM: Copy lambda data into cache. Since the cache uses
// pointers, it must be done this way.
for (std::size_t m = 0; m < nmat; ++m) {
if (!variadic_utils::is_nullptr(lambda[m])) {
for (int l = 0; l < eos[m].nlambda(); ++l) {
Cache[m][l] = lambda[m][l];
}
}
}
}

PORTABLE_INLINE_FUNCTION
Expand Down
6 changes: 3 additions & 3 deletions singularity-eos/eos/get_sg_eos_rho_p.cpp
View file
Original file line number Diff line number Diff line change
@@ -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
Expand Down Expand Up @@ -61,9 +61,9 @@ void get_sg_eos_rho_p(const char *name, int ncell, indirection_v &offsets_v,
Real *ptemp_pte = &temp_pte(tid, 0);
Real *ppress_pte = &press_pte(tid, 0);
Real *pscratch = &solver_scratch(tid, 0);
PTESolverFixedP<singularity::EOSAccessor_, Real *, Real *> method(
PTESolverFixedP<singularity::EOSAccessor_, Real *, Real **> method(
npte, eos_inx, vfrac_sum, press_pte(tid, 0), prho_pte, pvfrac_pte, psie_pte,
ptemp_pte, ppress_pte, cache[0], pscratch);
ptemp_pte, ppress_pte, cache, pscratch);
auto status = PTESolver(method);
pte_converged = status.converged;
// calculate total sie
Expand Down
13 changes: 8 additions & 5 deletions test/test_pte.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -46,10 +46,11 @@ using atomic_view = Kokkos::MemoryTraits<Kokkos::Atomic>;
#endif

template <template <typename... Types> class Method_t>
auto TestPTE(const std::string name, const std::size_t nscratch_vars) {
auto TestPTE(const std::string name, const std::size_t nscratch_vars,
std::size_t &nsuccess) {
constexpr Real EPS = 1e-5;
Real time;
std::size_t nsuccess = 0;
nsuccess = 0;

// EOS
#ifdef PORTABILITY_STRATEGY_KOKKOS
Expand Down Expand Up @@ -231,7 +232,7 @@ auto TestPTE(const std::string name, const std::size_t nscratch_vars) {
}
std::cout << std::endl;

return std::make_pair(nsuccess, rho_d);
return rho_d;
}

int main(int argc, char *argv[]) {
Expand All @@ -243,13 +244,15 @@ int main(int argc, char *argv[]) {
srand(time(NULL));

// scratch required for PTE solver
std::size_t ns_rt;
auto nscratch_vars_rt = PTESolverRhoTRequiredScratch(NMAT);
auto [ns_rt, rho_rt] = TestPTE<PTESolverRhoT>("PTESolverRhoT", nscratch_vars_rt);
auto rho_rt = TestPTE<PTESolverRhoT>("PTESolverRhoT", nscratch_vars_rt, ns_rt);
nsuccess += ns_rt;

// // scratch required for PTE solver
std::size_t ns_pt;
auto nscratch_vars_pt = PTESolverPTRequiredScratch(NMAT);
auto [ns_pt, rho_pt] = TestPTE<PTESolverPT>("PTESolverPT", nscratch_vars_pt);
auto rho_pt = TestPTE<PTESolverPT>("PTESolverPT", nscratch_vars_pt, ns_pt);
nsuccess += ns_pt;

int nmatch = 0;
Expand Down
Loading

0 comments on commit b2c6539

Please sign in to comment.