use multiple angle formula in sneut to speed up sin/cos (#1379)

Co-authored-by: Eric T. Johnson <yut23@users.noreply.github.com>

neutrinos/sneut5.H

@@ -3,6 +3,7 @@
 
 #include <AMReX_REAL.H>
 #include <AMReX_Array.H>
+#include <AMReX_Math.H>
 
 using namespace amrex;
 
@@ -824,21 +825,28 @@ void nu_photo(const sneutf_t& sf,
 
     Real taudt = nu_constants::iln10 * sf.tempi;
 
-    // equation 3.7, compute the expensive trig functions only one time
+    // equation 3.7
 
-    Real cos1 = std::cos(nu_constants::fac1 * tau);
-    Real cos2 = std::cos(nu_constants::fac1 * 2.0e0_rt * tau);
-    Real cos3 = std::cos(nu_constants::fac1 * 3.0e0_rt * tau);
-    Real cos4 = std::cos(nu_constants::fac1 * 4.0e0_rt * tau);
-    Real cos5 = std::cos(nu_constants::fac1 * 5.0e0_rt * tau);
-    Real last = std::cos(nu_constants::fac2 * tau);
+    const auto [cos1, sin1] = amrex::Math::sincos(nu_constants::fac1 * tau);
 
-    Real sin1 = std::sin(nu_constants::fac1*tau);
-    Real sin2 = std::sin(nu_constants::fac1*2.0e0_rt*tau);
-    Real sin3 = std::sin(nu_constants::fac1*3.0e0_rt*tau);
-    Real sin4 = std::sin(nu_constants::fac1*4.0e0_rt*tau);
-    Real sin5 = std::sin(nu_constants::fac1*5.0e0_rt*tau);
-    Real xast = std::sin(nu_constants::fac2*tau);
+    // double, triple, etc. angle formulas
+    // sin/cos (2 fac1 tau)
+    const Real sin2 = 2.0_rt * sin1 * cos1;
+    const Real cos2 = 2.0_rt * cos1 * cos1 - 1.0_rt;
+
+    // sin/cos (3 fac1 tau)
+    const Real sin3 = sin1 * (3.0_rt - 4.0_rt * sin1 * sin1);
+    const Real cos3 = cos1 * (4.0_rt * cos1 * cos1 - 3.0_rt);
+
+    // sin/cos (4 fac1 tau) -- use double angle on sin2/cos2
+    const Real sin4 = 2.0_rt * sin2 * cos2;
+    const Real cos4 = 2.0_rt * cos2 * cos2 - 1.0_rt;
+
+    // sin/cos (5 fac1 tau)
+    const Real sin5 = sin1 * (5.0_rt - sin1 * sin1 * (20.0_rt - 16.0_rt * sin1 * sin1));
+    const Real cos5 = cos1 * (cos1 * cos1 * (16.0_rt * cos1 * cos1 - 20.0_rt) + 5.0_rt);
+
+    const auto [xast, last] = amrex::Math::sincos(nu_constants::fac2 * tau);
 
     Real a0 = 0.5e0_rt * c00
         + c01 * cos1 + dd01 * sin1 + c02 * cos2 + dd02 * sin2
@@ -1166,17 +1174,25 @@ void nu_brem(const sneutf_t& sf,
         //a0    = iln10*fac3*sf.deni;
 
         // compute the expensive trig functions of equation 5.21 only once
-        Real cos1 = std::cos(u);
-        Real cos2 = std::cos(2.0e0_rt*u);
-        Real cos3 = std::cos(3.0e0_rt*u);
-        Real cos4 = std::cos(4.0e0_rt*u);
-        Real cos5 = std::cos(5.0e0_rt*u);
-
-        Real sin1 = std::sin(u);
-        Real sin2 = std::sin(2.0e0_rt*u);
-        Real sin3 = std::sin(3.0e0_rt*u);
-        Real sin4 = std::sin(4.0e0_rt*u);
-        //sin5 = std::sin(5.0e0_rt*u);
+
+        const auto [cos1, sin1] = amrex::Math::sincos(u);
+
+        // double, triple, etc. angle formulas
+        // sin/cos (2 fac1 tau)
+        const Real sin2 = 2.0_rt * sin1 * cos1;
+        const Real cos2 = 2.0_rt * cos1 * cos1 - 1.0_rt;
+
+        // sin/cos (3 fac1 tau)
+        const Real sin3 = sin1 * (3.0_rt - 4.0_rt * sin1 * sin1);
+        const Real cos3 = cos1 * (4.0_rt * cos1 * cos1 - 3.0_rt);
+
+        // sin/cos (4 fac1 tau) -- use double angle on sin2/cos2
+        const Real sin4 = 2.0_rt * sin2 * cos2;
+        const Real cos4 = 2.0_rt * cos2 * cos2 - 1.0_rt;
+
+        // sin/cos (5 fac1 tau)
+        //const Real sin5 = sin1 * (5.0_rt - sin1 * sin1 * (20.0_rt - 16.0_rt * sin1 * sin1));
+        const Real cos5 = cos1 * (cos1 * cos1 * (16.0_rt * cos1 * cos1 - 20.0_rt) + 5.0_rt);
 
         // equation 5.21
         Real fb = 0.5e0_rt * 0.17946e0_rt + 0.00945e0_rt * u + 0.34529e0_rt

unit_test/burn_cell/ci-benchmarks/aprox13_QSS_unit_test.out

@@ -1,4 +1,4 @@
-AMReX (22.11-1-gc4a4811c373d) initialized
+AMReX (23.10-13-gecaa14456e3f) initialized
 starting the single zone burn...
 Maximum Time (s): 0.01
 State Density (g/cm^3): 1000000
@@ -32,41 +32,41 @@ RHS at t = 0
   Ni56 4.07828597e-29
 ------------------------------------
 successful? 1
- - Hnuc = 3.448624864e+19
- - added e = 3.448624864e+17
- - final T = 3173085674
+ - Hnuc = 3.448624858e+19
+ - added e = 3.448624858e+17
+ - final T = 3173085673
 ------------------------------------
 e initial = 1.389440599e+18
-e final =   1.734303086e+18
+e final =   1.734303085e+18
 ------------------------------------
 new mass fractions: 
-He4 0.7718339061
-C12 0.0001227479966
-O16 8.071653239e-06
-Ne20 3.560425108e-07
-Mg24 6.162605681e-07
-Si28 5.774651995e-06
-S32 5.291893286e-06
-Ar36 8.531081471e-06
-Ca40 0.0001668187948
-Ti44 0.0004124607789
-Cr48 0.000710075351
-Fe52 0.003115481694
+He4 0.771833906
+C12 0.0001227479967
+O16 8.071653242e-06
+Ne20 3.56042511e-07
+Mg24 6.162605686e-07
+Si28 5.774651999e-06
+S32 5.291893291e-06
+Ar36 8.531081479e-06
+Ca40 0.000166818795
+Ti44 0.0004124607795
+Cr48 0.0007100753521
+Fe52 0.003115481699
 Ni56 0.2236098678
 ------------------------------------
 species creation rates: 
-omegadot(He4): -22.81660939
-omegadot(C12): 0.01227479966
-omegadot(O16): 0.0008071653239
-omegadot(Ne20): 3.560425108e-05
-omegadot(Mg24): 6.162605681e-05
-omegadot(Si28): 0.0005774651995
-omegadot(S32): 0.0005291893286
-omegadot(Ar36): 0.0008531081471
-omegadot(Ca40): 0.01668187948
-omegadot(Ti44): 0.04124607789
-omegadot(Cr48): 0.0710075351
-omegadot(Fe52): 0.3115481694
+omegadot(He4): -22.8166094
+omegadot(C12): 0.01227479967
+omegadot(O16): 0.0008071653242
+omegadot(Ne20): 3.56042511e-05
+omegadot(Mg24): 6.162605686e-05
+omegadot(Si28): 0.0005774651999
+omegadot(S32): 0.0005291893291
+omegadot(Ar36): 0.0008531081479
+omegadot(Ca40): 0.0166818795
+omegadot(Ti44): 0.04124607795
+omegadot(Cr48): 0.07100753521
+omegadot(Fe52): 0.3115481699
 omegadot(Ni56): 22.36098678
 number of steps taken: 3695
-AMReX (22.11-1-gc4a4811c373d) finalized
+AMReX (23.10-13-gecaa14456e3f) finalized

unit_test/burn_cell/ci-benchmarks/aprox13_RKC_unit_test.out

@@ -1,4 +1,4 @@
-AMReX (23.06-6-g1e1c433b401c) initialized
+AMReX (23.10-13-gecaa14456e3f) initialized
 starting the single zone burn...
 Maximum Time (s): 0.01
 State Density (g/cm^3): 1000000
@@ -32,8 +32,8 @@ RHS at t = 0
   Ni56 3.938787868e-40
 ------------------------------------
 successful? 1
- - Hnuc = 8.558390239e+18
- - added e = 8.558390239e+16
+ - Hnuc = 8.558390237e+18
+ - added e = 8.558390237e+16
  - final T = 1516425860
 ------------------------------------
 e initial = 1.284393683e+17
@@ -42,31 +42,31 @@ e final =   2.140232707e+17
 new mass fractions: 
 He4 0.8760723967
 C12 0.1064099566
-O16 0.0001403204362
+O16 0.0001403204361
 Ne20 8.129701233e-05
 Mg24 0.0002972369574
 Si28 0.01113151728
-S32 0.005297014797
-Ar36 0.0005641483161
-Ca40 6.109847027e-06
-Ti44 2.004306621e-09
-Cr48 3.272104798e-14
-Fe52 7.12614193e-20
-Ni56 1.224820021e-26
+S32 0.005297014786
+Ar36 0.0005641483117
+Ca40 6.109846921e-06
+Ti44 2.004306558e-09
+Cr48 3.272104624e-14
+Fe52 7.126141319e-20
+Ni56 1.224819862e-26
 ------------------------------------
 species creation rates: 
 omegadot(He4): -12.39276033
 omegadot(C12): 10.64099566
-omegadot(O16): 0.01403204362
+omegadot(O16): 0.01403204361
 omegadot(Ne20): 0.008129701233
 omegadot(Mg24): 0.02972369574
 omegadot(Si28): 1.113151728
-omegadot(S32): 0.5297014797
-omegadot(Ar36): 0.05641483161
-omegadot(Ca40): 0.0006109847027
-omegadot(Ti44): 2.004306621e-07
-omegadot(Cr48): 3.272104798e-12
-omegadot(Fe52): 7.12614193e-18
-omegadot(Ni56): 1.224720021e-24
+omegadot(S32): 0.5297014786
+omegadot(Ar36): 0.05641483117
+omegadot(Ca40): 0.0006109846921
+omegadot(Ti44): 2.004306558e-07
+omegadot(Cr48): 3.272104624e-12
+omegadot(Fe52): 7.126141319e-18
+omegadot(Ni56): 1.224719862e-24
 number of steps taken: 255
-AMReX (23.06-6-g1e1c433b401c) finalized
+AMReX (23.10-13-gecaa14456e3f) finalized

unit_test/burn_cell/ci-benchmarks/aprox19_nse_unit_test.out

@@ -1,4 +1,4 @@
-AMReX (23.01-12-ga05384c0b7c4) initialized
+AMReX (23.10-13-gecaa14456e3f) initialized
 starting the single zone burn...
 reading the NSE table (C++) ...
 Maximum Time (s): 0.01
@@ -43,57 +43,57 @@ RHS at t = 0
   Ni56 -4.649371884e-24
      n 6.323463259e-06
      p 6.323463259e-06
-burn entered NSE during integration (after 280 steps), zone = (0, 0, 0)
+burn entered NSE during integration (after 407 steps), zone = (0, 0, 0)
 recovering burn failure in NSE, zone = (0, 0, 0)
 ------------------------------------
 successful? 1
- - Hnuc = 2.552804555e+21
- - added e = 2.552804555e+19
- - final T = 1.060964458e+10
+ - Hnuc = 2.577064219e+21
+ - added e = 2.577064219e+19
+ - final T = 1.063350508e+10
 ------------------------------------
 e initial = 1.516256085e+18
-e final =   2.704430164e+19
+e final =   2.728689827e+19
 ------------------------------------
 new mass fractions: 
-H1 0.0001057720428
-He3 1.320599831e-10
-He4 0.06412077048
-C12 1.016297715e-11
-N14 9.999991592e-31
-O16 9.999991592e-31
-Ne20 3.890552725e-17
-Mg24 9.999991592e-31
-Si28 1.194568144e-21
-S32 1.956720348e-25
-Ar36 4.189728368e-29
-Ca40 9.999991592e-31
-Ti44 3.531480561e-30
-Cr48 9.999991592e-31
-Fe52 9.999991592e-31
-Fe54 9.999991592e-31
-Ni56 9.999991592e-31
-n 0.4665061227
-p 0.4692673346
+H1 6.505608868e-05
+He3 1.038927825e-10
+He4 0.06022070651
+C12 9.99999158e-31
+N14 1.731352613e-13
+O16 2.450459695e-15
+Ne20 9.99999158e-31
+Mg24 6.790701604e-20
+Si28 3.428830286e-22
+S32 9.99999158e-31
+Ar36 9.99999158e-31
+Ca40 2.030520937e-29
+Ti44 9.99999158e-31
+Cr48 5.120465763e-30
+Fe52 9.99999158e-31
+Fe54 3.385938672e-28
+Ni56 9.99999158e-31
+n 0.4684052388
+p 0.4713089985
 ------------------------------------
 species creation rates: 
-omegadot(H1): -9.989422796
-omegadot(He3): -2.499999987
-omegadot(He4): -73.58792295
-omegadot(C12): -2.499999999
+omegadot(H1): -9.993494391
+omegadot(He3): -2.49999999
+omegadot(He4): -73.97792935
+omegadot(C12): -2.5
 omegadot(N14): -2.5
 omegadot(O16): -2.5
-omegadot(Ne20): 3.890552725e-15
-omegadot(Mg24): -8.408373504e-35
-omegadot(Si28): 1.194568143e-19
-omegadot(S32): 1.956710348e-23
-omegadot(Ar36): 4.089728368e-27
-omegadot(Ca40): -8.408373504e-35
-omegadot(Ti44): 2.531480561e-28
-omegadot(Cr48): -8.408373504e-35
-omegadot(Fe52): -8.408373504e-35
-omegadot(Fe54): -8.408373504e-35
-omegadot(Ni56): -8.408373504e-35
-omegadot(n): 46.65061227
-omegadot(p): 46.92673346
-number of steps taken: 27384
-AMReX (23.01-12-ga05384c0b7c4) finalized
+omegadot(Ne20): -8.420134623e-35
+omegadot(Mg24): 6.790701604e-18
+omegadot(Si28): 3.428830276e-20
+omegadot(S32): -8.420134623e-35
+omegadot(Ar36): -8.420134623e-35
+omegadot(Ca40): 1.930520937e-27
+omegadot(Ti44): -8.420134623e-35
+omegadot(Cr48): 4.120465763e-28
+omegadot(Fe52): -8.420134623e-35
+omegadot(Fe54): 3.375938672e-26
+omegadot(Ni56): -8.420134623e-35
+omegadot(n): 46.84052388
+omegadot(p): 47.13089985
+number of steps taken: 22341
+AMReX (23.10-13-gecaa14456e3f) finalized