pow(x, 1/3) to cbrt(x) (#1705)

AMReX-Astro · Jan 11, 2025 · 8106adc · 8106adc
1 parent 22c1543
commit 8106adc
Show file tree

Hide file tree

Showing 12 changed files with 195 additions and 197 deletions.
diff --git a/EOS/helmholtz/actual_eos.H b/EOS/helmholtz/actual_eos.H
@@ -686,8 +686,8 @@ void apply_coulomb_corrections (T& state)
     amrex::Real dsdd     = z * dxnidd;
     [[maybe_unused]] amrex::Real dsda     = z * dxnida;
 
-    amrex::Real lami     = 1.0e0_rt / std::pow(s, onethird);
-    amrex::Real inv_lami = 1.0e0_rt / lami;
+    amrex::Real inv_lami = std::cbrt(s);
+    amrex::Real lami     = 1.0e0_rt / inv_lami;
     z             = -onethird * lami;
     amrex::Real lamidd   = z * dsdd / s;
     [[maybe_unused]] amrex::Real lamida   = z * dsda / s;

diff --git a/EOS/ztwd/actual_eos.H b/EOS/ztwd/actual_eos.H
@@ -176,7 +176,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the pressure and energy:
 
-        x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+        x = std::cbrt(dens / B);
         pres = pressure(x);
         eint = enth - pres / dens;
 
@@ -188,7 +188,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the pressure, energy and enthalpy:
 
-        x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+        x = std::cbrt(dens / B);
         pres = pressure(x);
         enth = enthalpy(x, B);
         eint = enth - pres / dens;
@@ -203,7 +203,7 @@ void actual_eos (I input, T& state)
 
         pres_iter(pres, dens, B);
 
-        x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+        x = std::cbrt(dens / B);
         enth = enthalpy(x, B);
         eint = enth - pres / dens;
 
@@ -215,7 +215,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the enthalpy and energy:
 
-        x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+        x = std::cbrt(dens / B);
         enth = enthalpy(x, B);
         eint = enth - pres / dens;
 
@@ -227,7 +227,7 @@ void actual_eos (I input, T& state)
 
         // Solve for the pressure and enthalpy:
 
-        x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+        x = std::cbrt(dens / B);
         pres = pressure(x);
         enth = enthalpy(x, B);
 
@@ -241,7 +241,7 @@ void actual_eos (I input, T& state)
 
         pres_iter(pres, dens, B);
 
-        x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+        x = std::cbrt(dens / B);
         enth = enthalpy(x, B);
         eint = enth - pres / dens;
 
@@ -318,7 +318,7 @@ void actual_eos (I input, T& state)
 
     // Density derivatives are computed using the chain rule, e.g. dpdr = dpdx * dxdr.
 
-    x = std::pow(dens / B, 1.0_rt / 3.0_rt);
+    x = std::cbrt(dens / B);
     dxdr = (1.0_rt / 3.0_rt) * x / dens;
 
     Real dpdr = dxdr * dpdx(x);

diff --git a/conductivity/stellar/actual_conductivity.H b/conductivity/stellar/actual_conductivity.H
@@ -299,7 +299,7 @@ actual_conductivity (T& state)
   // from yakovlev & urpin soviet astro 1980 24 303 and
   // potekhin et al. 1997 aa 323 415 for degenerate regimes
   if (dlog10 > drel10) {
-    amrex::Real xmas = meff * std::pow(state.xne, third);
+    amrex::Real xmas = meff * std::cbrt(state.xne);
     amrex::Real ymas = std::sqrt(1.0_rt + xmas*xmas);
     amrex::Real wfac = weid * state.T/ymas * state.xne;
     amrex::Real cint = 1.0_rt;
@@ -311,10 +311,10 @@ actual_conductivity (T& state)
     // electron-electron collision frequency and thermal conductivity
     amrex::Real tpe = xec * std::sqrt(state.xne/ymas);
     amrex::Real yg = rt3 * tpe/state.T;
-    amrex::Real xrel = 1.009_rt * std::pow(zbar/abar * state.rho * 1.0e-6_rt, third);
+    amrex::Real xrel = 1.009_rt * std::cbrt(zbar/abar * state.rho * 1.0e-6_rt);
     amrex::Real beta2 = xrel*xrel/(1.0_rt + xrel*xrel);
-    amrex::Real jy = (1.0_rt + 6.0_rt/(5.0_rt*xrel*xrel) + 2.0_rt/(5.0_rt*xrel*xrel*xrel*xrel))
-      * ( yg*yg*yg / (3.0_rt * amrex::Math::powi<3>(1.0_rt + 0.07414_rt * yg))
+    amrex::Real jy = (1.0_rt + 6.0_rt/(5.0_rt*xrel*xrel) + 2.0_rt/(5.0_rt*amrex::Math::powi<4>(xrel)))
+        * ( amrex::Math::powi<3>(yg) / (3.0_rt * amrex::Math::powi<3>(1.0_rt + 0.07414_rt * yg))
           * std::log((2.81_rt - 0.810_rt*beta2 + yg)/yg)
           + std::pow(M_PI, 5/6.0_rt) * amrex::Math::powi<4>(yg/(13.91_rt + yg)));
     amrex::Real vee = 0.511_rt * state.T*state.T * xmas/(ymas*ymas) * std::sqrt(xmas/ymas) * jy;

diff --git a/integration/BackwardEuler/be_integrator.H b/integration/BackwardEuler/be_integrator.H
@@ -271,7 +271,7 @@ int be_integrator (BurnT& state, BeT& be)
             // can we potentially increase the timestep?
             // backward-Euler has a local truncation error of dt**2
 
-            amrex::Real dt_new = dt_sub * std::pow(1.0_rt / rel_error, 0.5_rt);
+            amrex::Real dt_new = dt_sub * std::sqrt(1.0_rt / rel_error);
             dt_sub = std::clamp(dt_new, dt_sub / 2.0, 2.0 * dt_sub);
 
         } else {

diff --git a/integration/RKC/rkc.H b/integration/RKC/rkc.H
@@ -126,8 +126,6 @@ int rkclow (BurnT& state, RkcT& rstate)
     const amrex::Real onep1 = 1.1_rt;
     const amrex::Real p4 = 0.4_rt;
     const amrex::Real p8 = 0.8_rt;
-    const amrex::Real one3rd = 1.0_rt/3.0_rt;
-    const amrex::Real two3rd = 2.0_rt/3.0_rt;
 
     // Initialize on the first call.
 
@@ -265,7 +263,7 @@ int rkclow (BurnT& state, RkcT& rstate)
         if (err > 1.0_rt) {
             // Step is rejected.
             rstate.nrejct++;
-            absh = p8 * absh / std::pow(err, one3rd);
+            absh = p8 * absh / std::cbrt(err);
             if (absh < hmin) {
                 return IERR_DT_UNDERFLOW;
             } else {
@@ -296,13 +294,13 @@ int rkclow (BurnT& state, RkcT& rstate)
         }
         amrex::Real fac = 10.0_rt;
         if (rstate.naccpt == 1) {
-            amrex::Real temp2 = std::pow(err, one3rd);
+            amrex::Real temp2 = std::cbrt(err);
             if (p8 < fac * temp2) {
                 fac = p8/temp2;
             }
         } else {
-            amrex::Real temp1 = p8 * absh * std::pow(errold, one3rd);
-            amrex::Real temp2 = std::abs(hold) * std::pow(err, two3rd);
+            amrex::Real temp1 = p8 * absh * std::cbrt(errold);
+            amrex::Real temp2 = std::abs(hold) * amrex::Math::powi<2>(std::cbrt(err));
             if (temp1 < fac * temp2) {
                 fac = temp1 / temp2;
             }

diff --git a/unit_test/burn_cell/ci-benchmarks/aprox13_RKC_unit_test.out b/unit_test/burn_cell/ci-benchmarks/aprox13_RKC_unit_test.out
@@ -1,5 +1,5 @@
-Initializing AMReX (23.08-206-g798ca3356d16)...
-AMReX (23.08-206-g798ca3356d16) initialized
+Initializing AMReX (23.05-658-g0165b6743355)...
+AMReX (23.05-658-g0165b6743355) initialized
 starting the single zone burn...
 Maximum Time (s): 0.01
 State Density (g/cm^3): 1000000
@@ -33,8 +33,8 @@ RHS at t = 0
   Ni56 3.938787868e-40
 ------------------------------------
 successful? 1
- - Hnuc = 8.55839024e+18
- - added e = 8.55839024e+16
+ - Hnuc = 8.558390237e+18
+ - added e = 8.558390237e+16
  - final T = 1516425860
 ------------------------------------
 e initial = 1.284393683e+17
@@ -43,31 +43,31 @@ e final =   2.140232707e+17
 new mass fractions: 
 He4 0.8760723967
 C12 0.1064099566
-O16 0.0001403204362
+O16 0.0001403204361
 Ne20 8.129701234e-05
-Mg24 0.0002972369573
+Mg24 0.0002972369574
 Si28 0.01113151728
-S32 0.005297014805
-Ar36 0.0005641483217
-Ca40 6.109847214e-06
-Ti44 2.004306753e-09
-Cr48 3.2721052e-14
-Fe52 7.126143425e-20
-Ni56 1.224820426e-26
+S32 0.005297014784
+Ar36 0.0005641483118
+Ca40 6.109846937e-06
+Ti44 2.004306574e-09
+Cr48 3.272104685e-14
+Fe52 7.126141581e-20
+Ni56 1.22481994e-26
 ------------------------------------
 species creation rates: 
 omegadot(He4): -12.39276033
 omegadot(C12): 10.64099566
-omegadot(O16): 0.01403204362
+omegadot(O16): 0.01403204361
 omegadot(Ne20): 0.008129701234
-omegadot(Mg24): 0.02972369573
+omegadot(Mg24): 0.02972369574
 omegadot(Si28): 1.113151728
-omegadot(S32): 0.5297014805
-omegadot(Ar36): 0.05641483217
-omegadot(Ca40): 0.0006109847214
-omegadot(Ti44): 2.004306753e-07
-omegadot(Cr48): 3.2721052e-12
-omegadot(Fe52): 7.126143424e-18
-omegadot(Ni56): 1.224720426e-24
+omegadot(S32): 0.5297014784
+omegadot(Ar36): 0.05641483118
+omegadot(Ca40): 0.0006109846937
+omegadot(Ti44): 2.004306574e-07
+omegadot(Cr48): 3.272104685e-12
+omegadot(Fe52): 7.126141581e-18
+omegadot(Ni56): 1.22471994e-24
 number of steps taken: 255
-AMReX (23.08-206-g798ca3356d16) finalized
+AMReX (23.05-658-g0165b6743355) finalized
diff --git a/unit_test/burn_cell/ci-benchmarks/chamulak_VODE_unit_test.out b/unit_test/burn_cell/ci-benchmarks/chamulak_VODE_unit_test.out
@@ -1,5 +1,5 @@
-Initializing AMReX (24.07-16-gdcb9cc0383dc)...
-AMReX (24.07-16-gdcb9cc0383dc) initialized
+Initializing AMReX (23.05-658-g0165b6743355)...
+AMReX (23.05-658-g0165b6743355) initialized
 starting the single zone burn...
 Maximum Time (s): 0.01585
 State Density (g/cm^3): 1000000000
@@ -13,21 +13,21 @@ RHS at t = 0
    ash 0.01230280576
 ------------------------------------
 successful? 1
- - Hnuc = 5.277400893e+17
- - added e = 8.364680415e+15
- - final T = 1433712612
+ - Hnuc = 5.277036926e+17
+ - added e = 8.364103527e+15
+ - final T = 1433684587
 ------------------------------------
 e initial = 1.253426044e+18
-e final =   1.261790725e+18
+e final =   1.261790148e+18
 ------------------------------------
 new mass fractions: 
-C12 0.9657895158
+C12 0.9657918752
 O16 1e-30
-ash 0.03421048417
+ash 0.03420812477
 ------------------------------------
 species creation rates: 
-omegadot(C12): -2.158390168
-omegadot(O16): 8.840999775e-44
-omegadot(ash): 2.158390168
+omegadot(C12): -2.15824131
+omegadot(O16): 7.735874803e-44
+omegadot(ash): 2.15824131
 number of steps taken: 381
-AMReX (24.07-16-gdcb9cc0383dc) finalized
+AMReX (23.05-658-g0165b6743355) finalized
diff --git a/unit_test/burn_cell/ci-benchmarks/ecsn_unit_test.out b/unit_test/burn_cell/ci-benchmarks/ecsn_unit_test.out
@@ -1,5 +1,5 @@
-Initializing AMReX (23.11-4-ga7afcba3cffd)...
-AMReX (23.11-4-ga7afcba3cffd) initialized
+Initializing AMReX (23.05-658-g0165b6743355)...
+AMReX (23.05-658-g0165b6743355) initialized
 starting the single zone burn...
 reading in network electron-capture / beta-decay tables...
 Maximum Time (s): 0.01
@@ -30,37 +30,37 @@ RHS at t = 0
    S32 1323519.174
 ------------------------------------
 successful? 1
- - Hnuc = 4.111466233e+19
- - added e = 4.111466233e+17
- - final T = 6536235789
+ - Hnuc = 4.520819292e+19
+ - added e = 4.520819292e+17
+ - final T = 6714483735
 ------------------------------------
 e initial = 5.995956082e+17
-e final =   1.010742231e+18
+e final =   1.051677537e+18
 ------------------------------------
 new mass fractions: 
-H1 0.007813581548
-He4 9.986703417e-31
-O16 7.812412175e-07
-O20 0.009524501615
-F20 0.00952319448
-Ne20 9.986703417e-31
-Mg24 8.489446236e-12
-Al27 9.986703417e-31
-Si28 0.2115906158
-P31 9.986703417e-31
-S32 0.7615473253
+H1 0.009195917106
+He4 9.999999989e-31
+O16 5.670739329e-07
+O20 0.009949514518
+F20 0.009948148298
+Ne20 9.076266252e-14
+Mg24 1.62767428e-09
+Al27 9.999999989e-31
+Si28 0.2969684506
+P31 9.999999989e-31
+S32 0.6739374007
 ------------------------------------
 species creation rates: 
-omegadot(H1): -0.2186418452
+omegadot(H1): -0.08040828943
 omegadot(He4): -3
-omegadot(O16): -49.99992188
-omegadot(O20): -0.04754983851
-omegadot(F20): -0.04768055199
+omegadot(O16): -49.99994329
+omegadot(O20): -0.005048548176
+omegadot(F20): -0.005185170248
 omegadot(Ne20): -30
-omegadot(Mg24): -9.999999999
+omegadot(Mg24): -9.999999837
 omegadot(Al27): -1
-omegadot(Si28): 20.15906158
+omegadot(Si28): 28.69684506
 omegadot(P31): -1
-omegadot(S32): 75.15473253
-number of steps taken: 9639
-AMReX (23.11-4-ga7afcba3cffd) finalized
+omegadot(S32): 66.39374007
+number of steps taken: 45595
+AMReX (23.05-658-g0165b6743355) finalized