diff --git a/.clang-tidy b/.clang-tidy
index dae91c5e27..c2a01354d5 100644
--- a/.clang-tidy
+++ b/.clang-tidy
@@ -21,6 +21,7 @@ Checks: >
        -misc-use-internal-linkage,
      modernize-*,
        -modernize-avoid-c-arrays,
+       -modernize-use-designated-initializers,
        -modernize-use-trailing-return-type,
      performance-*,
        -performance-avoid-endl,
diff --git a/.github/workflows/castro-development.yml b/.github/workflows/castro-development.yml
index 6df70b8a35..7c595834af 100644
--- a/.github/workflows/castro-development.yml
+++ b/.github/workflows/castro-development.yml
@@ -47,7 +47,7 @@ jobs:
         export MICROPHYSICS_HOME=${PWD}
 
         cd Castro/Exec/science/flame_wave/
-        make -j 4 CCACHE=ccache USE_MPI=FALSE
+        make -j 4 CXXSTD=c++20 CCACHE=ccache USE_MPI=FALSE
 
         ccache -s
         du -hs ~/.cache/ccache
@@ -64,7 +64,7 @@ jobs:
         export MICROPHYSICS_HOME=${PWD}
 
         cd Castro/Exec/science/subchandra/
-        make -j 4 CCACHE=ccache USE_MPI=FALSE NETWORK_DIR=he-burn/he-burn-19am
+        make -j 4 CXXSTD=c++20 CCACHE=ccache USE_MPI=FALSE NETWORK_DIR=he-burn/he-burn-19am
 
         ccache -s
         du -hs ~/.cache/ccache
diff --git a/.github/workflows/castro.yml b/.github/workflows/castro.yml
index 2d1b1776cc..65580766f7 100644
--- a/.github/workflows/castro.yml
+++ b/.github/workflows/castro.yml
@@ -51,7 +51,7 @@ jobs:
         export MICROPHYSICS_HOME=${PWD}
 
         cd Castro/Exec/science/flame_wave/
-        make -j 4 CCACHE=ccache USE_MPI=FALSE
+        make -j 4 CXXSTD=c++20 CCACHE=ccache USE_MPI=FALSE
 
         ccache -s
         du -hs ~/.cache/ccache
@@ -68,7 +68,7 @@ jobs:
         export MICROPHYSICS_HOME=${PWD}
 
         cd Castro/Exec/science/subchandra/
-        make -j 4 CCACHE=ccache USE_MPI=FALSE NETWORK_DIR=he-burn/he-burn-19am
+        make -j 4 CXXSTD=c++20 CCACHE=ccache USE_MPI=FALSE NETWORK_DIR=he-burn/he-burn-19am
 
         ccache -s
         du -hs ~/.cache/ccache
diff --git a/.github/workflows/cuda.yml b/.github/workflows/cuda.yml
index 1a5054065d..3b73be7da8 100644
--- a/.github/workflows/cuda.yml
+++ b/.github/workflows/cuda.yml
@@ -8,7 +8,7 @@ concurrency:
 
 jobs:
   cuda-compile:
-    name: CUDA@11.7 GCC
+    name: CUDA@12.9 GCC
     runs-on: ubuntu-22.04
     steps:
       - uses: actions/checkout@v6
@@ -26,7 +26,7 @@ jobs:
           cd ../..
 
       - name: Dependencies
-        run: .github/workflows/dependencies/dependencies_nvcc.sh 11.7
+        run: .github/workflows/dependencies/dependencies_nvcc.sh 12.9
 
       - name: compile test_react (aprox13)
         run: |
diff --git a/Docs/source/getting_started.rst b/Docs/source/getting_started.rst
index b5912ca912..663960de7d 100644
--- a/Docs/source/getting_started.rst
+++ b/Docs/source/getting_started.rst
@@ -14,7 +14,7 @@ Microphysics requires
 
 optional dependencies are:
 
-* CUDA (≥ 11)
+* CUDA (≥ 12)
 * ROCm (≥ 6.3.1 --- earlier versions have register allocation bugs)
 
 Usage Modes
diff --git a/EOS/breakout/actual_eos.H b/EOS/breakout/actual_eos.H
index cef9e7629e..e7293ca407 100644
--- a/EOS/breakout/actual_eos.H
+++ b/EOS/breakout/actual_eos.H
@@ -67,14 +67,14 @@ void actual_eos (I input, T& state)
         // dens, temp and xmass are inputs
         amrex::Real cv = R / (state.mu * (gamma_const - 1.0_rt));
         amrex::Real e = cv * state.T;
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             state.cv = cv;
             state.e = e;
         }
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             state.p = (gamma_const - 1.0_rt) * state.rho * e;
         }
-        if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+        if constexpr (has_pressure<T> && has_energy<T>) {
             state.gam1 = gamma_const;
         }
 
@@ -102,10 +102,10 @@ void actual_eos (I input, T& state)
     {
         // dens, pres, and xmass are inputs
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             amrex::Real poverrho = state.p / state.rho;
             state.T = poverrho * state.mu * (1.0_rt / R);
-            if constexpr (has_energy<T>::value) {
+            if constexpr (has_energy<T>) {
                 state.e = poverrho * (1.0_rt / (gamma_const - 1.0_rt));
                 state.gam1 = gamma_const;
             }
@@ -117,11 +117,11 @@ void actual_eos (I input, T& state)
     {
         // dens, energy, and xmass are inputs
 
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             amrex::Real poverrho = (gamma_const - 1.0_rt) * state.e;
             state.T = poverrho * state.mu * (1.0_rt / R);
 
-            if constexpr (has_pressure<T>::value) {
+            if constexpr (has_pressure<T>) {
                 state.p = poverrho * state.rho;
                 state.gam1 = gamma_const;
 
@@ -135,7 +135,7 @@ void actual_eos (I input, T& state)
 
         // Try to avoid the expensive log function.  Since we don't need entropy
         // in hydro solver, set it to an invalid but "nice" value for the plotfile.
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             state.s = 1.0_rt;
         }
 
diff --git a/EOS/eos_composition.H b/EOS/eos_composition.H
index 9e4d909688..3a808ae8a0 100644
--- a/EOS/eos_composition.H
+++ b/EOS/eos_composition.H
@@ -85,7 +85,7 @@ eos_xderivs_t composition_derivatives (const T& state)
   eos_xderivs_t state_xderivs;
 
   for (int n = 0; n < NumSpec; n++) {
-      if constexpr (has_dpdA<T>::value && has_dpdZ<T>::value) {
+      if constexpr (has_dpdA<T> && has_dpdZ<T>) {
           state_xderivs.dpdX[n] =
               state.dpdA * (state.abar * aion_inv[n]) * (aion[n] - state.abar) +
               state.dpdZ * (state.abar * aion_inv[n]) * (zion[n] - state.zbar);
@@ -94,7 +94,7 @@ eos_xderivs_t composition_derivatives (const T& state)
           state_xderivs.dpdX[n] = 0.0_rt;
       }
 
-      if constexpr (has_dedA<T>::value && has_dedZ<T>::value) {
+      if constexpr (has_dedA<T> && has_dedZ<T>) {
           state_xderivs.dedX[n] =
               state.dedA * (state.abar * aion_inv[n]) * (aion[n] - state.abar) +
               state.dedZ * (state.abar * aion_inv[n]) * (zion[n] - state.zbar);
@@ -103,7 +103,7 @@ eos_xderivs_t composition_derivatives (const T& state)
           state_xderivs.dedX[n] = 0.0_rt;
       }
 
-      if constexpr (has_pressure<T>::value) {
+      if constexpr (has_pressure<T>) {
           if (state.dpdr != 0.0) {
 
               state_xderivs.dhdX[n] = state_xderivs.dedX[n]
diff --git a/EOS/gamma_law/actual_eos.H b/EOS/gamma_law/actual_eos.H
index 515ad674ad..902e88beae 100644
--- a/EOS/gamma_law/actual_eos.H
+++ b/EOS/gamma_law/actual_eos.H
@@ -54,7 +54,7 @@ void actual_eos (I input, T& state)
     // Get the mass of a nucleon from m_u.
     const amrex::Real m_nucleon = C::m_u;
 
-    if constexpr (has_xn<T>::value) {
+    if constexpr (has_xn<T>) {
         if (eos_rp::eos_assume_neutral) {
             state.mu = state.abar;
         } else {
@@ -86,7 +86,7 @@ void actual_eos (I input, T& state)
         // Solve for the temperature:
         // h = e + p/rho = (p/rho)*[1 + 1/(gamma-1)] = (p/rho)*gamma/(gamma-1)
 
-          if constexpr (has_enthalpy<T>::value) {
+          if constexpr (has_enthalpy<T>) {
               state.T = (state.h * state.mu * m_nucleon / C::k_B)*(eos_rp::eos_gamma - 1.0)/eos_rp::eos_gamma;
           }
 
@@ -99,7 +99,7 @@ void actual_eos (I input, T& state)
           // Solve for the density:
           // p = rho k T / (mu m_nucleon)
 
-          if constexpr (has_pressure<T>::value) {
+          if constexpr (has_pressure<T>) {
               state.rho =  state.p * state.mu * m_nucleon / (C::k_B * state.T);
           }
 
@@ -112,7 +112,7 @@ void actual_eos (I input, T& state)
           // Solve for the temperature:
           // p = rho k T / (mu m_nucleon)
 
-          if constexpr (has_pressure<T>::value) {
+          if constexpr (has_pressure<T>) {
               state.T = state.p * state.mu * m_nucleon / (C::k_B * state.rho);
           }
 
@@ -125,7 +125,7 @@ void actual_eos (I input, T& state)
           // Solve for the temperature
           // e = k T / [(mu m_nucleon)*(gamma-1)]
 
-          if constexpr (has_energy<T>::value) {
+          if constexpr (has_energy<T>) {
               state.T = state.e * state.mu * m_nucleon * (eos_rp::eos_gamma - 1.0) / C::k_B;
           }
 
@@ -139,7 +139,7 @@ void actual_eos (I input, T& state)
           // Invert Sackur-Tetrode eqn (below) using
           // rho = p mu m_nucleon / (k T)
 
-          if constexpr (has_pressure<T>::value && has_entropy<T>::value) {
+          if constexpr (has_pressure<T> && has_entropy<T>) {
               state.T = std::pow(state.p, 2.0/5.0) *
                         std::pow(2.0 * M_PI * C::hbar * C::hbar / (state.mu * m_nucleon), 3.0/5.0) *
                         std::exp(2.0 * state.mu * m_nucleon * state.s / (5.0 * C::k_B) - 1.0) / C::k_B;
@@ -148,7 +148,7 @@ void actual_eos (I input, T& state)
           // Solve for the density
           // rho = p mu m_nucleon / (k T)
 
-          if constexpr (has_pressure<T>::value) {
+          if constexpr (has_pressure<T>) {
               state.rho = state.p * state.mu * m_nucleon / (C::k_B * state.T);
           }
 
@@ -160,7 +160,7 @@ void actual_eos (I input, T& state)
 
           // Solve for temperature and density
 
-          if constexpr (has_pressure<T>::value && has_enthalpy<T>::value) {
+          if constexpr (has_pressure<T> && has_enthalpy<T>) {
               state.rho = state.p / state.h * eos_rp::eos_gamma / (eos_rp::eos_gamma - 1.0);
               state.T = state.p * state.mu * m_nucleon / (C::k_B * state.rho);
           }
@@ -199,21 +199,21 @@ void actual_eos (I input, T& state)
     // specific internal energy using the gamma-law EOS relation.
     amrex::Real pressure = state.rho * state.T * C::k_B / (state.mu * m_nucleon);
     amrex::Real energy = pressure / (eos_rp::eos_gamma - 1.0) * rhoinv;
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p = pressure;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e = energy;
     }
 
     // enthalpy is h = e + p/rho
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.h = energy + pressure * rhoinv;
     }
 
     // entropy (per gram) of an ideal monoatomic gas (the Sackur-Tetrode equation)
     // NOTE: this expression is only valid for gamma = 5/3.
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         const amrex::Real fac = 1.0 / std::pow(2.0 * M_PI * C::hbar * C::hbar, 1.5);
 
         state.s = (C::k_B / (state.mu * m_nucleon)) *
@@ -222,39 +222,39 @@ void actual_eos (I input, T& state)
     }
 
     // Compute the thermodynamic derivatives and specific heats
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpdT = state.p * Tinv;
         state.dpdr = state.p * rhoinv;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedT = state.e * Tinv;
         state.dedr = 0.0;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.dsdT = 1.5 * (C::k_B / (state.mu * m_nucleon)) * Tinv;
         state.dsdr = - (C::k_B / (state.mu * m_nucleon)) * rhoinv;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.dhdT = state.dedT + state.dpdT * rhoinv;
         state.dhdr = 0.0;
     }
 
-    if constexpr (has_xne_xnp<T>::value) {
+    if constexpr (has_xne_xnp<T>) {
         state.xne = 0.0;
         state.xnp = 0.0;
     }
-    if constexpr (has_eta<T>::value) {
+    if constexpr (has_eta<T>) {
         state.eta = 0.0;
     }
-    if constexpr (has_pele_ppos<T>::value) {
+    if constexpr (has_pele_ppos<T>) {
         state.pele = 0.0;
         state.ppos = 0.0;
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.cv = state.dedT;
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             state.cp = eos_rp::eos_gamma * state.cv;
 
             state.gam1 = eos_rp::eos_gamma;
@@ -264,31 +264,31 @@ void actual_eos (I input, T& state)
 
             // sound speed
             state.cs = std::sqrt(eos_rp::eos_gamma * state.p * rhoinv);
-            if constexpr (has_G<T>::value) {
+            if constexpr (has_G<T>) {
                 state.G = 0.5 * (1.0 + eos_rp::eos_gamma);
             }
         }
     }
 
-    if constexpr (has_dpdA<T>::value) {
+    if constexpr (has_dpdA<T>) {
         state.dpdA = - state.p * (1.0 / state.abar);
     }
-    if constexpr (has_dedA<T>::value) {
+    if constexpr (has_dedA<T>) {
         state.dedA = - state.e * (1.0 / state.abar);
     }
 
     if (eos_rp::eos_assume_neutral) {
-        if constexpr (has_dpdZ<T>::value) {
+        if constexpr (has_dpdZ<T>) {
             state.dpdZ = 0.0;
         }
-        if constexpr (has_dedZ<T>::value) {
+        if constexpr (has_dedZ<T>) {
             state.dedZ = 0.0;
         }
     } else {
-        if constexpr (has_dpdZ<T>::value) {
+        if constexpr (has_dpdZ<T>) {
             state.dpdZ = state.p * (1.0 / (1.0 + state.zbar));
         }
-        if constexpr (has_dedZ<T>::value) {
+        if constexpr (has_dedZ<T>) {
             state.dedZ = state.e * (1.0/(1.0 + state.zbar));
         }
     }
diff --git a/EOS/helmholtz/actual_eos.H b/EOS/helmholtz/actual_eos.H
index 4ae178b0b1..118c029896 100644
--- a/EOS/helmholtz/actual_eos.H
+++ b/EOS/helmholtz/actual_eos.H
@@ -322,7 +322,7 @@ void apply_electrons (T& state)
 
     // pressure derivative with density
     [[maybe_unused]] amrex::Real dpepdd = 0.0e0_rt;
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         amrex::constexpr_for<0, 16>([&] (auto K) {
             constexpr int k = K;
             constexpr auto m = electron_table_indexing::map[k];
@@ -334,7 +334,7 @@ void apply_electrons (T& state)
 
     // electron chemical potential etaele
     [[maybe_unused]] amrex::Real etaele{};
-    if constexpr (has_eta<T>::value) {
+    if constexpr (has_eta<T>) {
         amrex::constexpr_for<0, 16>([&] (auto K) {
             constexpr int k = K;
             constexpr auto m = electron_table_indexing::map[k];
@@ -345,7 +345,7 @@ void apply_electrons (T& state)
 
     // electron + positron number densities
     [[maybe_unused]] amrex::Real xnefer{};
-    if constexpr (has_xne_xnp<T>::value) {
+    if constexpr (has_xne_xnp<T>) {
         amrex::constexpr_for<0, 16>([&] (auto K) {
             constexpr int k = K;
             constexpr auto m = electron_table_indexing::map[k];
@@ -365,13 +365,13 @@ void apply_electrons (T& state)
     amrex::Real dpepdt  = x * df_dt;
     [[maybe_unused]] amrex::Real dpepda{};
     [[maybe_unused]] amrex::Real dpepdz{};
-    if constexpr (has_pressure<T>::value) {
-        if constexpr (has_dpdA<T>::value || has_dpdZ<T>::value) {
+    if constexpr (has_pressure<T>) {
+        if constexpr (has_dpdA<T> || has_dpdZ<T>) {
             amrex::Real s       = dpepdd/state.y_e - 2.0e0_rt * din * df_d;
-            if constexpr (has_dpdA<T>::value) {
+            if constexpr (has_dpdA<T>) {
                 dpepda  = -ytot1 * (2.0e0_rt * pele + s * din);
             }
-            if constexpr (has_dpdZ<T>::value) {
+            if constexpr (has_dpdZ<T>) {
                 dpepdz  = state.rho*ytot1*(2.0e0_rt * din * df_d  +  s);
             }
         }
@@ -387,57 +387,57 @@ void apply_electrons (T& state)
     amrex::Real deepdd  = x * df_d + state.T * dsepdd;
 
     [[maybe_unused]] amrex::Real deepda{};
-    if constexpr (has_dedA<T>::value) {
+    if constexpr (has_dedA<T>) {
         amrex::Real dsepda = ytot1 * (state.y_e * df_dt * din - sele);
         deepda = -state.y_e * ytot1 * (free +  df_d * din) + state.T * dsepda;
     }
 
     [[maybe_unused]] amrex::Real deepdz{};
-    if constexpr (has_dedZ<T>::value) {
+    if constexpr (has_dedZ<T>) {
         amrex::Real dsepdz = -ytot1 * (state.y_e * df_dt * state.rho  + df_t);
         deepdz = ytot1* (free + state.y_e * df_d * state.rho) + state.T * dsepdz;
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p    = state.p + pele;
         state.dpdT = state.dpdT + dpepdt;
         state.dpdr = state.dpdr + dpepdd;
-        if constexpr (has_dpdA<T>::value) {
+        if constexpr (has_dpdA<T>) {
             state.dpdA = state.dpdA + dpepda;
         }
-        if constexpr (has_dpdZ<T>::value) {
+        if constexpr (has_dpdZ<T>) {
             state.dpdZ = state.dpdZ + dpepdz;
         }
     }
 
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s    = state.s + sele;
         state.dsdT = state.dsdT + dsepdt;
         state.dsdr = state.dsdr + dsepdd;
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e    = state.e + eele;
         state.dedT = state.dedT + deepdt;
         state.dedr = state.dedr + deepdd;
-        if constexpr (has_dedA<T>::value) {
+        if constexpr (has_dedA<T>) {
             state.dedA = state.dedA + deepda;
         }
-        if constexpr (has_dedZ<T>::value) {
+        if constexpr (has_dedZ<T>) {
             state.dedZ = state.dedZ + deepdz;
         }
     }
 
-    if constexpr (has_eta<T>::value) {
+    if constexpr (has_eta<T>) {
         state.eta = etaele;
     }
 
-    if constexpr (has_xne_xnp<T>::value) {
+    if constexpr (has_xne_xnp<T>) {
         state.xne = xnefer;
         state.xnp = 0.0e0_rt;
     }
 
-    if constexpr (has_pele_ppos<T>::value) {
+    if constexpr (has_pele_ppos<T>) {
         state.pele = pele;
         state.ppos = 0.0e0_rt;
     }
@@ -451,7 +451,7 @@ void apply_ions (T& state)
 {
     using namespace helmholtz;
 
-    constexpr amrex::Real pi      = 3.1415926535897932384e0_rt;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
     constexpr amrex::Real sioncon = (2.0e0_rt * pi * amu * kerg)/(h*h);
     constexpr amrex::Real kergavo = kerg * avo_eos;
 
@@ -480,7 +480,7 @@ void apply_ions (T& state)
     [[maybe_unused]] amrex::Real sion{};
     [[maybe_unused]] amrex::Real dsiondd{};
     [[maybe_unused]] amrex::Real dsiondt{};
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         amrex::Real x       = state.abar * state.abar * std::sqrt(state.abar) * deni / avo_eos;
         amrex::Real s       = sioncon * state.T;
         amrex::Real z       = x * s * std::sqrt(s);
@@ -493,31 +493,31 @@ void apply_ions (T& state)
             1.5e0_rt * kergavo * tempi * ytot1;
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p    = state.p + pion;
         state.dpdT = state.dpdT + dpiondt;
         state.dpdr = state.dpdr + dpiondd;
-        if constexpr (has_dpdA<T>::value) {
+        if constexpr (has_dpdA<T>) {
             state.dpdA = state.dpdA + dpionda;
         }
-        if constexpr (has_dpdZ<T>::value) {
+        if constexpr (has_dpdZ<T>) {
             state.dpdZ = state.dpdZ + dpiondz;
         }
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e    = state.e + eion;
         state.dedT = state.dedT + deiondt;
         state.dedr = state.dedr + deiondd;
-        if constexpr (has_dedA<T>::value) {
+        if constexpr (has_dedA<T>) {
             state.dedA = state.dedA + deionda;
         }
-        if constexpr (has_dedZ<T>::value) {
+        if constexpr (has_dedZ<T>) {
             state.dedZ = state.dedZ + deiondz;
         }
     }
 
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s    = state.s + sion;
         state.dsdT = state.dsdT + dsiondt;
         state.dsdr = state.dsdr + dsiondd;
@@ -568,7 +568,7 @@ void apply_radiation (T& state)
     [[maybe_unused]] amrex::Real srad{};
     [[maybe_unused]] amrex::Real dsraddd{};
     [[maybe_unused]] amrex::Real dsraddt{};
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         srad    = (prad * deni + erad) * tempi;
         dsraddd = (dpraddd * deni - prad * deni * deni + deraddd) * tempi;
         dsraddt = (dpraddt * deni + deraddt - srad) * tempi;
@@ -578,31 +578,31 @@ void apply_radiation (T& state)
     // sets these terms instead of adding to them,
     // since it comes first.
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p    = prad;
         state.dpdr = dpraddd;
         state.dpdT = dpraddt;
-        if constexpr (has_dpdA<T>::value) {
+        if constexpr (has_dpdA<T>) {
             state.dpdA = 0.0_rt;
         }
-        if constexpr (has_dpdZ<T>::value) {
+        if constexpr (has_dpdZ<T>) {
             state.dpdZ = 0.0_rt;
         }
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e    = erad;
         state.dedr = deraddd;
         state.dedT = deraddt;
-        if constexpr (has_dedA<T>::value) {
+        if constexpr (has_dedA<T>) {
             state.dedA = 0.0_rt;
         }
-        if constexpr (has_dedZ<T>::value) {
+        if constexpr (has_dedZ<T>) {
             state.dedZ = 0.0_rt;
         }
     }
 
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s    = srad;
         state.dsdr = dsraddd;
         state.dsdT = dsraddt;
@@ -630,7 +630,7 @@ void apply_coulomb_corrections (T& state)
     constexpr amrex::Real esqu = qe * qe;
     constexpr amrex::Real onethird = 1.0e0_rt/3.0e0_rt;
     constexpr amrex::Real forth = 4.0e0_rt/3.0e0_rt;
-    constexpr amrex::Real pi    = 3.1415926535897932384e0_rt;
+    constexpr amrex::Real pi = std::numbers::pi_v<amrex::Real>;
 
     [[maybe_unused]] amrex::Real pcoul    = 0.e0_rt;
     [[maybe_unused]] amrex::Real dpcouldd = 0.e0_rt;
@@ -688,7 +688,7 @@ void apply_coulomb_corrections (T& state)
         ecoul    = y * state.T * (a1 * plasg + b1 * x + c1 / x + d1);
         pcoul    = onethird * state.rho * ecoul;
 
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             scoul  = -y * (3.0e0_rt * b1 * x - 5.0e0_rt*c1 / x +
                            d1 * (std::log(plasg) - 1.0e0_rt) - e1);
         }
@@ -707,7 +707,7 @@ void apply_coulomb_corrections (T& state)
         dpcoulda = y * decoulda;
         dpcouldz = y * decouldz;
 
-       if constexpr (has_entropy<T>::value) {
+       if constexpr (has_entropy<T>) {
            y = -avo_eos * kerg / (state.abar * plasg) *
                (0.75e0_rt * b1 * x + 1.25e0_rt * c1 / x + d1);
            dscouldd = y * plasgdd;
@@ -730,7 +730,7 @@ void apply_coulomb_corrections (T& state)
         z        = c2 * x - onethird * a2 * y;
         pcoul    = -pion * z;
         ecoul    = 3.0e0_rt * pcoul / state.rho;
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             scoul = -avo_eos * ytot1 * kerg * (c2 * x -a2 * (b2 - 1.0e0_rt) / b2 * y);
         }
 
@@ -746,7 +746,7 @@ void apply_coulomb_corrections (T& state)
         decoulda = s * dpcoulda;
         decouldz = s * dpcouldz;
 
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             s = -avo_eos * kerg / (state.abar * plasg) *
                 (1.5e0_rt * c2 * x - a2 * (b2 - 1.0e0_rt) * y);
             dscouldd = s * plasgdd;
@@ -760,10 +760,10 @@ void apply_coulomb_corrections (T& state)
     amrex::Real p_temp = std::numeric_limits<amrex::Real>::max();
     amrex::Real e_temp = std::numeric_limits<amrex::Real>::max();
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         p_temp = state.p + pcoul;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         e_temp = state.e + ecoul;
     }
 
@@ -775,7 +775,7 @@ void apply_coulomb_corrections (T& state)
         ecoul    = 0.0e0_rt;
         decouldd = 0.0e0_rt;
         decouldt = 0.0e0_rt;
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             scoul    = 0.0e0_rt;
             dscouldd = 0.0e0_rt;
             dscouldt = 0.0e0_rt;
@@ -787,31 +787,31 @@ void apply_coulomb_corrections (T& state)
         decouldz = 0.0e0_rt;
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p    = state.p + pcoul;
         state.dpdr = state.dpdr + dpcouldd;
         state.dpdT = state.dpdT + dpcouldt;
-        if constexpr (has_dpdA<T>::value) {
+        if constexpr (has_dpdA<T>) {
             state.dpdA = state.dpdA + dpcoulda;
         }
-        if constexpr (has_dpdZ<T>::value) {
+        if constexpr (has_dpdZ<T>) {
             state.dpdZ = state.dpdZ + dpcouldz;
         }
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e    = state.e + ecoul;
         state.dedr = state.dedr + decouldd;
         state.dedT = state.dedT + decouldt;
-        if constexpr (has_dedA<T>::value) {
+        if constexpr (has_dedA<T>) {
             state.dedA = state.dedA + decoulda;
         }
-        if constexpr (has_dedZ<T>::value) {
+        if constexpr (has_dedZ<T>) {
             state.dedZ = state.dedZ + decouldz;
         }
     }
 
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s    = state.s + scoul;
         state.dsdr = state.dsdr + dscouldd;
         state.dsdT = state.dsdT + dscouldt;
@@ -838,7 +838,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_rh) {
 
-        if constexpr (has_enthalpy<T>::value) {
+        if constexpr (has_enthalpy<T>) {
             v_want = state.h;
             var  = ienth;
             dvar = itemp;
@@ -847,7 +847,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_tp) {
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             v_want = state.p;
             var  = ipres;
             dvar = idens;
@@ -856,7 +856,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_rp) {
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             v_want = state.p;
             var  = ipres;
             dvar = itemp;
@@ -865,7 +865,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_re) {
 
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             v_want = state.e;
             var  = iener;
             dvar = itemp;
@@ -874,7 +874,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_ps) {
 
-        if constexpr (has_pressure<T>::value && has_entropy<T>::value) {
+        if constexpr (has_pressure<T> && has_entropy<T>) {
             single_iter = false;
             v1_want = state.p;
             v2_want = state.s;
@@ -885,7 +885,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_ph) {
 
-        if constexpr (has_pressure<T>::value && has_enthalpy<T>::value) {
+        if constexpr (has_pressure<T> && has_enthalpy<T>) {
             single_iter = false;
             v1_want = state.p;
             v2_want = state.h;
@@ -896,7 +896,7 @@ void prepare_for_iterations (I input, const T& state,
     }
     else if (input == eos_input_th) {
 
-        if constexpr (has_enthalpy<T>::value) {
+        if constexpr (has_enthalpy<T>) {
             v_want = state.h;
             var  = ienth;
             dvar = idens;
@@ -932,25 +932,25 @@ void single_iter_update (T& state, int var, int dvar,
         xtol = ttol;
 
         if (var == ipres) {
-            if constexpr (has_pressure<T>::value) {
+            if constexpr (has_pressure<T>) {
                 v    = state.p;
                 dvdx = state.dpdT;
             }
         }
         else if (var == iener) {
-            if constexpr (has_energy<T>::value) {
+            if constexpr (has_energy<T>) {
                 v    = state.e;
                 dvdx = state.dedT;
             }
         }
         else if (var == ientr) {
-            if constexpr (has_entropy<T>::value) {
+            if constexpr (has_entropy<T>) {
                 v    = state.s;
                 dvdx = state.dsdT;
             }
         }
         else if (var == ienth) {
-            if constexpr (has_enthalpy<T>::value) {
+            if constexpr (has_enthalpy<T>) {
                 v    = state.h;
                 dvdx = state.dhdT;
             }
@@ -965,25 +965,25 @@ void single_iter_update (T& state, int var, int dvar,
         xtol = dtol;
 
         if (var == ipres) {
-            if constexpr (has_pressure<T>::value) {
+            if constexpr (has_pressure<T>) {
                 v    = state.p;
                 dvdx = state.dpdr;
             }
         }
         else if (var == iener) {
-            if constexpr (has_energy<T>::value) {
+            if constexpr (has_energy<T>) {
                 v    = state.e;
                 dvdx = state.dedr;
             }
         }
         else if (var == ientr) {
-            if constexpr (has_entropy<T>::value) {
+            if constexpr (has_entropy<T>) {
                 v    = state.s;
                 dvdx = state.dsdr;
             }
         }
         else if (var == ienth) {
-            if constexpr (has_enthalpy<T>::value) {
+            if constexpr (has_enthalpy<T>) {
                 v    = state.h;
                 dvdx = state.dhdr;
             }
@@ -1035,28 +1035,28 @@ void double_iter_update (T& state, int var1, int var2,
     amrex::Real rold = state.rho;
 
     if (var1 == ipres) {
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             v1    = state.p;
             dv1dt = state.dpdT;
             dv1dr = state.dpdr;
         }
     }
     else if (var1 == iener) {
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             v1    = state.e;
             dv1dt = state.dedT;
             dv1dr = state.dedr;
         }
     }
     else if (var1 == ientr) {
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             v1    = state.s;
             dv1dt = state.dsdT;
             dv1dr = state.dsdr;
         }
     }
     else if (var1 == ienth) {
-        if constexpr (has_enthalpy<T>::value) {
+        if constexpr (has_enthalpy<T>) {
             v1    = state.h;
             dv1dt = state.dhdT;
             dv1dr = state.dhdr;
@@ -1064,28 +1064,28 @@ void double_iter_update (T& state, int var1, int var2,
     }
 
     if (var2 == ipres) {
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             v2    = state.p;
             dv2dt = state.dpdT;
             dv2dr = state.dpdr;
         }
     }
     else if (var2 == iener) {
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             v2    = state.e;
             dv2dt = state.dedT;
             dv2dr = state.dedr;
         }
     }
     else if (var2 == ientr) {
-        if constexpr (has_entropy<T>::value) {
+        if constexpr (has_entropy<T>) {
             v2    = state.s;
             dv2dt = state.dsdT;
             dv2dr = state.dsdr;
         }
     }
     else if (var2 == ienth) {
-        if constexpr (has_enthalpy<T>::value) {
+        if constexpr (has_enthalpy<T>) {
             v2    = state.h;
             dv2dt = state.dhdT;
             dv2dr = state.dhdr;
@@ -1141,16 +1141,16 @@ void finalize_state (I input, T& state,
     state.mu = 1.0e0_rt / (1.0e0_rt / state.abar + 1.0e0_rt / state.mu_e);
 
     // Calculate some remaining derivatives
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpde = state.dpdT / state.dedT;
         state.dpdr_e = state.dpdr - state.dpdT * state.dedr / state.dedT;
     }
 
     // Specific heats and Gamma_1
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.cv = state.dedT;
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             amrex::Real chit = state.T / state.p * state.dpdT;
             amrex::Real chid = state.dpdr * state.rho / state.p;
 
@@ -1161,63 +1161,63 @@ void finalize_state (I input, T& state,
 
     // Use the non-relativistic version of the sound speed, cs = sqrt(gam_1 * P / rho).
     // This replaces the relativistic version that comes out of helmeos.
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.cs = std::sqrt(state.gam1 * state.p / state.rho);
     }
 
     if (input_is_constant) {
 
-       if (input == eos_input_rh) {
+       if (input == eos_input_rh) {  // NOLINT(bugprone-branch-clone)
 
-           if constexpr (has_enthalpy<T>::value) {
+           if constexpr (has_enthalpy<T>) {
                state.h = v_want;
            }
 
        }
-       else if (input == eos_input_tp) {
+       else if (input == eos_input_tp) {  // NOLINT(bugprone-branch-clone)
 
-           if constexpr (has_pressure<T>::value) {
+           if constexpr (has_pressure<T>) {
                state.p = v_want;
            }
 
        }
        else if (input == eos_input_rp) {
 
-           if constexpr (has_pressure<T>::value) {
+           if constexpr (has_pressure<T>) {
                state.p = v_want;
            }
 
        }
        else if (input == eos_input_re) {
 
-           if constexpr (has_energy<T>::value) {
+           if constexpr (has_energy<T>) {
                state.e = v_want;
            }
 
        }
        else if (input == eos_input_ps) {
 
-           if constexpr (has_pressure<T>::value) {
+           if constexpr (has_pressure<T>) {
                state.p = v1_want;
            }
-           if constexpr (has_entropy<T>::value) {
+           if constexpr (has_entropy<T>) {
                state.s = v2_want;
            }
 
        }
        else if (input == eos_input_ph) {
 
-           if constexpr (has_pressure<T>::value) {
+           if constexpr (has_pressure<T>) {
                state.p = v1_want;
            }
-           if constexpr (has_enthalpy<T>::value) {
+           if constexpr (has_enthalpy<T>) {
                state.h = v2_want;
            }
 
        }
        else if (input == eos_input_th) {
 
-           if constexpr (has_enthalpy<T>::value) {
+           if constexpr (has_enthalpy<T>) {
                state.h = v_want;
            }
 
@@ -1270,7 +1270,7 @@ void actual_eos (I input, T& state)
 
         // Calculate enthalpy the usual way, h = e + p / rho.
 
-        if constexpr (has_enthalpy<T>::value) {
+        if constexpr (has_enthalpy<T>) {
             state.h = state.e + state.p / state.rho;
             state.dhdr = state.dedr + state.dpdr / state.rho - state.p / (state.rho * state.rho);
             state.dhdT = state.dedT + state.dpdT / state.rho;
diff --git a/EOS/metal_chem/actual_eos.H b/EOS/metal_chem/actual_eos.H
index bf903302df..09326500e4 100644
--- a/EOS/metal_chem/actual_eos.H
+++ b/EOS/metal_chem/actual_eos.H
@@ -144,7 +144,7 @@ void actual_eos (I input, T& state)
         // dens, energy, and xmass are inputs
         // Solve for the temperature
 
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             dens = state.rho;
 
             // stop the integration if the internal energy < 0
@@ -174,7 +174,7 @@ void actual_eos (I input, T& state)
     case eos_input_rp:
         // dens, pressure, and xmass are inputs
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             dens = state.rho;
 
             // stop the integration if the pressure < 0
@@ -226,28 +226,28 @@ void actual_eos (I input, T& state)
     state.T   = temp;
     state.rho = dens;
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e = eint;
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         amrex::Real pressure = state.rho * eint / sum_gammasinv;
         state.p = pressure;
 
         state.dpdT = pressure / temp;
         state.dpdr = pressure / dens;
         state.cs = std::sqrt((1.0 + 1.0/sum_gammasinv) * state.p /state.rho);
-        if constexpr (has_G<T>::value) {
+        if constexpr (has_G<T>) {
             state.G = 0.5 * (1.0 + (1.0 + 1.0/sum_gammasinv));
         }
     }
 
     amrex::Real dedT = sum_gammasinv * sum_Abarinv * gasconstant;
     amrex::Real dedr = 0.0_rt;
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedT = dedT;
         state.dedr = dedr;
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             state.dpde = state.dpdT / state.dedT;
         }
     }
diff --git a/EOS/multigamma/actual_eos.H b/EOS/multigamma/actual_eos.H
index 61bb9f9f17..5d712596ca 100644
--- a/EOS/multigamma/actual_eos.H
+++ b/EOS/multigamma/actual_eos.H
@@ -107,7 +107,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the temperature:
         // h = e + p/rho = (p/rho)*[1 + 1/(gamma-1)] = (p/rho)*gamma/(gamma-1)
-        if constexpr (has_enthalpy<T>::value) {
+        if constexpr (has_enthalpy<T>) {
             dens = state.rho;
             temp = state.h * (m_nucleon / C::k_B) / sumYg_gm1;
         }
@@ -120,7 +120,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the density:
         // p = rho k T / (abar m_nucleon)
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             dens = state.p * state.abar * (m_nucleon / C::k_B) / state.T;
             temp = state.T;
         }
@@ -133,7 +133,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the temperature:
         // p = rho k T / (mu m_nucleon)
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             dens = state.rho;
             temp = state.p * state.abar * (m_nucleon / C::k_B) / state.rho;
         }
@@ -146,7 +146,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the temperature
         // e = k T / [(mu m_nucleon)*(gamma-1)]
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             dens = state.rho;
             temp = state.e * (m_nucleon / C::k_B) / sumY_gm1;
         }
@@ -167,7 +167,7 @@ void actual_eos (I input, T& state)
         // pressure, enthalpy and xmass are inputs
 
         // Solve for temperature and density
-        if constexpr (has_pressure<T>::value && has_enthalpy<T>::value) {
+        if constexpr (has_pressure<T> && has_enthalpy<T>) {
             dens = state.p * state.abar / state.h * sumYg_gm1;
             temp = state.p * state.abar * (m_nucleon / C::k_B) / dens;
         }
@@ -206,21 +206,21 @@ void actual_eos (I input, T& state)
 
     // Compute the pressure simply from the ideal gas law, and the
     // specific internal energy using the gamma-law EOS relation.
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p = pres;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e = ener;
     }
 
     // Enthalpy is h = e + p/rho
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.h = enth;
     }
 
     // entropy (per gram) -- this is wrong. Not sure what the expression
     // is for a multigamma gas
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s = ((C::k_B / m_nucleon) / state.abar) *
                   (2.5_rt + std::log((std::pow(state.abar * m_nucleon, 2.5_rt) / dens) *
                                      std::pow(C::k_B * temp, 1.5_rt) /
@@ -230,36 +230,36 @@ void actual_eos (I input, T& state)
     // Compute the thermodynamic derivatives and specific heats
     amrex::Real dpdT = pres / temp;
     amrex::Real dpdr = pres / dens;
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpdT = dpdT;
         state.dpdr = dpdr;
     }
 
     amrex::Real dedT = ener / temp;
     amrex::Real dedr = 0.0_rt;
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedT = dedT;
         state.dedr = dedr;
     }
 
     amrex::Real dsdT = 0.0_rt;
     amrex::Real dsdr = 0.0_rt;
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.dsdT = dsdT;
         state.dsdr = dsdr;
     }
 
     amrex::Real dhdT = dedT + dpdT / dens;
     amrex::Real dhdr = 0.0_rt;
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.dhdT = dhdT;
         state.dhdr = dhdr;
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.cv = state.dedT;
     }
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.cp = enth / state.T;
 
         state.gam1 = state.cp / state.cv;
@@ -269,22 +269,22 @@ void actual_eos (I input, T& state)
     }
 
     // These need to be worked out.
-    if constexpr (has_dpdA<T>::value) {
+    if constexpr (has_dpdA<T>) {
         state.dpdA = 0.0_rt;
     }
-    if constexpr (has_dpdZ<T>::value) {
+    if constexpr (has_dpdZ<T>) {
         state.dpdZ = 0.0_rt;
     }
 
-    if constexpr (has_dedA<T>::value) {
+    if constexpr (has_dedA<T>) {
         state.dedA = 0.0_rt;
     }
-    if constexpr (has_dedZ<T>::value) {
+    if constexpr (has_dedZ<T>) {
         state.dedZ = 0.0_rt;
     }
 
     // Sound speed
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.cs = std::sqrt(state.gam1 * state.p / dens);
     }
 }
diff --git a/EOS/polytrope/actual_eos.H b/EOS/polytrope/actual_eos.H
index 1db5a9f61b..5ee49f3f2e 100644
--- a/EOS/polytrope/actual_eos.H
+++ b/EOS/polytrope/actual_eos.H
@@ -127,22 +127,22 @@ void actual_eos (I input, T& state)
     amrex::Real temp = state.T;
 
     amrex::Real pres = 1.0_rt;
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         pres = state.p;
     }
 
     amrex::Real enth = 1.0_rt;
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         enth = state.h;
     }
 
     amrex::Real eint = 1.0_rt;
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         eint = state.e;
     }
 
     amrex::Real entr = 1.0_rt;
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         entr = state.s;
     }
 
@@ -266,46 +266,46 @@ void actual_eos (I input, T& state)
 
     state.T   = temp;
     state.rho = dens;
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.h   = enth;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s   = entr;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e   = eint;
     }
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p   = pres;
     }
 
     // Compute the thermodynamic derivatives and specific heats
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpdT = 0.0_rt;
         state.dpdr = gamma_const * pres / dens;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedT = 0.0_rt;
         state.dedr = pres / (dens * dens);
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.dsdT = 0.0_rt;
         state.dsdr = 0.0_rt;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.dhdT = 0.0_rt;
         state.dhdr = state.dedr + gm1 * pres / (dens * dens);
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpde = 0.0_rt;
         state.dpdr_e = gamma_const * pres / dens;
     }
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.cv = state.dedT;
     }
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.cp = gamma_const * state.cv;
 
         state.gam1 = gamma_const;
@@ -313,22 +313,22 @@ void actual_eos (I input, T& state)
 
     // Compute dpdX, dedX, dhdX.
 
-    if constexpr (has_dpdA<T>::value) {
+    if constexpr (has_dpdA<T>) {
         state.dpdA = - state.p / state.abar;
     }
-    if constexpr (has_dpdZ<T>::value) {
+    if constexpr (has_dpdZ<T>) {
         state.dpdZ =   state.p / (1.0_rt + state.zbar);
     }
 
-    if constexpr (has_dedA<T>::value) {
+    if constexpr (has_dedA<T>) {
         state.dedA = - state.e / state.abar;
     }
-    if constexpr (has_dedZ<T>::value) {
+    if constexpr (has_dedZ<T>) {
         state.dedZ =   state.e / (1.0_rt + state.zbar);
     }
 
     // sound speed
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.cs = std::sqrt(gamma_const * pres / dens);
     }
 }
diff --git a/EOS/primordial_chem/actual_eos.H b/EOS/primordial_chem/actual_eos.H
index f5eecca717..9eb74b3236 100644
--- a/EOS/primordial_chem/actual_eos.H
+++ b/EOS/primordial_chem/actual_eos.H
@@ -125,7 +125,7 @@ void actual_eos (I input, T& state)
         // dens, energy, and xmass are inputs
         // Solve for the temperature
 
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             dens = state.rho;
 
             // stop the integration if the internal energy < 0
@@ -155,7 +155,7 @@ void actual_eos (I input, T& state)
     case eos_input_rp:
         // dens, pressure, and xmass are inputs
 
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             dens = state.rho;
 
             // stop the integration if the pressure < 0
@@ -207,28 +207,28 @@ void actual_eos (I input, T& state)
     state.T   = temp;
     state.rho = dens;
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e = eint;
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         amrex::Real pressure = state.rho * eint / sum_gammasinv;
         state.p = pressure;
 
         state.dpdT = pressure / temp;
         state.dpdr = pressure / dens;
         state.cs = std::sqrt((1.0 + 1.0/sum_gammasinv) * state.p /state.rho);
-        if constexpr (has_G<T>::value) {
+        if constexpr (has_G<T>) {
             state.G = 0.5 * (1.0 + (1.0 + 1.0/sum_gammasinv));
         }
     }
 
     amrex::Real dedT = sum_gammasinv * sum_Abarinv * gasconstant;
     amrex::Real dedr = 0.0_rt;
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedT = dedT;
         state.dedr = dedr;
-        if constexpr (has_pressure<T>::value) {
+        if constexpr (has_pressure<T>) {
             state.dpde = state.dpdT / state.dedT;
         }
     }
diff --git a/EOS/rad_power_law/actual_eos.H b/EOS/rad_power_law/actual_eos.H
index c259fe5d31..fbe2d9e73b 100644
--- a/EOS/rad_power_law/actual_eos.H
+++ b/EOS/rad_power_law/actual_eos.H
@@ -75,7 +75,7 @@ void actual_eos (I input, T& state)
 
     case eos_input_rt:
 
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             state.cv = eos_rp::eos_const_c_v * std::pow(state.rho, eos_rp::eos_c_v_exp_m) * std::pow(state.T, -eos_rp::eos_c_v_exp_n);
             state.e  = eos_rp::eos_const_c_v * std::pow(state.rho, eos_rp::eos_c_v_exp_m) * std::pow(state.T, 1 - eos_rp::eos_c_v_exp_n) / (1 - eos_rp::eos_c_v_exp_n);
         }
@@ -84,7 +84,7 @@ void actual_eos (I input, T& state)
 
     case eos_input_re:
 
-        if constexpr (has_energy<T>::value) {
+        if constexpr (has_energy<T>) {
             state.T  = std::pow((1 - eos_rp::eos_c_v_exp_n) * state.e * std::pow(state.rho, -eos_rp::eos_c_v_exp_m) / eos_rp::eos_const_c_v, 1.0_rt / (1.0_rt - eos_rp::eos_c_v_exp_n));
             state.cv = eos_rp::eos_const_c_v * std::pow(state.rho, eos_rp::eos_c_v_exp_m) * std::pow(state.T, -eos_rp::eos_c_v_exp_n);
         }
@@ -104,23 +104,23 @@ void actual_eos (I input, T& state)
     // Set some data to nonsense values so that things intentionally go wrong
     // if this EOS is somehow used for hydro.
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p    = -1.e0_rt;
     }
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.gam1 = -1.e0_rt;
         state.cs   = -1.e0_rt;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s    = -1.e0_rt;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.h    = -1.e0_rt;
     }
 
     // Give dpdr a value for the purposes of the composition_derivatives routine.
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpdr = 0.e0_rt;
     }
 }
diff --git a/EOS/tillotson/actual_eos.H b/EOS/tillotson/actual_eos.H
index 883f5fd010..7693e7886a 100644
--- a/EOS/tillotson/actual_eos.H
+++ b/EOS/tillotson/actual_eos.H
@@ -126,7 +126,7 @@ void actual_eos (I input, T& state)
 
     // Pressure and pressure derivatives
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p = pres;
         state.dpdr_e = dpdr_e;
         state.dpde = dpde;
@@ -141,13 +141,13 @@ void actual_eos (I input, T& state)
 
     // Sound speed
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.cs = std::sqrt(state.dpdr_e);
     }
 
     // Compute gamma from sound speed
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.gam1 = state.cs * state.cs * state.rho / state.p;
     }
 
@@ -155,7 +155,7 @@ void actual_eos (I input, T& state)
 
     state.cv = eos_rp::eos_c_v;
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.cp = state.gam1 * eos_rp::eos_c_v;
     }
 
@@ -169,13 +169,13 @@ void actual_eos (I input, T& state)
 
     // Usual formula for enthalpy
 
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.h = state.e + state.p / state.rho;
     }
 
     // Arbitrary value for entropy
 
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s = 1.0_rt;
     }
 }
diff --git a/EOS/ztwd/actual_eos.H b/EOS/ztwd/actual_eos.H
index a2d6ff3d6b..264174952e 100644
--- a/EOS/ztwd/actual_eos.H
+++ b/EOS/ztwd/actual_eos.H
@@ -136,22 +136,22 @@ void actual_eos (I input, T& state)
     amrex::Real temp = state.T;
 
     amrex::Real pres = 1.0_rt;
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         pres = state.p;
     }
 
     amrex::Real enth = 1.0_rt;
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         enth = state.h;
     }
 
     amrex::Real eint = 1.0_rt;
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         eint = state.e;
     }
 
     amrex::Real entr = 1.0_rt;
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         entr = state.s;
     }
 
@@ -289,31 +289,31 @@ void actual_eos (I input, T& state)
 
     state.T   = temp;
     state.rho = dens;
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.h   = enth;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.s   = entr;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.e   = eint;
     }
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.p   = pres;
     }
 
     // All temperature derivatives are zero since the gas is temperature-independent.
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpdT = 0.0_rt;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.dhdT = 0.0_rt;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedT = 0.0_rt;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.dsdT = 0.0_rt;
     }
 
@@ -325,53 +325,53 @@ void actual_eos (I input, T& state)
     amrex::Real dpdr = dxdr * dpdx(x);
     amrex::Real dhdr = dxdr * dhdx(x, B);
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         state.dpdr = dpdr;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         state.dhdr = dhdr;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.dedr = dhdr - dpdr / state.rho + pres / (state.rho * state.rho);
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         state.dsdr = 0.0_rt;
     }
 
     // Heat capacities are zero: the gas properties don't change when the temperature changes.
 
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         state.cv = 0.0_rt;
     }
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.cp = 0.0_rt;
     }
 
     // Adiabatic gamma_1 == d(log p) / d(log rho) |_s.
 
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.gam1 = state.dpdr * (state.rho / state.p);
     }
 
     // Derivatives with respect to A and Z.
 
-    if constexpr (has_dpdA<T>::value) {
+    if constexpr (has_dpdA<T>) {
         state.dpdA = -state.p / state.abar;
     }
-    if constexpr (has_dpdZ<T>::value) {
+    if constexpr (has_dpdZ<T>) {
         state.dpdZ =  state.p / (1.0_rt + state.zbar);
     }
 
-    if constexpr (has_dedA<T>::value) {
+    if constexpr (has_dedA<T>) {
         state.dedA = -state.e / state.abar;
     }
-    if constexpr (has_dedZ<T>::value) {
+    if constexpr (has_dedZ<T>) {
         state.dedZ =  state.e / (1.0_rt + state.zbar);
     }
 
     // Sound speed.
 
-    if constexpr (has_pressure<T>::value && has_energy<T>::value) {
+    if constexpr (has_pressure<T> && has_energy<T>) {
         state.cs = std::sqrt(state.dpdr);
     }
 }
diff --git a/Make.Microphysics b/Make.Microphysics
index 5dfd0a9737..0254bc7fd3 100644
--- a/Make.Microphysics
+++ b/Make.Microphysics
@@ -32,8 +32,8 @@ ifeq ($(USE_MICROPHYSICS_DEBUG),TRUE)
   DEFINES += -DMICROPHYSICS_DEBUG
 endif
 
-# Require C++17
-CXXSTD := c++17
+# Require C++20
+CXXSTD := c++20
 
 ifeq ($(USE_GPU_PRINTF),TRUE)
   DEFINES += -DALLOW_GPU_PRINTF
diff --git a/integration/integrator_data.H b/integration/integrator_data.H
index 5fb00aed4d..e39d369746 100644
--- a/integration/integrator_data.H
+++ b/integration/integrator_data.H
@@ -42,11 +42,11 @@ constexpr int integrator_neqs ()
     int int_neqs = 0;
 
 #ifdef STRANG
-    if constexpr (has_xn<BurnT>::value) {
+    if constexpr (has_xn<BurnT>) {
         int_neqs += NumSpec;
     }
 
-    if constexpr (has_energy<BurnT>::value) {
+    if constexpr (has_energy<BurnT>) {
         int_neqs += 1;
     }
 #else
diff --git a/interfaces/eos.H b/interfaces/eos.H
index 3abccc2c77..96340e65b5 100644
--- a/interfaces/eos.H
+++ b/interfaces/eos.H
@@ -98,7 +98,7 @@ template <typename T>
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void reset_e (T& state, bool& has_been_reset)
 {
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         if (state.e < EOSData::mine || state.e > EOSData::maxe) {
             eos_reset(state, has_been_reset);
         }
@@ -109,7 +109,7 @@ template <typename T>
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void reset_h (T& state, bool& has_been_reset)
 {
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         if (state.h < EOSData::minh || state.h > EOSData::maxh) {
             eos_reset(state, has_been_reset);
         }
@@ -120,7 +120,7 @@ template <typename T>
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void reset_s (T& state, bool& has_been_reset)
 {
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         if (state.s < EOSData::mins || state.s > EOSData::maxs) {
             eos_reset(state, has_been_reset);
         }
@@ -131,7 +131,7 @@ template <typename T>
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
 void reset_p (T& state, bool& has_been_reset)
 {
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         if (state.p < EOSData::minp || state.p > EOSData::maxp) {
             eos_reset(state, has_been_reset);
         }
@@ -252,7 +252,7 @@ template <typename T>
 AMREX_INLINE
 void check_e (T& state)
 {
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         if (state.e < EOSData::mine) {
             print_state(state);
             amrex::Error("EOS: e smaller than mine.");
@@ -268,7 +268,7 @@ template <typename T>
 AMREX_INLINE
 void check_h (T& state)
 {
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         if (state.h < EOSData::minh) {
             print_state(state);
             amrex::Error("EOS: h smaller than minh.");
@@ -285,7 +285,7 @@ template <typename T>
 AMREX_INLINE
 void check_s (T& state)
 {
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         if (state.s < EOSData::mins) {
             print_state(state);
             amrex::Error("EOS: s smaller than mins.");
@@ -301,7 +301,7 @@ template <typename T>
 AMREX_INLINE
 void check_p (T& state)
 {
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         if (state.p < EOSData::minp) {
             print_state(state);
             amrex::Error("EOS: p smaller than minp.");
@@ -410,7 +410,7 @@ void eos (const I input, T& state)
   }
 #endif
 
-  if constexpr (has_xn<T>::value) {
+  if constexpr (has_xn<T>) {
       // Get abar, zbar, etc.
       composition(state);
   }
diff --git a/interfaces/eos_type.H b/interfaces/eos_type.H
index 7db2db24a6..fcb0280817 100644
--- a/interfaces/eos_type.H
+++ b/interfaces/eos_type.H
@@ -164,131 +164,56 @@ struct chem_eos_t:eos_base_t {
 };
 
 
-template <typename T, typename Enable = void>
-struct has_energy
-    : std::false_type {};
-
 template <typename T>
-struct has_energy<T, decltype((void)T::e, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_enthalpy
-    : std::false_type {};
+concept has_energy = requires { std::declval<T&>().e; };
 
 template <typename T>
-struct has_enthalpy<T, decltype((void)T::h, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_entropy
-    : std::false_type {};
+concept has_enthalpy = requires { std::declval<T&>().h; };
 
 template <typename T>
-struct has_entropy<T, decltype((void)T::s, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_pressure
-    : std::false_type {};
+concept has_entropy = requires { std::declval<T&>().s; };
 
 template <typename T>
-struct has_pressure<T, decltype((void)T::p, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_dpdA
-    : std::false_type {};
+concept has_pressure = requires { std::declval<T&>().p; };
 
 template <typename T>
-struct has_dpdA<T, decltype((void)T::dpdA, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_G
-    : std::false_type {};
+concept has_dpdA = requires { std::declval<T&>().dpdA; };
 
 template <typename T>
-struct has_G<T, decltype((void)T::G, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_dpdZ
-    : std::false_type {};
+concept has_dpdZ = requires { std::declval<T&>().dpdZ; };
 
 template <typename T>
-struct has_dpdZ<T, decltype((void)T::dpdZ, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_dedA
-    : std::false_type {};
+concept has_dedA = requires { std::declval<T&>().dedA; };
 
 template <typename T>
-struct has_dedA<T, decltype((void)T::dedA, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_dedZ
-    : std::false_type {};
+concept has_dedZ = requires { std::declval<T&>().dedZ; };
 
 template <typename T>
-struct has_dedZ<T, decltype((void)T::dedZ, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_pele_ppos
-    : std::false_type {};
+concept has_G = requires { std::declval<T&>().G; };
 
 template <typename T>
-struct has_pele_ppos<T, decltype((void)T::pele, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_xne_xnp
-    : std::false_type {};
+concept has_pele_ppos = requires { std::declval<T&>().pele; };
 
 template <typename T>
-struct has_xne_xnp<T, decltype((void)T::xne, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_eta
-    : std::false_type {};
+concept has_xne_xnp = requires { std::declval<T&>().xne; };
 
 template <typename T>
-struct has_eta<T, decltype((void)T::eta, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_conductivity
-    : std::false_type {};
+concept has_eta = requires { std::declval<T&>().eta; };
 
 template <typename T>
-struct has_conductivity<T, decltype((void)T::conductivity, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_xn
-    : std::false_type {};
+concept has_conductivity = requires { std::declval<T&>().conductivity; };
 
 template <typename T>
-struct has_xn<T, decltype((void)T::xn, void())>
-    : std::true_type {};
-
-template <typename T, typename Enable = void>
-struct has_base_variables
-    : std::false_type {};
+concept has_xn = requires { std::declval<T&>().xn; };
 
 template <typename T>
-struct has_base_variables<T, decltype((void)T::rho, void())>
-    : std::true_type {};
+concept has_base_variables = requires { std::declval<T&>().rho; };
 
 template <typename T>
 inline
 std::ostream& print_state (std::ostream& o, T const& eos_state)
 {
-    if constexpr (has_base_variables<T>::value) {
+    if constexpr (has_base_variables<T>) {
         o << "rho = " << eos_state.rho << std::endl;
         o << "T =   " << eos_state.T << std::endl;
         o << "xn = ";
@@ -304,86 +229,86 @@ std::ostream& print_state (std::ostream& o, T const& eos_state)
         o << std::endl;
 #endif
     }
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         o << "p = " << eos_state.p << std::endl;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         o << "e = " << eos_state.e << std::endl;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         o << "h = " << eos_state.h << std::endl;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         o << "s = " << eos_state.s << std::endl;
     }
 
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         o << "dpdT = " << eos_state.dpdT << std::endl;
         o << "dpdr = " << eos_state.dpdr << std::endl;
     }
-    if constexpr (has_energy<T>::value) {
+    if constexpr (has_energy<T>) {
         o << "dedT = " << eos_state.dedT << std::endl;
         o << "dedr = " << eos_state.dedr << std::endl;
     }
-    if constexpr (has_enthalpy<T>::value) {
+    if constexpr (has_enthalpy<T>) {
         o << "dhdT = " << eos_state.dhdT << std::endl;
         o << "dhdr = " << eos_state.dhdr << std::endl;
     }
-    if constexpr (has_entropy<T>::value) {
+    if constexpr (has_entropy<T>) {
         o << "dsdT = " << eos_state.dsdT << std::endl;
         o << "dsdr = " << eos_state.dsdr << std::endl;
     }
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         o << "dpde = " << eos_state.dpde << std::endl;
         o << "dpdr_e = " << eos_state.dpdr_e  << std::endl;
     }
 
-    if constexpr (has_base_variables<T>::value) {
+    if constexpr (has_base_variables<T>) {
         o << "cv = " << eos_state.cv << std::endl;
     }
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         o << "cp = " << eos_state.cp << std::endl;
     }
-    if constexpr (has_xne_xnp<T>::value) {
+    if constexpr (has_xne_xnp<T>) {
         o << "xne = " << eos_state.xne << std::endl;
         o << "xnp = " << eos_state.xnp << std::endl;
     }
-    if constexpr (has_eta<T>::value) {
+    if constexpr (has_eta<T>) {
         o << "eta = " << eos_state.eta << std::endl;
     }
-    if constexpr (has_pele_ppos<T>::value) {
+    if constexpr (has_pele_ppos<T>) {
         o << "pele = " << eos_state.pele << std::endl;
         o << "ppos = " << eos_state.ppos << std::endl;
     }
-    if constexpr (has_base_variables<T>::value) {
+    if constexpr (has_base_variables<T>) {
         o << "mu = " << eos_state.mu << std::endl;
         o << "mu_e = " << eos_state.mu_e << std::endl;
         o << "y_e = " << eos_state.y_e << std::endl;
     }
-    if constexpr (has_pressure<T>::value) {
+    if constexpr (has_pressure<T>) {
         o << "gam1 = " << eos_state.gam1 << std::endl;
         o << "cs = " << eos_state.cs << std::endl;
     }
 
-    if constexpr (has_base_variables<T>::value) {
+    if constexpr (has_base_variables<T>) {
         o << "abar = " << eos_state.abar << std::endl;
         o << "zbar = " << eos_state.zbar << std::endl;
     }
 
-    if constexpr (has_dpdA<T>::value) {
+    if constexpr (has_dpdA<T>) {
         o << "dpdA = " << eos_state.dpdA << std::endl;
     }
-    if constexpr (has_dpdZ<T>::value) {
+    if constexpr (has_dpdZ<T>) {
         o << "dpdZ = " << eos_state.dpdZ << std::endl;
     }
-    if constexpr (has_dedA<T>::value) {
+    if constexpr (has_dedA<T>) {
         o << "dedA = " << eos_state.dedA << std::endl;
     }
-    if constexpr (has_dedZ<T>::value) {
+    if constexpr (has_dedZ<T>) {
         o << "dedZ = " << eos_state.dedZ << std::endl;
     }
 #ifdef CONDUCTIVITY
-    if constexpr (has_conductivity<T>::value) {
+    if constexpr (has_conductivity<T>) {
         o << "conductivity = " << eos_state.conductivity << std::endl;
     }
 #endif
diff --git a/rates/aprox_rates.H b/rates/aprox_rates.H
index 1193a56017..68feb1310a 100644
--- a/rates/aprox_rates.H
+++ b/rates/aprox_rates.H
@@ -184,7 +184,7 @@ void rate_triplealf(const tf_t<number_t>& tf, const amrex::Real den, number_t& f
     // this is a(a,g)be8
     const number_t aa    = 7.40e5_rt * tf.t9i32 * admath::exp(-1.0663_rt*tf.t9i);
     const number_t bb    = 4.164e9_rt * tf.t9i23 * admath::exp(-13.49_rt*tf.t9i13 - tf.t92*q1);
-    const number_t cc    = 1.0e0_rt + 0.031_rt*tf.t913 + 8.009_rt*tf.t923 + 1.732_rt*tf.t9
+    const number_t cc    = 1.0e0_rt + 0.031_rt*tf.t913 + 8.009_rt*tf.t923 + 1.732_rt*tf.t9  // NOLINT(modernize-use-std-numbers)
           + 49.883_rt*tf.t943 + 27.426_rt*tf.t953;
     const number_t r2abe    = aa + bb * cc;
 
@@ -1384,13 +1384,13 @@ void ecapnuc(const amrex::Real etakep, const amrex::Real temp,
 {
     constexpr amrex::Real qn1    = -2.0716446e-06_rt;
     constexpr amrex::Real ftinv  = 1.0e0_rt/1083.9269e0_rt;
-    constexpr amrex::Real twoln  = 0.6931472e0_rt;
+    constexpr amrex::Real twoln  = 0.6931472e0_rt;  // NOLINT(modernize-use-std-numbers)
     constexpr amrex::Real cmk5   = 1.3635675e-49_rt;
     constexpr amrex::Real cmk6   = 2.2993864e-59_rt;
     constexpr amrex::Real bk     = 1.38062e-16_rt;
     constexpr amrex::Real qn2    = 2.0716446e-06_rt;
     constexpr amrex::Real c2me   = 8.1872665e-07_rt;
-    constexpr amrex::Real pi     = 3.1415927e0_rt;
+    constexpr amrex::Real pi     = 3.1415927e0_rt;  // NOLINT(modernize-use-std-numbers)
     constexpr amrex::Real pi2    = pi * pi;
 
     // c2me is the constant used to convert the neutrino energy
diff --git a/unit_test/CMakeLists.txt b/unit_test/CMakeLists.txt
index 7310d8ace1..c5d88877d4 100644
--- a/unit_test/CMakeLists.txt
+++ b/unit_test/CMakeLists.txt
@@ -1,4 +1,6 @@
 
+set(CMAKE_CXX_STANDARD 20)
+
 if(BUILD_UNIT_TEST_PC)
   #add primordial_chem test
   add_subdirectory(burn_cell_primordial_chem)
diff --git a/util/approx_math/approx_math.H b/util/approx_math/approx_math.H
index 97fa6b43b3..9ad26284c6 100644
--- a/util/approx_math/approx_math.H
+++ b/util/approx_math/approx_math.H
@@ -8,6 +8,12 @@
 #include <cstdint>
 #include <cstring>
 
+// some of these algorithms need constants in particular precision, so
+// we use the hardcoded values instead of getting them from
+// std::numbers
+
+// NOLINTBEGIN(modernize-use-std-numbers)
+
 using namespace amrex::literals;
 
 AMREX_GPU_HOST_DEVICE AMREX_INLINE
@@ -474,4 +480,7 @@ amrex::Real fast_pow(const amrex::Real x, const amrex::Real y) {
 
     return frac_ans * int_pow(x, n);
 }
+
+// NOLINTEND(modernize-use-std-numbers)
+
 #endif
diff --git a/util/microphysics_math.H b/util/microphysics_math.H
index 606981baea..c1646285f6 100644
--- a/util/microphysics_math.H
+++ b/util/microphysics_math.H
@@ -7,6 +7,8 @@
 // Override definitions of constants in GCEM to work around lack of support
 // on Power9 for long double. We don't use arithmetic > 64 bit anyway.
 
+// NOLINTBEGIN(modernize-use-std-numbers)
+
 #define GCEM_LOG_2 0.6931471805599453094172321214581765680755
 #define GCEM_LOG_10 2.3025850929940456840179914546843642076011
 #define GCEM_PI 3.1415926535897932384626433832795028841972
@@ -19,6 +21,8 @@
 #define GCEM_SQRT_HALF_PI 1.2533141373155002512078826424055226265035
 #define GCEM_E 2.7182818284590452353602874713526624977572
 
+// NOLINTEND(modernize-use-std-numbers)
+
 #include <gcem.hpp>
 
 #endif