From 448132f96834e1958e5fe38a10f08610ec665981 Mon Sep 17 00:00:00 2001 From: Zhi Chen <62574124+zhichen3@users.noreply.github.com> Date: Fri, 12 Jul 2024 07:30:45 -0400 Subject: [PATCH 1/5] remove test_nse and rename test_ase -> test_nse_net (#1610) remove redundant test_nse and rename test_ase to test_nse_net for better clarity --- .github/workflows/cuda.yml | 4 +- .github/workflows/nse_net.yml | 32 +++---- .github/workflows/nse_test.yml | 48 ----------- unit_test/test_nse/GNUmakefile | 42 ---------- unit_test/test_nse/Make.package | 2 - unit_test/test_nse/README.md | 7 -- unit_test/test_nse/_parameters | 17 ---- .../test_nse/ci-benchmarks/aprox21_ci.out | 54 ------------ unit_test/test_nse/inputs_aprox13 | 21 ----- unit_test/test_nse/inputs_aprox19 | 30 ------- unit_test/test_nse/inputs_aprox21 | 32 ------- unit_test/test_nse/inputs_iso7 | 16 ---- unit_test/test_nse/inputs_subch_full | 25 ------ unit_test/test_nse/main.cpp | 34 -------- unit_test/test_nse/nse_example.H | 83 ------------------- .../{test_ase => test_nse_net}/GNUmakefile | 0 .../{test_ase => test_nse_net}/Make.package | 0 .../{test_ase => test_nse_net}/README.md | 4 +- .../{test_ase => test_nse_net}/_parameters | 0 .../{test_ase => test_nse_net}/burn_cell.H | 0 .../ci-benchmarks/nse_net_unit_test.out} | 0 .../{test_ase => test_nse_net}/inputs_ase | 0 unit_test/{test_ase => test_nse_net}/main.cpp | 0 .../make_table/GNUmakefile | 0 .../make_table/Make.package | 0 .../make_table/README.md | 0 .../make_table/_parameters | 0 .../make_table/burn_cell.H | 0 .../nse_net_make_table_unit_test.out} | 0 .../make_table/main.cpp | 0 30 files changed, 20 insertions(+), 431 deletions(-) delete mode 100644 .github/workflows/nse_test.yml delete mode 100644 unit_test/test_nse/GNUmakefile delete mode 100644 unit_test/test_nse/Make.package delete mode 100644 unit_test/test_nse/README.md delete mode 100644 unit_test/test_nse/_parameters delete mode 100644 unit_test/test_nse/ci-benchmarks/aprox21_ci.out delete mode 100644 unit_test/test_nse/inputs_aprox13 delete mode 100644 unit_test/test_nse/inputs_aprox19 delete mode 100644 unit_test/test_nse/inputs_aprox21 delete mode 100644 unit_test/test_nse/inputs_iso7 delete mode 100644 unit_test/test_nse/inputs_subch_full delete mode 100644 unit_test/test_nse/main.cpp delete mode 100644 unit_test/test_nse/nse_example.H rename unit_test/{test_ase => test_nse_net}/GNUmakefile (100%) rename unit_test/{test_ase => test_nse_net}/Make.package (100%) rename unit_test/{test_ase => test_nse_net}/README.md (82%) rename unit_test/{test_ase => test_nse_net}/_parameters (100%) rename unit_test/{test_ase => test_nse_net}/burn_cell.H (100%) rename unit_test/{test_ase/ci-benchmarks/ase_nse_net_unit_test.out => test_nse_net/ci-benchmarks/nse_net_unit_test.out} (100%) rename unit_test/{test_ase => test_nse_net}/inputs_ase (100%) rename unit_test/{test_ase => test_nse_net}/main.cpp (100%) rename unit_test/{test_ase => test_nse_net}/make_table/GNUmakefile (100%) rename unit_test/{test_ase => test_nse_net}/make_table/Make.package (100%) rename unit_test/{test_ase => test_nse_net}/make_table/README.md (100%) rename unit_test/{test_ase => test_nse_net}/make_table/_parameters (100%) rename unit_test/{test_ase => test_nse_net}/make_table/burn_cell.H (100%) rename unit_test/{test_ase/make_table/ci-benchmarks/ase_nse_net_make_table_unit_test.out => test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out} (100%) rename unit_test/{test_ase => test_nse_net}/make_table/main.cpp (100%) diff --git a/.github/workflows/cuda.yml b/.github/workflows/cuda.yml index 60aff9d0e8..7e418249a8 100644 --- a/.github/workflows/cuda.yml +++ b/.github/workflows/cuda.yml @@ -41,8 +41,8 @@ jobs: make realclean make NETWORK_DIR=ignition_reaclib/URCA-simple USE_CUDA=TRUE COMP=gnu USE_MPI=FALSE -j 2 - - name: compile test_ase (ase) + - name: compile test_nse_net (ase) run: | export PATH=/usr/local/nvidia/bin:/usr/local/cuda/bin:${PATH} - cd unit_test/test_ase + cd unit_test/test_nse_net make USE_CUDA=TRUE COMP=gnu USE_MPI=FALSE -j 2 diff --git a/.github/workflows/nse_net.yml b/.github/workflows/nse_net.yml index e6e8741724..fc1124fe61 100644 --- a/.github/workflows/nse_net.yml +++ b/.github/workflows/nse_net.yml @@ -26,42 +26,42 @@ jobs: sudo apt-get update -y -qq sudo apt-get -qq -y install curl cmake jq clang g++>=9.3.0 - - name: Compile, test_ase (NSE_NET, ase) + - name: Compile, test_nse_net (NSE_NET, ase) run: | - cd unit_test/test_ase + cd unit_test/test_nse_net make realclean make -j 4 - - name: Run test_ase (NSE_NET, ase) + - name: Run test_nse_net (NSE_NET, ase) run: | - cd unit_test/test_ase + cd unit_test/test_nse_net ./main3d.gnu.ex inputs_ase amrex.fpe_trap_{invalid,zero,overflow}=1 > test.out - name: Print backtrace - if: ${{ failure() && hashFiles('unit_test/test_ase/Backtrace.0') != '' }} - run: cat unit_test/test_ase/Backtrace.0 + if: ${{ failure() && hashFiles('unit_test/test_nse_net/Backtrace.0') != '' }} + run: cat unit_test/test_nse_net/Backtrace.0 - name: Compare to stored output (NSE_NET, ase) run: | - cd unit_test/test_ase - diff -I "^Initializing AMReX" -I "^AMReX" -I "^reading in reaclib rates" test.out ci-benchmarks/ase_nse_net_unit_test.out + cd unit_test/test_nse_net + diff -I "^Initializing AMReX" -I "^AMReX" -I "^reading in reaclib rates" test.out ci-benchmarks/nse_net_unit_test.out - - name: Compile, test_ase/make_table (NSE_NET, ase, make_table) + - name: Compile, test_nse_net/make_table (NSE_NET, ase, make_table) run: | - cd unit_test/test_ase/make_table + cd unit_test/test_nse_net/make_table make realclean make -j 4 - - name: Run, test_ase/make_table (NSE_NET, ase, make_table) + - name: Run, test_nse_net/make_table (NSE_NET, ase, make_table) run: | - cd unit_test/test_ase/make_table + cd unit_test/test_nse_net/make_table ./main3d.gnu.ex amrex.fpe_trap_{invalid,zero,overflow}=1 > test.out - name: Print backtrace - if: ${{ failure() && hashFiles('unit_test/test_ase/make_table/Backtrace.0') != '' }} - run: cat unit_test/test_ase/make_table/Backtrace.0 + if: ${{ failure() && hashFiles('unit_test/test_nse_net/make_table/Backtrace.0') != '' }} + run: cat unit_test/test_nse_net/make_table/Backtrace.0 - name: Compare to stored output (NSE_NET, ase, make_table) run: | - cd unit_test/test_ase/make_table - diff -I "^Initializing AMReX" -I "^AMReX" -I "^reading in reaclib rates" test.out ci-benchmarks/ase_nse_net_make_table_unit_test.out + cd unit_test/test_nse_net/make_table + diff -I "^Initializing AMReX" -I "^AMReX" -I "^reading in reaclib rates" test.out ci-benchmarks/nse_net_make_table_unit_test.out diff --git a/.github/workflows/nse_test.yml b/.github/workflows/nse_test.yml deleted file mode 100644 index 1b898f0e8c..0000000000 --- a/.github/workflows/nse_test.yml +++ /dev/null @@ -1,48 +0,0 @@ -name: nse_test - -on: [pull_request] -jobs: - nse_test: - runs-on: ubuntu-latest - steps: - - uses: actions/checkout@v4 - with: - fetch-depth: 0 - - - name: Get AMReX - run: | - mkdir external - cd external - git clone https://github.com/AMReX-Codes/amrex.git - cd amrex - git checkout development - echo 'AMREX_HOME=$(GITHUB_WORKSPACE)/external/amrex' >> $GITHUB_ENV - echo $AMREX_HOME - if [[ -n "${AMREX_HOME}" ]]; then exit 1; fi - cd ../.. - - - name: Install dependencies - run: | - sudo apt-get update -y -qq - sudo apt-get -qq -y install curl cmake jq clang g++>=9.3.0 - - - name: Compile, test_nse (aprox21) - run: | - cd unit_test/test_nse - make realclean - make NETWORK_DIR=aprox21 -j 4 - - - name: Run test_nse (aprox21) - run: | - cd unit_test/test_nse - ./main3d.gnu.ex inputs_aprox21 amrex.fpe_trap_{invalid,zero,overflow}=1 > test.out - - - name: Print backtrace - if: ${{ failure() && hashFiles('unit_test/test_nse/Backtrace.0') != '' }} - run: cat unit_test/test_nse/Backtrace.0 - - - name: Compare to stored output (aprox21) - run: | - cd unit_test/test_nse - diff -I "^Initializing AMReX" -I "^AMReX" -I "^reading in reaclib rates" -I "^0x" test.out ci-benchmarks/aprox21_ci.out - diff --git a/unit_test/test_nse/GNUmakefile b/unit_test/test_nse/GNUmakefile deleted file mode 100644 index be47b8580a..0000000000 --- a/unit_test/test_nse/GNUmakefile +++ /dev/null @@ -1,42 +0,0 @@ -PRECISION = DOUBLE -PROFILE = FALSE - -DEBUG = FALSE - -DIM = 3 - -COMP = gnu - -USE_MPI = FALSE -USE_OMP = FALSE - -USE_REACT = TRUE - -USE_NSE_NET = TRUE -EBASE = main - -BL_NO_FORT = TRUE - -# define the location of the Microphysics top directory -MICROPHYSICS_HOME := ../.. - -# This sets the EOS directory -EOS_DIR := helmholtz - -# This sets the network directory -NETWORK_DIR := aprox21 - -SCREEN_METHOD := chabrier1998 - -USE_SIMPLIFIED_SDC = TRUE - -CONDUCTIVITY_DIR := stellar - -INTEGRATOR_DIR = VODE - -EXTERN_SEARCH += . .. - -Bpack := ./Make.package -Blocs := . - -include $(MICROPHYSICS_HOME)/unit_test/Make.unit_test diff --git a/unit_test/test_nse/Make.package b/unit_test/test_nse/Make.package deleted file mode 100644 index 6bc6cd127c..0000000000 --- a/unit_test/test_nse/Make.package +++ /dev/null @@ -1,2 +0,0 @@ -CEXE_sources += main.cpp -CEXE_headers += nse_example.H diff --git a/unit_test/test_nse/README.md b/unit_test/test_nse/README.md deleted file mode 100644 index 27f809bb53..0000000000 --- a/unit_test/test_nse/README.md +++ /dev/null @@ -1,7 +0,0 @@ -# nse_example - -`nse_example` finds the NSE mass fractions of a given state. -The density, temperature, and composition are set in the -inputs file. - -Upon completion, the new state is printed to the screen. \ No newline at end of file diff --git a/unit_test/test_nse/_parameters b/unit_test/test_nse/_parameters deleted file mode 100644 index 355960df5e..0000000000 --- a/unit_test/test_nse/_parameters +++ /dev/null @@ -1,17 +0,0 @@ -@namespace: unit_test - -run_prefix string "" - -small_temp real 1.e5 -small_dens real 1.e5 - -density real 1.e7 - -temperature real 3.e9 - -mu_p real -3.0 -mu_n real -12.0 - -y_e real 0.5 - - diff --git a/unit_test/test_nse/ci-benchmarks/aprox21_ci.out b/unit_test/test_nse/ci-benchmarks/aprox21_ci.out deleted file mode 100644 index bb3bd62e2e..0000000000 --- a/unit_test/test_nse/ci-benchmarks/aprox21_ci.out +++ /dev/null @@ -1,54 +0,0 @@ -Initializing AMReX (23.07-395-ge6c93bf22695)... -AMReX (23.07-395-ge6c93bf22695) initialized -starting the single zone burn... -0x7ffcbd574c90 -State Density (g/cm^3): 277355338.4 -State Temperature (K): 5197769252 -Mass Fraction (H1): 0.4 -Mass Fraction (He3): 0.2 -Mass Fraction (He4): 0.2 -Mass Fraction (C12): 0.01111111111 -Mass Fraction (N14): 0.01111111111 -Mass Fraction (O16): 0.01111111111 -Mass Fraction (Ne20): 0.01111111111 -Mass Fraction (Mg24): 0.01111111111 -Mass Fraction (Si28): 0.01111111111 -Mass Fraction (S32): 0.01111111111 -Mass Fraction (Ar36): 0.01111111111 -Mass Fraction (Ca40): 0.01111111111 -Mass Fraction (Ti44): 0.01111111111 -Mass Fraction (Cr48): 0.01111111111 -Mass Fraction (Cr56): 0.01111111111 -Mass Fraction (Fe52): 0.01111111111 -Mass Fraction (Fe54): 0.01111111111 -Mass Fraction (Fe56): 0.01111111111 -Mass Fraction (Ni56): 0.01111111111 -Mass Fraction (n): 0.01111111111 -Mass Fraction (p): 0.01111111111 -electron fraction is 0.5 -After solving: -chemical potential of proton is -5.581982209 -chemical potential of neutron is -12.0247581 -NSE state: -H1 : 0 -He3 : 6.726596351e-12 -He4 : 0.002680193661 -C12 : 2.088656717e-08 -N14 : 3.63235298e-13 -O16 : 8.472970972e-08 -Ne20 : 1.346294667e-09 -Mg24 : 5.596482124e-07 -Si28 : 0.0010080397 -S32 : 0.002046314471 -Ar36 : 0.00209146282 -Ca40 : 0.005238622021 -Ti44 : 0.0001839172934 -Cr48 : 0.002037708642 -Cr56 : 7.745054245e-22 -Fe52 : 0.03983693988 -Fe54 : 0.04496220141 -Fe56 : 1.685809845e-05 -Ni56 : 0.8982306026 -n : 9.441707272e-10 -p : 0.001666471811 -AMReX (23.07-395-ge6c93bf22695) finalized diff --git a/unit_test/test_nse/inputs_aprox13 b/unit_test/test_nse/inputs_aprox13 deleted file mode 100644 index 1f6020da45..0000000000 --- a/unit_test/test_nse/inputs_aprox13 +++ /dev/null @@ -1,21 +0,0 @@ -unit_test.small_temp = 1.e5 -unit_test.small_dens = 1.e5 - -unit_test.density = 1.e6 -unit_test.temperature = 3.e9 - -unit_test.X1 = 1.0 -unit_test.X2 = 0.0 -unit_test.X3 = 0.0 -unit_test.X4 = 0.0 -unit_test.X5 = 0.0 -unit_test.X6 = 0.0 -unit_test.X7 = 0.0 -unit_test.X8 = 0.0 -unit_test.X9 = 0.0 -unit_test.X10 = 0.0 -unit_test.X11 = 0.0 -unit_test.X12 = 0.0 -unit_test.X13 = 0.0 - -unit_test.y_e = 0.5 diff --git a/unit_test/test_nse/inputs_aprox19 b/unit_test/test_nse/inputs_aprox19 deleted file mode 100644 index 3ac62abc70..0000000000 --- a/unit_test/test_nse/inputs_aprox19 +++ /dev/null @@ -1,30 +0,0 @@ -unit_test.small_temp = 1.e5 -unit_test.small_dens = 1.e5 - -unit_test.density = 1.e6 -unit_test.temperature = 3.e9 - -unit_test.X1 = 0.7 -unit_test.X2 = 0.025 -unit_test.X3 = 0.2 -unit_test.X4 = 0.025 -unit_test.X5 = 0.025 -unit_test.X6 = 0.025 -unit_test.X7 = 0.0 -unit_test.X8 = 0.0 -unit_test.X9 = 0.0 -unit_test.X10 = 0.0 -unit_test.X11 = 0.0 -unit_test.X12 = 0.0 -unit_test.X13 = 0.0 -unit_test.X14 = 0.0 -unit_test.X15 = 0.0 -unit_test.X16 = 0.0 -unit_test.X17 = 0.0 -unit_test.X18 = 0.0 -unit_test.X19 = 0.0 - -unit_test.mu_p = -3.0 -unit_test.mu_n = -12.0 - -unit_test.y_e = 0.5 diff --git a/unit_test/test_nse/inputs_aprox21 b/unit_test/test_nse/inputs_aprox21 deleted file mode 100644 index 480b74ca54..0000000000 --- a/unit_test/test_nse/inputs_aprox21 +++ /dev/null @@ -1,32 +0,0 @@ -unit_test.small_temp = 1e5 -unit_test.small_dens = 1e5 - -unit_test.density = 277355338.38976681 -unit_test.temperature = 5197769252.1761198 - -unit_test.X1 = 0.4 -unit_test.X2 = 0.2 -unit_test.X3 = 0.2 -unit_test.X4 = 0.01111111111 -unit_test.X5 = 0.01111111111 -unit_test.X6 = 0.01111111111 -unit_test.X7 = 0.01111111111 -unit_test.X8 = 0.01111111111 -unit_test.X9 = 0.01111111111 -unit_test.X10 = 0.01111111111 -unit_test.X11 = 0.01111111111 -unit_test.X12 = 0.01111111111 -unit_test.X13 = 0.01111111111 -unit_test.X14 = 0.01111111111 -unit_test.X15 = 0.01111111111 -unit_test.X16 = 0.01111111111 -unit_test.X17 = 0.01111111111 -unit_test.X18 = 0.01111111111 -unit_test.X19 = 0.01111111111 -unit_test.X20 = 0.01111111111 -unit_test.X21 = 0.01111111111 - -unit_test.mu_p = -3.0 -unit_test.mu_n = -12.0 - -unit_test.y_e = 0.5 diff --git a/unit_test/test_nse/inputs_iso7 b/unit_test/test_nse/inputs_iso7 deleted file mode 100644 index 45d1a6c7cc..0000000000 --- a/unit_test/test_nse/inputs_iso7 +++ /dev/null @@ -1,16 +0,0 @@ -unit_test.small_temp = 1e5 -unit_test.small_dens = 1e5 - -unit_test.density = 1.e6 -unit_test.temperature = 3.e9 - -unit_test.X1 = 1.0 -unit_test.X2 = 0.0 -unit_test.X3 = 0.0 -unit_test.X4 = 0.0 -unit_test.X5 = 0.0 -unit_test.X6 = 0.0 -unit_test.X7 = 0.0 - -unit_test.y_e = 0.5 - diff --git a/unit_test/test_nse/inputs_subch_full b/unit_test/test_nse/inputs_subch_full deleted file mode 100644 index 171b663f4e..0000000000 --- a/unit_test/test_nse/inputs_subch_full +++ /dev/null @@ -1,25 +0,0 @@ -unit_test.small_temp = 1.e5 -unit_test.small_dens = 1.e5 - -unit_test.density = 1.e6 -unit_test.temperature = 3.e9 - -unit_test.X1 = 0.5 -unit_test.X2 = 0.1 -unit_test.X3 = 0.1 -unit_test.X4 = 0.1 -unit_test.X5 = 0.1 -unit_test.X6 = 0.1 -unit_test.X7 = 0.0 -unit_test.X8 = 0.0 -unit_test.X9 = 0.0 -unit_test.X10 = 0.0 -unit_test.X11 = 0.0 -unit_test.X12 = 0.0 -unit_test.X13 = 0.0 - -unit_test.mu_p = -3.0 -unit_test.mu_n = -12.0 - -unit_test.y_e = 0.5 - diff --git a/unit_test/test_nse/main.cpp b/unit_test/test_nse/main.cpp deleted file mode 100644 index 04ed4eec82..0000000000 --- a/unit_test/test_nse/main.cpp +++ /dev/null @@ -1,34 +0,0 @@ -#include -#include -#include - -#include -#include -using namespace amrex; - -#include -#include -#include -#include -#include - -int main(int argc, char *argv[]) { - - amrex::Initialize(argc, argv); - - std::cout << "starting the single zone burn..." << std::endl; - - ParmParse ppa("amr"); - - init_unit_test(); - - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) - eos_init(small_temp, small_dens); - - // C++ Network, RHS, screening, rates initialization - network_init(); - - nse_example_c(); - - amrex::Finalize(); -} diff --git a/unit_test/test_nse/nse_example.H b/unit_test/test_nse/nse_example.H deleted file mode 100644 index 5a712485d2..0000000000 --- a/unit_test/test_nse/nse_example.H +++ /dev/null @@ -1,83 +0,0 @@ -#ifndef NSE_EXAMPLE_H -#define NSE_EXAMPLE_H - -#include -#include -#include -#include -#include -#include -#include -#include - -using namespace unit_test_rp; - -AMREX_INLINE -void nse_example_c() -{ - - burn_t state; - - amrex::Real massfractions[NumSpec] = {-1.0}; - std::cout << massfractions << std::endl; - // Make sure user set all the mass fractions to values in the interval [0, 1] - for (int n = 1; n <= NumSpec; ++n) { - - massfractions[n-1] = get_xn(n); - - if (massfractions[n-1] < 0 || massfractions[n-1] > 1) { - amrex::Error("mass fraction for " + short_spec_names_cxx[n-1] + " not initialized in the interval [0,1]!"); - } - - } - - // Echo initial conditions at burn and fill burn state input - - std::cout << "State Density (g/cm^3): " << density << std::endl; - std::cout << "State Temperature (K): " << temperature << std::endl; - for (int n = 0; n < NumSpec; ++n) { - std::cout << "Mass Fraction (" << short_spec_names_cxx[n] << "): " << massfractions[n] << std::endl; - } - - state.T = temperature; - state.rho = density; - for (int n = 0; n < NumSpec; ++n) { - state.xn[n] = massfractions[n]; - } - - - // normalize -- just in case - - normalize_abundances_burn(state); - - // compute the initial Ye - - state.y_e = y_e; - - // composition(state); - - std::cout << "electron fraction is " << state.y_e << std::endl; - - // set initial chemical potential of proton and neutron - - state.mu_p = mu_p; - state.mu_n = mu_n; - - const bool assume_ye_valid = true; - amrex::Real eps = 1.0e-10_rt; - - // find the nse state - use_hybrid_solver = 1; - - auto NSE_STATE = get_actual_nse_state(state, eps, assume_ye_valid); - - std::cout << "After solving: " << std::endl; - std::cout << "chemical potential of proton is " << state.mu_p << std::endl; - std::cout << "chemical potential of neutron is " << state.mu_n << std::endl; - - std::cout << "NSE state: " << std::endl; - for (int n = 0; n < NumSpec; ++n) { - std::cout << short_spec_names_cxx[n] << " : " << NSE_STATE.xn[n]<< std::endl; - } -} -#endif diff --git a/unit_test/test_ase/GNUmakefile b/unit_test/test_nse_net/GNUmakefile similarity index 100% rename from unit_test/test_ase/GNUmakefile rename to unit_test/test_nse_net/GNUmakefile diff --git a/unit_test/test_ase/Make.package b/unit_test/test_nse_net/Make.package similarity index 100% rename from unit_test/test_ase/Make.package rename to unit_test/test_nse_net/Make.package diff --git a/unit_test/test_ase/README.md b/unit_test/test_nse_net/README.md similarity index 82% rename from unit_test/test_ase/README.md rename to unit_test/test_nse_net/README.md index f7069ddcdd..098db51823 100644 --- a/unit_test/test_ase/README.md +++ b/unit_test/test_nse_net/README.md @@ -1,6 +1,6 @@ -# test_ase +# test_nse_net -`test_ase` finds the NSE mass fractions of a given state. +`test_nse_net` finds the NSE mass fractions of a given state. The density, temperature, and composition are set in the inputs file. Then we update the NSE mass fraction to the current state to make sure we're in NSE. Then we diff --git a/unit_test/test_ase/_parameters b/unit_test/test_nse_net/_parameters similarity index 100% rename from unit_test/test_ase/_parameters rename to unit_test/test_nse_net/_parameters diff --git a/unit_test/test_ase/burn_cell.H b/unit_test/test_nse_net/burn_cell.H similarity index 100% rename from unit_test/test_ase/burn_cell.H rename to unit_test/test_nse_net/burn_cell.H diff --git a/unit_test/test_ase/ci-benchmarks/ase_nse_net_unit_test.out b/unit_test/test_nse_net/ci-benchmarks/nse_net_unit_test.out similarity index 100% rename from unit_test/test_ase/ci-benchmarks/ase_nse_net_unit_test.out rename to unit_test/test_nse_net/ci-benchmarks/nse_net_unit_test.out diff --git a/unit_test/test_ase/inputs_ase b/unit_test/test_nse_net/inputs_ase similarity index 100% rename from unit_test/test_ase/inputs_ase rename to unit_test/test_nse_net/inputs_ase diff --git a/unit_test/test_ase/main.cpp b/unit_test/test_nse_net/main.cpp similarity index 100% rename from unit_test/test_ase/main.cpp rename to unit_test/test_nse_net/main.cpp diff --git a/unit_test/test_ase/make_table/GNUmakefile b/unit_test/test_nse_net/make_table/GNUmakefile similarity index 100% rename from unit_test/test_ase/make_table/GNUmakefile rename to unit_test/test_nse_net/make_table/GNUmakefile diff --git a/unit_test/test_ase/make_table/Make.package b/unit_test/test_nse_net/make_table/Make.package similarity index 100% rename from unit_test/test_ase/make_table/Make.package rename to unit_test/test_nse_net/make_table/Make.package diff --git a/unit_test/test_ase/make_table/README.md b/unit_test/test_nse_net/make_table/README.md similarity index 100% rename from unit_test/test_ase/make_table/README.md rename to unit_test/test_nse_net/make_table/README.md diff --git a/unit_test/test_ase/make_table/_parameters b/unit_test/test_nse_net/make_table/_parameters similarity index 100% rename from unit_test/test_ase/make_table/_parameters rename to unit_test/test_nse_net/make_table/_parameters diff --git a/unit_test/test_ase/make_table/burn_cell.H b/unit_test/test_nse_net/make_table/burn_cell.H similarity index 100% rename from unit_test/test_ase/make_table/burn_cell.H rename to unit_test/test_nse_net/make_table/burn_cell.H diff --git a/unit_test/test_ase/make_table/ci-benchmarks/ase_nse_net_make_table_unit_test.out b/unit_test/test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out similarity index 100% rename from unit_test/test_ase/make_table/ci-benchmarks/ase_nse_net_make_table_unit_test.out rename to unit_test/test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out diff --git a/unit_test/test_ase/make_table/main.cpp b/unit_test/test_nse_net/make_table/main.cpp similarity index 100% rename from unit_test/test_ase/make_table/main.cpp rename to unit_test/test_nse_net/make_table/main.cpp From b413290def9550c8fac7035ec614f90f68c12ba4 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Fri, 12 Jul 2024 11:34:44 -0400 Subject: [PATCH 2/5] remove some more Fortran comments (#1609) --- Make.Microphysics | 4 ---- nse_solver/make_table/main.cpp | 2 +- nse_solver/nse_compatibility/main.cpp | 2 +- unit_test/burn_cell/main.cpp | 2 +- unit_test/burn_cell_primordial_chem/main.cpp | 2 +- unit_test/burn_cell_sdc/main.cpp | 2 +- unit_test/eos_cell/main.cpp | 2 +- unit_test/nse_table_cell/main.cpp | 2 +- unit_test/test_aprox_rates/README.md | 3 +-- unit_test/test_eos/README.md | 4 +--- unit_test/test_jac/main.cpp | 2 +- unit_test/test_linear_algebra/main.cpp | 2 +- unit_test/test_nse_interp/main.cpp | 2 +- unit_test/test_nse_net/main.cpp | 2 +- unit_test/test_nse_net/make_table/main.cpp | 2 +- unit_test/test_part_func/main.cpp | 2 +- unit_test/test_react/main.cpp | 2 +- unit_test/test_rhs/main.cpp | 2 +- unit_test/test_sdc/main.cpp | 2 +- unit_test/test_sdc_vode_rhs/main.cpp | 2 +- 20 files changed, 19 insertions(+), 26 deletions(-) diff --git a/Make.Microphysics b/Make.Microphysics index a6b7ece9de..9b6fa9e526 100644 --- a/Make.Microphysics +++ b/Make.Microphysics @@ -118,10 +118,6 @@ Blocs += $(foreach dir, $(EXTERN_CORE), $(dir)) include $(Bpack) -# this is a safety from the mega-fortran attempts -f90EXE_sources += $(ca_f90EXE_sources) -F90EXE_sources += $(ca_F90EXE_sources) - INCLUDE_LOCATIONS += $(Blocs) VPATH_LOCATIONS += $(Blocs) diff --git a/nse_solver/make_table/main.cpp b/nse_solver/make_table/main.cpp index 7295b50c85..4ac7112bf1 100644 --- a/nse_solver/make_table/main.cpp +++ b/nse_solver/make_table/main.cpp @@ -22,7 +22,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(); // C++ Network, RHS, screening, rates initialization diff --git a/nse_solver/nse_compatibility/main.cpp b/nse_solver/nse_compatibility/main.cpp index 2581910048..a5982d52bf 100644 --- a/nse_solver/nse_compatibility/main.cpp +++ b/nse_solver/nse_compatibility/main.cpp @@ -22,7 +22,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/burn_cell/main.cpp b/unit_test/burn_cell/main.cpp index 3b401b12ea..faffbc210b 100644 --- a/unit_test/burn_cell/main.cpp +++ b/unit_test/burn_cell/main.cpp @@ -16,7 +16,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/burn_cell_primordial_chem/main.cpp b/unit_test/burn_cell_primordial_chem/main.cpp index 9a08bec7e5..5e458f8936 100644 --- a/unit_test/burn_cell_primordial_chem/main.cpp +++ b/unit_test/burn_cell_primordial_chem/main.cpp @@ -32,7 +32,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and + // C++ EOS initialization (must be done after // init_extern_parameters) eos_init(unit_test_rp::small_temp, unit_test_rp::small_dens); diff --git a/unit_test/burn_cell_sdc/main.cpp b/unit_test/burn_cell_sdc/main.cpp index 3b401b12ea..faffbc210b 100644 --- a/unit_test/burn_cell_sdc/main.cpp +++ b/unit_test/burn_cell_sdc/main.cpp @@ -16,7 +16,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/eos_cell/main.cpp b/unit_test/eos_cell/main.cpp index c7b0a86fa5..b66f0c18d0 100644 --- a/unit_test/eos_cell/main.cpp +++ b/unit_test/eos_cell/main.cpp @@ -16,7 +16,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/nse_table_cell/main.cpp b/unit_test/nse_table_cell/main.cpp index b97e628012..166a2c2b4a 100644 --- a/unit_test/nse_table_cell/main.cpp +++ b/unit_test/nse_table_cell/main.cpp @@ -14,7 +14,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_aprox_rates/README.md b/unit_test/test_aprox_rates/README.md index c22ea9d000..be323d5d6d 100644 --- a/unit_test/test_aprox_rates/README.md +++ b/unit_test/test_aprox_rates/README.md @@ -1,5 +1,4 @@ # test_aprox_rates -This simply tests the aprox rates (without integration). It tests -both the Fortran and C++ implementations (via the do_cxx parameters). +This simply tests the aprox rates (without integration). diff --git a/unit_test/test_eos/README.md b/unit_test/test_eos/README.md index 2a4294bafd..3e97ea194c 100644 --- a/unit_test/test_eos/README.md +++ b/unit_test/test_eos/README.md @@ -1,8 +1,6 @@ # test_eos -This is the unit test for testing the equation of state. This tests -both the Fortran and C++ implementations (according to the parameter -`do_cxx`). +This is the unit test for testing the equation of state. The test sets up a cube of data (density on one axis, temperature on another, and composition on the third) and calls the EOS in various diff --git a/unit_test/test_jac/main.cpp b/unit_test/test_jac/main.cpp index f43ff9e9a5..e9ca42f556 100644 --- a/unit_test/test_jac/main.cpp +++ b/unit_test/test_jac/main.cpp @@ -98,7 +98,7 @@ void main_main () init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_linear_algebra/main.cpp b/unit_test/test_linear_algebra/main.cpp index f1941303e0..e16897053f 100644 --- a/unit_test/test_linear_algebra/main.cpp +++ b/unit_test/test_linear_algebra/main.cpp @@ -14,7 +14,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_nse_interp/main.cpp b/unit_test/test_nse_interp/main.cpp index 56ad8e7488..9de0af0faf 100644 --- a/unit_test/test_nse_interp/main.cpp +++ b/unit_test/test_nse_interp/main.cpp @@ -14,7 +14,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(unit_test_rp::small_temp, unit_test_rp::small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_nse_net/main.cpp b/unit_test/test_nse_net/main.cpp index f172e2fe10..134b08b579 100644 --- a/unit_test/test_nse_net/main.cpp +++ b/unit_test/test_nse_net/main.cpp @@ -22,7 +22,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_nse_net/make_table/main.cpp b/unit_test/test_nse_net/make_table/main.cpp index f172e2fe10..134b08b579 100644 --- a/unit_test/test_nse_net/make_table/main.cpp +++ b/unit_test/test_nse_net/make_table/main.cpp @@ -22,7 +22,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_part_func/main.cpp b/unit_test/test_part_func/main.cpp index 4f16284455..9fa048d804 100644 --- a/unit_test/test_part_func/main.cpp +++ b/unit_test/test_part_func/main.cpp @@ -22,7 +22,7 @@ int main(int argc, char *argv[]) { init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_react/main.cpp b/unit_test/test_react/main.cpp index ef051bf564..0a47397c39 100644 --- a/unit_test/test_react/main.cpp +++ b/unit_test/test_react/main.cpp @@ -104,7 +104,7 @@ void main_main () init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_rhs/main.cpp b/unit_test/test_rhs/main.cpp index 1b04270667..4c08451bb3 100644 --- a/unit_test/test_rhs/main.cpp +++ b/unit_test/test_rhs/main.cpp @@ -106,7 +106,7 @@ void main_main () amrex::Error("test_rhs only works for jacobian = 1"); } - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_sdc/main.cpp b/unit_test/test_sdc/main.cpp index 7967b7ecb6..0a9886b46b 100644 --- a/unit_test/test_sdc/main.cpp +++ b/unit_test/test_sdc/main.cpp @@ -102,7 +102,7 @@ void main_main () init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization diff --git a/unit_test/test_sdc_vode_rhs/main.cpp b/unit_test/test_sdc_vode_rhs/main.cpp index 1db72c42b6..8369c89429 100644 --- a/unit_test/test_sdc_vode_rhs/main.cpp +++ b/unit_test/test_sdc_vode_rhs/main.cpp @@ -41,7 +41,7 @@ void main_main () init_unit_test(); - // C++ EOS initialization (must be done after Fortran eos_init and init_extern_parameters) + // C++ EOS initialization (must be done after init_extern_parameters) eos_init(small_temp, small_dens); // C++ Network, RHS, screening, rates initialization From a5408bad445e49e4b8c0622843905a4c08b1ed6d Mon Sep 17 00:00:00 2001 From: Zhi Chen <62574124+zhichen3@users.noreply.github.com> Date: Sun, 14 Jul 2024 15:02:00 -0400 Subject: [PATCH 3/5] update self-consistent NSE SDC update with Tabular NSE (#1569) the self-consistent and tabular code paths are now consistent -- both do an update that should be second-order accurate. --- integration/nse_update_sdc.H | 363 +++++++++++------- networks/CNO_extras/actual_rhs.H | 80 +++- networks/ECSN/actual_rhs.H | 83 +++- networks/He-C-Fe-group/actual_rhs.H | 183 ++++++++- networks/ase/actual_rhs.H | 82 +++- .../C-burn-simple/actual_rhs.H | 52 ++- .../ignition_reaclib/URCA-medium/actual_rhs.H | 92 ++++- .../ignition_reaclib/URCA-simple/actual_rhs.H | 81 +++- networks/nova/actual_rhs.H | 62 ++- networks/nova2/actual_rhs.H | 70 +++- networks/partition_test/actual_rhs.H | 46 ++- networks/sn160/actual_rhs.H | 356 ++++++++++++++++- networks/subch_base/actual_rhs.H | 72 +++- networks/subch_simple/actual_rhs.H | 80 +++- nse_solver/Make.package | 1 + nse_solver/nse_eos.H | 174 +++++++++ 16 files changed, 1734 insertions(+), 143 deletions(-) create mode 100644 nse_solver/nse_eos.H diff --git a/integration/nse_update_sdc.H b/integration/nse_update_sdc.H index bfc24ec230..ed4913bb27 100644 --- a/integration/nse_update_sdc.H +++ b/integration/nse_update_sdc.H @@ -21,6 +21,7 @@ #endif #ifdef NSE_NET #include +#include #endif using namespace amrex::literals; @@ -301,207 +302,311 @@ void sdc_nse_burn(BurnT& state, const amrex::Real dt) { } #else +/// Self-Consistent NSE VERSION + /// -/// update the state due to NSE changes for the simplified-SDC case -/// this version works with the self-consistent NSE +/// this acts as an explicit Euler step for the system (rho e, rho) +/// on input, *_source are the reactive sources at time t0 and on output +/// they are the sources at time t0+dt /// -template AMREX_GPU_HOST_DEVICE AMREX_INLINE -void sdc_nse_burn(BurnT& state, const amrex::Real dt) { +void nse_derivs(const amrex::Real rho0, const amrex::Real rhoe0, + const amrex::Real rhoYe0, const amrex::Real T_old, + const amrex::Real dt, amrex::Real &mu_p, amrex::Real &mu_n, + amrex::Real &drhoyedt_weak, const amrex::Real drhoyedt_a, + amrex::Array1D &ydot_weak, + const amrex::Real *ydot_a, amrex::Real& drhoedt, + const amrex::Real T_fixed) { - state.success = true; - state.n_rhs = 0; - state.n_jac = 0; + // start with the current state and compute T - // store the initial mass fractions -- we will need these - // to compute the energy release. + amrex::Real T0; + amrex::Real Ye0 = rhoYe0 / rho0; + amrex::Real abar{0.0_rt}; - amrex::Real X_old[NumSpec]; + // Get initial temperature and abar for - for (int n = 0; n < NumSpec; ++n) { - X_old[n] = state.y[SFS+n] / state.y[SRHO]; + if (T_fixed > 0) { + T0 = T_fixed; + } else { + // initial guess for eos + T0 = T_old; + amrex::Real e0 = rhoe0 / rho0; + nse_T_abar_from_e(rho0, e0, Ye0, T0, abar, mu_p, mu_n); } - // density and momentum have no reactive sources - amrex::Real rho_old = state.y[SRHO]; + // initialize burn_state for NSE state - state.y[SRHO] += dt * state.ydot_a[SRHO]; - state.y[SMX] += dt * state.ydot_a[SMX]; - state.y[SMY] += dt * state.ydot_a[SMY]; - state.y[SMZ] += dt * state.ydot_a[SMZ]; + burn_t burn_state; - // if we are doing drive_initial_convection, we want to use - // the temperature that comes in through T_fixed + burn_state.T = T0; + burn_state.rho = rho0; + burn_state.y_e = Ye0; + burn_state.mu_p = mu_p; + burn_state.mu_n = mu_n; - amrex::Real T_in = state.T_fixed > 0.0_rt ? state.T_fixed : state.T; + // get NSE state at t0 - // get the neutrino loss term -- we want to use the state that we - // came in here with, so the original Abar and Zbar - amrex::Real snu{0.0}; - amrex::Real dsnudt{0.0}; - amrex::Real dsnudd{0.0}; - amrex::Real dsnuda{0.0}; - amrex::Real dsnudz{0.0}; + auto nse_state0 = get_actual_nse_state(burn_state, 1.0e-10_rt, true); - amrex::Real abar{0.0}; - amrex::Real zbar{0.0}; - for (int n = 0; n < NumSpec; ++n) { - abar += X_old[n] * aion_inv[n]; - zbar += X_old[n] * zion[n] * aion_inv[n]; + // Still need to get abar for T_fixed > 0 + + if (T_fixed > 0) { + for (int n = 0; n < NumSpec; ++n) { + abar += nse_state0.xn[n] * aion_inv[n]; + } + abar = 1.0_rt / abar; } - abar = 1.0 / abar; - zbar *= abar; + + // compute the plasma neutrino losses + // here abar should already be in-sync with NSE + + amrex::Real snu{0.0_rt}; + amrex::Real dsnudt{0.0_rt}; + amrex::Real dsnudd{0.0_rt}; + amrex::Real dsnuda{0.0_rt}; + amrex::Real dsnudz{0.0_rt}; + + amrex::Real zbar = abar * Ye0; #ifdef NEUTRINOS { constexpr int do_derivatives = 0; - sneut5(T_in, rho_old, abar, zbar, - snu, dsnudt, dsnudd, dsnuda, dsnudz); + sneut5(T0, rho0, abar, zbar, + snu, dsnudt, dsnudd, + dsnuda, dsnudz); } #endif - amrex::Real snu_old = snu; + // Get molar fractions, ydots and neutrino loss due to weak rates - // if our network could return the evolution of Ye due to the - // weak interactions, we would evaluate the NSE state here and - // get dYe/dt. + amrex::Array1D Y; + amrex::Real e_nu {0.0_rt}; - //amrex::Real dyedt_old = 0.0; + for (int n = 1; n <= NumSpec; ++n) { + Y(n) = nse_state0.xn[n-1] * aion_inv[n-1]; + } + // Fill in ydot with only weak rates contributing - // predict the U^{n+1,*} state with only estimates of the X - // updates with advection to dt and the neutrino loss term in - // energy + get_ydot_weak(nse_state0, ydot_weak, e_nu, Y); - BurnT burn_state; - burn_state.T = T_in; // initial guess - burn_state.mu_p = state.mu_p; - burn_state.mu_n = state.mu_n; + /// + /// construct initial sources at t0 + /// - amrex::Real rhoX_source[NumSpec] = {0.0_rt}; + // Find drhoyedt_weak from weak reaction - amrex::Real rhoX_new[NumSpec]; + drhoyedt_weak = 0.0_rt; + for (int n = 0; n < NumSpec; ++n) { + drhoyedt_weak += rho0 * zion[n] * ydot_weak(n+1); + } - amrex::Real rhoe_new; - amrex::Real rho_enucdot = -rho_old * snu; + // find rhoe_source - amrex::Real rho_half = 0.5_rt * (rho_old + state.y[SRHO]); + amrex::Real rhoe_source{0.0_rt}; + for (int n = 1; n <= NumSpec; ++n) { + rhoe_source += rho0 * ydot_weak(n) * network::mion(n); + } + rhoe_source *= C::Legacy::enuc_conv2; - // predict the new aux for the first iteration -- this is really - // just including the advection bits + // include plasma neutrino terms - for (int n = 0; n < NumSpec; n++) { - rhoX_new[n] = state.y[SFS+n] + dt * state.ydot_a[SFS+n] + dt * rhoX_source[n]; + rhoe_source -= rho0 * snu; + + // include neutrino loss terms from weak rates + // Note that e_nu here is already negative so we add here. + + if (integrator_rp::nse_include_enu_weak == 1) { + rhoe_source += rho0 * e_nu; } - burn_t nse_state; + // + // evolve for eps * dt; + // - for (int iter = 0; iter < integrator_rp::nse_iters; iter++) { + amrex::Real tau = integrator_rp::nse_deriv_dt_factor * dt; - // update (rho e)^{n+1} based on the new energy generation rate - rhoe_new = state.y[SEINT] + dt * state.ydot_a[SEINT] + dt * rho_enucdot; + amrex::Real rho1 = rho0 + tau * ydot_a[SRHO]; + amrex::Real rhoe1 = rhoe0 + tau * (ydot_a[SEINT] + rhoe_source); + amrex::Real rhoYe1 = rhoYe0 + tau * (drhoyedt_weak + drhoyedt_a); - // call the EOS to get the updated T* + // compute the temperature at t0 + tau - amrex::Real T_new; - burn_state.rho = state.y[SRHO]; - burn_state.e = rhoe_new / state.y[SRHO]; - for (int n = 0; n < NumSpec; n++) { - burn_state.xn[n] = rhoX_new[n] / state.y[SRHO]; - burn_state.y[SFS+n] = rhoX_new[n]; - } + amrex::Real T1; + amrex::Real Ye1 = rhoYe1 / rho1; - if (state.T_fixed > 0) { - T_new = state.T_fixed; - } else { - eos(eos_input_re, burn_state); - T_new = burn_state.T; - } - burn_state.T = T_new; + if (T_fixed > 0) { + T1 = T_fixed; + } else { + // initial guess for eos + T1 = T0; + amrex::Real e1 = rhoe1 / rho1; + amrex::Real abar1; + nse_T_abar_from_e(rho1, e1, Ye1, T1, abar1, mu_p, mu_n); + } - // solve for the NSE state for this network we need to call - // get_nse_state with a burn_t. We will have it compute - // Ye from the input X's + // update burn_state at t0 + tau - nse_state = get_actual_nse_state(burn_state); + burn_state.T = T1; + burn_state.rho = rho1; + burn_state.y_e = Ye1; - // compute the energy release. The mass fractions in nse_state.xn[] - // include the advective parts, so first we need to remove that. + // call NSE at t0 + tau with - amrex::Real rhoX_tilde[NumSpec]; - for (int n = 0; n < NumSpec; ++n) { - rhoX_tilde[n] = state.y[SRHO] * nse_state.xn[n] - dt * state.ydot_a[SFS+n]; - } + // Note this new NSE state has advection contribution - //amrex::Real dyedt = 0.0_rt; // we can update this in the future by calling actual_rhs() + auto nse_state1 = get_actual_nse_state(burn_state, 1.0e-10_rt, true); - // we want to compute (rho eps) = - N_A c^2 sum{m_i (rhoX_tilde - rhoX_old) / A_i} - rho_enucdot = 0.0; - for (int n = 0; n < NumSpec; ++n) { - rho_enucdot += (rhoX_tilde[n] - rho_old * X_old[n]) * - network::mion(n+1) * aion_inv[n]; - } - rho_enucdot *= C::Legacy::enuc_conv2; + // construct the finite-difference approximation to the derivatives + // subtract off advection contribution to get pure drhoedt from reaction - // now get the updated neutrino term - abar = 0.0; - zbar = 0.0; - for (int n = 0; n < NumSpec; ++n) { - abar += nse_state.xn[n] * aion_inv[n]; - zbar += nse_state.xn[n] * zion[n] * aion_inv[n]; - } - abar = 1.0 / abar; - zbar *= abar; + drhoedt = 0.0_rt; + for (int n = 0; n < NumSpec; ++n) { + drhoedt += (nse_state1.y[SFS+n] - tau * ydot_a[SFS+n] - + nse_state0.y[SFS+n]) * aion_inv[n] * network::mion(n+1); + } + drhoedt *= C::Legacy::enuc_conv2 / tau; -#ifdef NEUTRINOS - constexpr int do_derivatives = 0; - sneut5(T_new, state.y[SRHO], abar, zbar, - snu, dsnudt, dsnudd, dsnuda, dsnudz); -#endif + // include neutrino terms again + + drhoedt -= rho1 * snu; + if (integrator_rp::nse_include_enu_weak == 1) { + drhoedt += rho1 * e_nu; + } - rho_enucdot -= 0.5_rt * rho_half * (snu_old + snu); +} - // update the new state for the next pass -- this should - // already implicitly have the advective portion included, - // since it was there when we called the NSE state - for (int n = 0; n < NumSpec; n++) { - rhoX_new[n] = state.y[SRHO] * nse_state.xn[n]; - } +/// +/// update the state due to NSE changes for the simplified-SDC case +/// this version works with the self-consistent NSE +/// + +template +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void sdc_nse_burn(BurnT& state, const amrex::Real dt) { + + state.success = true; + state.n_rhs = 0; + state.n_jac = 0; + + // store the initial state + + amrex::Real rho_old = state.y[SRHO]; + amrex::Real rhoe_old = state.y[SEINT]; + amrex::Real rhoye_old = rho_old * state.y_e; + amrex::Real T_old = state.T; + amrex::Real mu_p = state.mu_p; + amrex::Real mu_n = state.mu_n; + // density and momentum have no reactive sources + + state.y[SRHO] += dt * state.ydot_a[SRHO]; + state.y[SMX] += dt * state.ydot_a[SMX]; + state.y[SMY] += dt * state.ydot_a[SMY]; + state.y[SMZ] += dt * state.ydot_a[SMZ]; + + // do an RK2 integration + + // get the derivatives, rho_dyedt, drhoedt, ydot_weak at t = t^n + + amrex::Array1D ydot_weak; + amrex::Real drhoedt{0.0_rt}; + amrex::Real drhoyedt_weak{0.0_rt}; + amrex::Real drhoyedt_a{0.0_rt}; + + // find the advection contribution: drhoyedt_a and drhoabardt_a + + for (int n = 0; n < NumSpec; ++n) { + drhoyedt_a += state.ydot_a[SFS+n] * zion[n] * aion_inv[n]; } - // now update the aux quantities + nse_derivs(rho_old, rhoe_old, + rhoye_old, T_old, + dt, mu_p, mu_n, + drhoyedt_weak, drhoyedt_a, + ydot_weak, state.ydot_a, + drhoedt, state.T_fixed); - // the new mass fractions are just those that come from the table + // evolve to the midpoint in time + + amrex::Real rho_tmp = rho_old + 0.5_rt * dt * state.ydot_a[SRHO]; + amrex::Real rhoe_tmp = rhoe_old + 0.5_rt * dt * (state.ydot_a[SEINT] + drhoedt); + amrex::Real rhoye_tmp = rhoye_old + 0.5_rt * dt * (drhoyedt_weak + drhoyedt_a); + + // compute the derivatives at the midpoint in time + + nse_derivs(rho_tmp, rhoe_tmp, + rhoye_tmp, T_old, + dt, mu_p, mu_n, + drhoyedt_weak, drhoyedt_a, + ydot_weak, state.ydot_a, + drhoedt, state.T_fixed); + + // evolve to the new time + + amrex::Real rho_new = rho_old + dt * state.ydot_a[SRHO]; + amrex::Real rhoe_new = rhoe_old + dt * (state.ydot_a[SEINT] + drhoedt); + amrex::Real rhoye_new = rhoye_old + dt * (drhoyedt_weak + drhoyedt_a); + + // Update the temperature with new internal energy + + amrex::Real T_new; + amrex::Real Ye_new = rhoye_new / rho_new; + + if (state.T_fixed > 0) { + T_new = state.T_fixed; + } else { + // initial eos guess + T_new = T_old; + amrex::Real e_new = rhoe_new / rho_new; + amrex::Real abar_new; + nse_T_abar_from_e(rho_new, e_new, Ye_new, T_new, abar_new, mu_p, mu_n); + } + + // With updated temp, rho, and Ye get final NSE state + + burn_t burn_state; + burn_state.T = T_new; + burn_state.rho = rho_new; + burn_state.y_e = Ye_new; + burn_state.mu_p = mu_p; + burn_state.mu_n = mu_n; + + auto nse_state = get_actual_nse_state(burn_state, 1.0e-10_rt, true); + + // store the new state // make sure they are normalized + amrex::Real sum_X{0.0_rt}; - amrex::Real X[NumSpec] = {0.0_rt}; - for (int n = 0; n < NumSpec; ++n) { - X[n] = amrex::max(small_x, amrex::min(1.0_rt, nse_state.xn[n])); - sum_X += X[n]; + for (auto & xn : nse_state.xn) { + xn = std::clamp(xn, small_x, 1.0_rt); + sum_X += xn; } - for (int n = 0; n < NumSpec; ++n) { - X[n] /= sum_X; + for (auto & xn : nse_state.xn) { + xn /= sum_X; } for (int n = 0; n < NumSpec; ++n) { - state.y[SFS+n] = state.y[SRHO] * X[n]; + state.y[SFS+n] = rho_new * nse_state.xn[n]; } - // density and momenta have already been updated - // update the total and internal energy now + // get the energy release from the change in rhoe over the timestep + // excluding any advection, and use that to update the total energy - state.y[SEINT] += dt * state.ydot_a[SEINT] + dt * rho_enucdot; + amrex::Real rho_enucdot = (rhoe_new - rhoe_old) / dt - state.ydot_a[SEINT]; + + state.y[SEINT] = rhoe_new; state.y[SEDEN] += dt * state.ydot_a[SEDEN] + dt * rho_enucdot; // store the chemical potentials + state.mu_p = nse_state.mu_p; state.mu_n = nse_state.mu_n; - } #endif diff --git a/networks/CNO_extras/actual_rhs.H b/networks/CNO_extras/actual_rhs.H index a9eaa3f3c0..b74804639b 100644 --- a/networks/CNO_extras/actual_rhs.H +++ b/networks/CNO_extras/actual_rhs.H @@ -693,11 +693,87 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(C13) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(N15) = 0.0_rt; + + ydot_nuc(O14) = 0.0_rt; + + ydot_nuc(O15) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(O17) = 0.0_rt; + + ydot_nuc(O18) = 0.0_rt; + + ydot_nuc(F17) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + ydot_nuc(F19) = 0.0_rt; + + ydot_nuc(Ne18) = 0.0_rt; + + ydot_nuc(Ne19) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Mg22) = 0.0_rt; + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Fe56) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, @@ -858,7 +934,7 @@ void rhs_nuc(const burn_t& state, (screened_rates(k_He4_Ne20_to_Mg24)*Y(He4)*Y(Ne20)*state.rho + -screened_rates(k_Mg24_to_He4_Ne20)*Y(Mg24)) + (screened_rates(k_C12_O16_to_He4_Mg24)*Y(C12)*Y(O16)*state.rho + -screened_rates(k_He4_Mg24_to_C12_O16)*Y(He4)*Y(Mg24)*state.rho); - ydot_nuc(Fe56) = 0.0; + ydot_nuc(Fe56) = 0.0_rt; } diff --git a/networks/ECSN/actual_rhs.H b/networks/ECSN/actual_rhs.H index 251dd5610e..0b404497ef 100644 --- a/networks/ECSN/actual_rhs.H +++ b/networks/ECSN/actual_rhs.H @@ -275,7 +275,7 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; tabular_evaluate(j_F20_O20_meta, j_F20_O20_rhoy, j_F20_O20_temp, j_F20_O20_data, rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); @@ -312,6 +312,87 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + tabular_evaluate(j_F20_O20_meta, j_F20_O20_rhoy, j_F20_O20_temp, j_F20_O20_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_F20_to_O20) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(F20) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Ne20_F20_meta, j_Ne20_F20_rhoy, j_Ne20_F20_temp, j_Ne20_F20_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Ne20_to_F20) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Ne20) * (edot_nu + edot_gamma); + + tabular_evaluate(j_O20_F20_meta, j_O20_F20_rhoy, j_O20_F20_temp, j_O20_F20_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_O20_to_F20) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(O20) * (edot_nu + edot_gamma); + + tabular_evaluate(j_F20_Ne20_meta, j_F20_Ne20_rhoy, j_F20_Ne20_temp, j_F20_Ne20_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_F20_to_Ne20) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(F20) * (edot_nu + edot_gamma); + + auto screened_rates = rate_eval.screened_rates; + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(O20) = + (-screened_rates(k_O20_to_F20)*Y(O20) + screened_rates(k_F20_to_O20)*Y(F20)); + + ydot_nuc(F20) = + (screened_rates(k_O20_to_F20)*Y(O20) + -screened_rates(k_F20_to_O20)*Y(F20)) + + (-screened_rates(k_F20_to_Ne20)*Y(F20) + screened_rates(k_Ne20_to_F20)*Y(Ne20)); + + ydot_nuc(Ne20) = + (screened_rates(k_F20_to_Ne20)*Y(F20) + -screened_rates(k_Ne20_to_F20)*Y(Ne20)); + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Al27) = 0.0_rt; + + ydot_nuc(Si28) = 0.0_rt; + + ydot_nuc(P31) = 0.0_rt; + + ydot_nuc(S32) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/He-C-Fe-group/actual_rhs.H b/networks/He-C-Fe-group/actual_rhs.H index 12c101e7ca..8736460ae3 100644 --- a/networks/He-C-Fe-group/actual_rhs.H +++ b/networks/He-C-Fe-group/actual_rhs.H @@ -1276,7 +1276,7 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; tabular_evaluate(j_Co55_Fe55_meta, j_Co55_Fe55_rhoy, j_Co55_Fe55_temp, j_Co55_Fe55_data, rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); @@ -1377,6 +1377,187 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + tabular_evaluate(j_Co55_Fe55_meta, j_Co55_Fe55_rhoy, j_Co55_Fe55_temp, j_Co55_Fe55_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Co55_to_Fe55) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Co55) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Co56_Fe56_meta, j_Co56_Fe56_rhoy, j_Co56_Fe56_temp, j_Co56_Fe56_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Co56_to_Fe56) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Co56) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Co56_Ni56_meta, j_Co56_Ni56_rhoy, j_Co56_Ni56_temp, j_Co56_Ni56_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Co56_to_Ni56) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Co56) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Co57_Ni57_meta, j_Co57_Ni57_rhoy, j_Co57_Ni57_temp, j_Co57_Ni57_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Co57_to_Ni57) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Co57) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Fe55_Co55_meta, j_Fe55_Co55_rhoy, j_Fe55_Co55_temp, j_Fe55_Co55_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Fe55_to_Co55) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Fe55) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Fe55_Mn55_meta, j_Fe55_Mn55_rhoy, j_Fe55_Mn55_temp, j_Fe55_Mn55_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Fe55_to_Mn55) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Fe55) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Fe56_Co56_meta, j_Fe56_Co56_rhoy, j_Fe56_Co56_temp, j_Fe56_Co56_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Fe56_to_Co56) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Fe56) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Mn55_Fe55_meta, j_Mn55_Fe55_rhoy, j_Mn55_Fe55_temp, j_Mn55_Fe55_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Mn55_to_Fe55) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Mn55) * (edot_nu + edot_gamma); + + tabular_evaluate(j_n_p_meta, j_n_p_rhoy, j_n_p_temp, j_n_p_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_n_to_p) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(N) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Ni56_Co56_meta, j_Ni56_Co56_rhoy, j_Ni56_Co56_temp, j_Ni56_Co56_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Ni56_to_Co56) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Ni56) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Ni57_Co57_meta, j_Ni57_Co57_rhoy, j_Ni57_Co57_temp, j_Ni57_Co57_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Ni57_to_Co57) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Ni57) * (edot_nu + edot_gamma); + + tabular_evaluate(j_p_n_meta, j_p_n_rhoy, j_p_n_temp, j_p_n_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_p_to_n) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(H1) * (edot_nu + edot_gamma); + + auto screened_rates = rate_eval.screened_rates; + + ydot_nuc(N) = + -screened_rates(k_n_to_p)*Y(N) + + screened_rates(k_p_to_n)*Y(H1); + + ydot_nuc(H1) = + screened_rates(k_n_to_p)*Y(N) + + -screened_rates(k_p_to_n)*Y(H1); + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Ne21) = 0.0_rt; + + ydot_nuc(Na22) = 0.0_rt; + + ydot_nuc(Na23) = 0.0_rt; + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Al27) = 0.0_rt; + + ydot_nuc(Si28) = 0.0_rt; + + ydot_nuc(P31) = 0.0_rt; + + ydot_nuc(S32) = 0.0_rt; + + ydot_nuc(Ar36) = 0.0_rt; + + ydot_nuc(Ca40) = 0.0_rt; + + ydot_nuc(Ti44) = 0.0_rt; + + ydot_nuc(Cr48) = 0.0_rt; + + ydot_nuc(Mn51) = 0.0_rt; + + ydot_nuc(Mn55) = + (screened_rates(k_Fe55_to_Mn55)*Y(Fe55) + -screened_rates(k_Mn55_to_Fe55)*Y(Mn55)); + + ydot_nuc(Fe52) = 0.0_rt; + + ydot_nuc(Fe53) = 0.0_rt; + + ydot_nuc(Fe54) = 0.0_rt; + + ydot_nuc(Fe55) = + (screened_rates(k_Co55_to_Fe55)*Y(Co55) + -screened_rates(k_Fe55_to_Co55)*Y(Fe55)) + + (-screened_rates(k_Fe55_to_Mn55)*Y(Fe55) + screened_rates(k_Mn55_to_Fe55)*Y(Mn55)); + + ydot_nuc(Fe56) = + (screened_rates(k_Co56_to_Fe56)*Y(Co56) + -screened_rates(k_Fe56_to_Co56)*Y(Fe56)); + + ydot_nuc(Co55) = + (-screened_rates(k_Co55_to_Fe55)*Y(Co55) + screened_rates(k_Fe55_to_Co55)*Y(Fe55)); + + ydot_nuc(Co56) = + (-screened_rates(k_Co56_to_Fe56)*Y(Co56) + screened_rates(k_Fe56_to_Co56)*Y(Fe56)) + + (screened_rates(k_Ni56_to_Co56)*Y(Ni56) + -screened_rates(k_Co56_to_Ni56)*Y(Co56)); + + ydot_nuc(Co57) = + (screened_rates(k_Ni57_to_Co57)*Y(Ni57) + -screened_rates(k_Co57_to_Ni57)*Y(Co57)); + + ydot_nuc(Ni56) = + (-screened_rates(k_Ni56_to_Co56)*Y(Ni56) + screened_rates(k_Co56_to_Ni56)*Y(Co56)); + + ydot_nuc(Ni57) = + (-screened_rates(k_Ni57_to_Co57)*Y(Ni57) + screened_rates(k_Co57_to_Ni57)*Y(Co57)); + + ydot_nuc(Ni58) = 0.0_rt; + + ydot_nuc(Cu59) = 0.0_rt; + + ydot_nuc(Zn60) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/ase/actual_rhs.H b/networks/ase/actual_rhs.H index f5925ec203..4d08b90d65 100644 --- a/networks/ase/actual_rhs.H +++ b/networks/ase/actual_rhs.H @@ -894,11 +894,91 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(N) = 0.0_rt; + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Ne21) = 0.0_rt; + + ydot_nuc(Na22) = 0.0_rt; + + ydot_nuc(Na23) = 0.0_rt; + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Al27) = 0.0_rt; + + ydot_nuc(Si28) = 0.0_rt; + + ydot_nuc(P31) = 0.0_rt; + + ydot_nuc(S32) = 0.0_rt; + + ydot_nuc(Ar36) = 0.0_rt; + + ydot_nuc(Ca40) = 0.0_rt; + + ydot_nuc(Ti44) = 0.0_rt; + + ydot_nuc(Cr48) = 0.0_rt; + + ydot_nuc(Fe52) = 0.0_rt; + + ydot_nuc(Ni56) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/ignition_reaclib/C-burn-simple/actual_rhs.H b/networks/ignition_reaclib/C-burn-simple/actual_rhs.H index ad623a6c06..d96207d283 100644 --- a/networks/ignition_reaclib/C-burn-simple/actual_rhs.H +++ b/networks/ignition_reaclib/C-burn-simple/actual_rhs.H @@ -140,11 +140,61 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(N) = 0.0_rt; + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Na23) = 0.0_rt; + + ydot_nuc(Mg23) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/ignition_reaclib/URCA-medium/actual_rhs.H b/networks/ignition_reaclib/URCA-medium/actual_rhs.H index 1b9f728b47..a6d35692e1 100644 --- a/networks/ignition_reaclib/URCA-medium/actual_rhs.H +++ b/networks/ignition_reaclib/URCA-medium/actual_rhs.H @@ -300,7 +300,7 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; tabular_evaluate(j_Na23_Ne23_meta, j_Na23_Ne23_rhoy, j_Na23_Ne23_temp, j_Na23_Ne23_data, rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); @@ -345,6 +345,96 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + tabular_evaluate(j_Na23_Ne23_meta, j_Na23_Ne23_rhoy, j_Na23_Ne23_temp, j_Na23_Ne23_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Na23_to_Ne23) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Na23) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Ne23_Na23_meta, j_Ne23_Na23_rhoy, j_Ne23_Na23_temp, j_Ne23_Na23_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Ne23_to_Na23) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Ne23) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Mg23_Na23_meta, j_Mg23_Na23_rhoy, j_Mg23_Na23_temp, j_Mg23_Na23_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Mg23_to_Na23) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Mg23) * (edot_nu + edot_gamma); + + tabular_evaluate(j_n_p_meta, j_n_p_rhoy, j_n_p_temp, j_n_p_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_n_to_p) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(N) * (edot_nu + edot_gamma); + + tabular_evaluate(j_p_n_meta, j_p_n_rhoy, j_p_n_temp, j_p_n_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_p_to_n) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(H1) * (edot_nu + edot_gamma); + + auto screened_rates = rate_eval.screened_rates; + + ydot_nuc(N) = + -screened_rates(k_n_to_p)*Y(N) + + screened_rates(k_p_to_n)*Y(H1); + + ydot_nuc(H1) = + screened_rates(k_n_to_p)*Y(N) + + -screened_rates(k_p_to_n)*Y(H1); + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Ne23) = + (-screened_rates(k_Ne23_to_Na23)*Y(Ne23) + screened_rates(k_Na23_to_Ne23)*Y(Na23)); + + ydot_nuc(Na23) = + (screened_rates(k_Ne23_to_Na23)*Y(Ne23) + -screened_rates(k_Na23_to_Ne23)*Y(Na23)) + + screened_rates(k_Mg23_to_Na23)*Y(Mg23); + + ydot_nuc(Mg23) = + -screened_rates(k_Mg23_to_Na23)*Y(Mg23); + + ydot_nuc(Mg24) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/ignition_reaclib/URCA-simple/actual_rhs.H b/networks/ignition_reaclib/URCA-simple/actual_rhs.H index 627b2a97f8..ff15e513fb 100644 --- a/networks/ignition_reaclib/URCA-simple/actual_rhs.H +++ b/networks/ignition_reaclib/URCA-simple/actual_rhs.H @@ -140,7 +140,7 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; tabular_evaluate(j_Na23_Ne23_meta, j_Na23_Ne23_rhoy, j_Na23_Ne23_temp, j_Na23_Ne23_data, rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); @@ -177,6 +177,85 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + tabular_evaluate(j_Na23_Ne23_meta, j_Na23_Ne23_rhoy, j_Na23_Ne23_temp, j_Na23_Ne23_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Na23_to_Ne23) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Na23) * (edot_nu + edot_gamma); + + tabular_evaluate(j_Ne23_Na23_meta, j_Ne23_Na23_rhoy, j_Ne23_Na23_temp, j_Ne23_Na23_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_Ne23_to_Na23) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(Ne23) * (edot_nu + edot_gamma); + + tabular_evaluate(j_n_p_meta, j_n_p_rhoy, j_n_p_temp, j_n_p_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_n_to_p) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(N) * (edot_nu + edot_gamma); + + tabular_evaluate(j_p_n_meta, j_p_n_rhoy, j_p_n_temp, j_p_n_data, + rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma); + rate_eval.screened_rates(k_p_to_n) = rate; + rate_eval.enuc_weak += C::Legacy::n_A * Y(H1) * (edot_nu + edot_gamma); + + auto screened_rates = rate_eval.screened_rates; + + ydot_nuc(N) = + -screened_rates(k_n_to_p)*Y(N) + + screened_rates(k_p_to_n)*Y(H1); + + ydot_nuc(H1) = + screened_rates(k_n_to_p)*Y(N) + + -screened_rates(k_p_to_n)*Y(H1); + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Ne23) = + (-screened_rates(k_Ne23_to_Na23)*Y(Ne23) + screened_rates(k_Na23_to_Ne23)*Y(Na23)); + + ydot_nuc(Na23) = + (screened_rates(k_Ne23_to_Na23)*Y(Ne23) + -screened_rates(k_Na23_to_Ne23)*Y(Na23)); + + ydot_nuc(Mg23) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/nova/actual_rhs.H b/networks/nova/actual_rhs.H index 0c2dd612d3..8f3e7c1c60 100644 --- a/networks/nova/actual_rhs.H +++ b/networks/nova/actual_rhs.H @@ -362,11 +362,71 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(C13) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(N15) = 0.0_rt; + + ydot_nuc(O14) = 0.0_rt; + + ydot_nuc(O15) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(O17) = 0.0_rt; + + ydot_nuc(F17) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/nova2/actual_rhs.H b/networks/nova2/actual_rhs.H index 1219595277..6cb140fd89 100644 --- a/networks/nova2/actual_rhs.H +++ b/networks/nova2/actual_rhs.H @@ -559,11 +559,79 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(H2) = 0.0_rt; + + ydot_nuc(He3) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(Be7) = 0.0_rt; + + ydot_nuc(B8) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(C13) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(N15) = 0.0_rt; + + ydot_nuc(O14) = 0.0_rt; + + ydot_nuc(O15) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(O17) = 0.0_rt; + + ydot_nuc(F17) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/partition_test/actual_rhs.H b/networks/partition_test/actual_rhs.H index ba170ca871..72296fab0e 100644 --- a/networks/partition_test/actual_rhs.H +++ b/networks/partition_test/actual_rhs.H @@ -140,11 +140,55 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(Fe52) = 0.0_rt; + + ydot_nuc(Co55) = 0.0_rt; + + ydot_nuc(Ni56) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/sn160/actual_rhs.H b/networks/sn160/actual_rhs.H index 9d14330199..3ff15b8115 100644 --- a/networks/sn160/actual_rhs.H +++ b/networks/sn160/actual_rhs.H @@ -8815,11 +8815,365 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(N) = 0.0_rt; + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(H2) = 0.0_rt; + + ydot_nuc(He3) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(Li6) = 0.0_rt; + + ydot_nuc(Li7) = 0.0_rt; + + ydot_nuc(Be7) = 0.0_rt; + + ydot_nuc(Be9) = 0.0_rt; + + ydot_nuc(B8) = 0.0_rt; + + ydot_nuc(B10) = 0.0_rt; + + ydot_nuc(B11) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(C13) = 0.0_rt; + + ydot_nuc(C14) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(N15) = 0.0_rt; + + ydot_nuc(O14) = 0.0_rt; + + ydot_nuc(O15) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(O17) = 0.0_rt; + + ydot_nuc(O18) = 0.0_rt; + + ydot_nuc(F17) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + ydot_nuc(F19) = 0.0_rt; + + ydot_nuc(Ne18) = 0.0_rt; + + ydot_nuc(Ne19) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Ne21) = 0.0_rt; + + ydot_nuc(Ne22) = 0.0_rt; + + ydot_nuc(Na21) = 0.0_rt; + + ydot_nuc(Na22) = 0.0_rt; + + ydot_nuc(Na23) = 0.0_rt; + + ydot_nuc(Mg23) = 0.0_rt; + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Mg25) = 0.0_rt; + + ydot_nuc(Mg26) = 0.0_rt; + + ydot_nuc(Al25) = 0.0_rt; + + ydot_nuc(Al26) = 0.0_rt; + + ydot_nuc(Al27) = 0.0_rt; + + ydot_nuc(Si28) = 0.0_rt; + + ydot_nuc(Si29) = 0.0_rt; + + ydot_nuc(Si30) = 0.0_rt; + + ydot_nuc(Si31) = 0.0_rt; + + ydot_nuc(Si32) = 0.0_rt; + + ydot_nuc(P29) = 0.0_rt; + + ydot_nuc(P30) = 0.0_rt; + + ydot_nuc(P31) = 0.0_rt; + + ydot_nuc(P32) = 0.0_rt; + + ydot_nuc(P33) = 0.0_rt; + + ydot_nuc(S32) = 0.0_rt; + + ydot_nuc(S33) = 0.0_rt; + + ydot_nuc(S34) = 0.0_rt; + + ydot_nuc(S35) = 0.0_rt; + + ydot_nuc(S36) = 0.0_rt; + + ydot_nuc(Cl33) = 0.0_rt; + + ydot_nuc(Cl34) = 0.0_rt; + + ydot_nuc(Cl35) = 0.0_rt; + + ydot_nuc(Cl36) = 0.0_rt; + + ydot_nuc(Cl37) = 0.0_rt; + + ydot_nuc(Ar36) = 0.0_rt; + + ydot_nuc(Ar37) = 0.0_rt; + + ydot_nuc(Ar38) = 0.0_rt; + + ydot_nuc(Ar39) = 0.0_rt; + + ydot_nuc(Ar40) = 0.0_rt; + + ydot_nuc(K37) = 0.0_rt; + + ydot_nuc(K38) = 0.0_rt; + + ydot_nuc(K39) = 0.0_rt; + + ydot_nuc(K40) = 0.0_rt; + + ydot_nuc(K41) = 0.0_rt; + + ydot_nuc(Ca40) = 0.0_rt; + + ydot_nuc(Ca41) = 0.0_rt; + + ydot_nuc(Ca42) = 0.0_rt; + + ydot_nuc(Ca43) = 0.0_rt; + + ydot_nuc(Ca44) = 0.0_rt; + + ydot_nuc(Ca45) = 0.0_rt; + + ydot_nuc(Ca46) = 0.0_rt; + + ydot_nuc(Ca47) = 0.0_rt; + + ydot_nuc(Ca48) = 0.0_rt; + + ydot_nuc(Sc43) = 0.0_rt; + + ydot_nuc(Sc44) = 0.0_rt; + + ydot_nuc(Sc45) = 0.0_rt; + + ydot_nuc(Sc46) = 0.0_rt; + + ydot_nuc(Sc47) = 0.0_rt; + + ydot_nuc(Sc48) = 0.0_rt; + + ydot_nuc(Sc49) = 0.0_rt; + + ydot_nuc(Ti44) = 0.0_rt; + + ydot_nuc(Ti45) = 0.0_rt; + + ydot_nuc(Ti46) = 0.0_rt; + + ydot_nuc(Ti47) = 0.0_rt; + + ydot_nuc(Ti48) = 0.0_rt; + + ydot_nuc(Ti49) = 0.0_rt; + + ydot_nuc(Ti50) = 0.0_rt; + + ydot_nuc(Ti51) = 0.0_rt; + + ydot_nuc(V46) = 0.0_rt; + + ydot_nuc(V47) = 0.0_rt; + + ydot_nuc(V48) = 0.0_rt; + + ydot_nuc(V49) = 0.0_rt; + + ydot_nuc(V50) = 0.0_rt; + + ydot_nuc(V51) = 0.0_rt; + + ydot_nuc(V52) = 0.0_rt; + + ydot_nuc(Cr48) = 0.0_rt; + + ydot_nuc(Cr49) = 0.0_rt; + + ydot_nuc(Cr50) = 0.0_rt; + + ydot_nuc(Cr51) = 0.0_rt; + + ydot_nuc(Cr52) = 0.0_rt; + + ydot_nuc(Cr53) = 0.0_rt; + + ydot_nuc(Cr54) = 0.0_rt; + + ydot_nuc(Mn50) = 0.0_rt; + + ydot_nuc(Mn51) = 0.0_rt; + + ydot_nuc(Mn52) = 0.0_rt; + + ydot_nuc(Mn53) = 0.0_rt; + + ydot_nuc(Mn54) = 0.0_rt; + + ydot_nuc(Mn55) = 0.0_rt; + + ydot_nuc(Fe52) = 0.0_rt; + + ydot_nuc(Fe53) = 0.0_rt; + + ydot_nuc(Fe54) = 0.0_rt; + + ydot_nuc(Fe55) = 0.0_rt; + + ydot_nuc(Fe56) = 0.0_rt; + + ydot_nuc(Fe57) = 0.0_rt; + + ydot_nuc(Fe58) = 0.0_rt; + + ydot_nuc(Co53) = 0.0_rt; + + ydot_nuc(Co54) = 0.0_rt; + + ydot_nuc(Co55) = 0.0_rt; + + ydot_nuc(Co56) = 0.0_rt; + + ydot_nuc(Co57) = 0.0_rt; + + ydot_nuc(Co58) = 0.0_rt; + + ydot_nuc(Co59) = 0.0_rt; + + ydot_nuc(Ni56) = 0.0_rt; + + ydot_nuc(Ni57) = 0.0_rt; + + ydot_nuc(Ni58) = 0.0_rt; + + ydot_nuc(Ni59) = 0.0_rt; + + ydot_nuc(Ni60) = 0.0_rt; + + ydot_nuc(Ni61) = 0.0_rt; + + ydot_nuc(Ni62) = 0.0_rt; + + ydot_nuc(Ni63) = 0.0_rt; + + ydot_nuc(Ni64) = 0.0_rt; + + ydot_nuc(Cu57) = 0.0_rt; + + ydot_nuc(Cu58) = 0.0_rt; + + ydot_nuc(Cu59) = 0.0_rt; + + ydot_nuc(Cu60) = 0.0_rt; + + ydot_nuc(Cu61) = 0.0_rt; + + ydot_nuc(Cu62) = 0.0_rt; + + ydot_nuc(Cu63) = 0.0_rt; + + ydot_nuc(Cu64) = 0.0_rt; + + ydot_nuc(Cu65) = 0.0_rt; + + ydot_nuc(Zn59) = 0.0_rt; + + ydot_nuc(Zn60) = 0.0_rt; + + ydot_nuc(Zn61) = 0.0_rt; + + ydot_nuc(Zn62) = 0.0_rt; + + ydot_nuc(Zn63) = 0.0_rt; + + ydot_nuc(Zn64) = 0.0_rt; + + ydot_nuc(Zn65) = 0.0_rt; + + ydot_nuc(Zn66) = 0.0_rt; + + ydot_nuc(Ga62) = 0.0_rt; + + ydot_nuc(Ga63) = 0.0_rt; + + ydot_nuc(Ga64) = 0.0_rt; + + ydot_nuc(Ge63) = 0.0_rt; + + ydot_nuc(Ge64) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/subch_base/actual_rhs.H b/networks/subch_base/actual_rhs.H index 6b083c0284..c1fae4ff7e 100644 --- a/networks/subch_base/actual_rhs.H +++ b/networks/subch_base/actual_rhs.H @@ -804,11 +804,81 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Na23) = 0.0_rt; + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Al27) = 0.0_rt; + + ydot_nuc(Si28) = 0.0_rt; + + ydot_nuc(P31) = 0.0_rt; + + ydot_nuc(S32) = 0.0_rt; + + ydot_nuc(Ar36) = 0.0_rt; + + ydot_nuc(Ca40) = 0.0_rt; + + ydot_nuc(Ti44) = 0.0_rt; + + ydot_nuc(Cr48) = 0.0_rt; + + ydot_nuc(Fe52) = 0.0_rt; + + ydot_nuc(Ni56) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/networks/subch_simple/actual_rhs.H b/networks/subch_simple/actual_rhs.H index 5a72edee3c..9e86ce9ddf 100644 --- a/networks/subch_simple/actual_rhs.H +++ b/networks/subch_simple/actual_rhs.H @@ -1038,11 +1038,89 @@ void evaluate_rates(const burn_t& state, T& rate_eval) { [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; - rate_eval.enuc_weak = 0.0; + rate_eval.enuc_weak = 0.0_rt; } +#ifdef NSE_NET +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void get_ydot_weak(const burn_t& state, + amrex::Array1D& ydot_nuc, + amrex::Real& enuc_weak, + [[maybe_unused]] const amrex::Array1D& Y) { + /// + /// Calculate Ydots contribute only from weak reactions. + /// This is used to calculate dyedt and energy generation from + /// weak reactions for self-consistent NSE + /// + + + // initialize ydot_nuc to 0 + + for (int i = 1; i <= neqs; ++i) { + ydot_nuc(i) = 0.0_rt; + } + + rate_t rate_eval; + + [[maybe_unused]] amrex::Real rate, drate_dt, edot_nu, edot_gamma; + [[maybe_unused]] amrex::Real rhoy = state.rho * state.y_e; + + rate_eval.enuc_weak = 0.0_rt; + + // Calculate tabular rates and get ydot_weak + + + ydot_nuc(H1) = 0.0_rt; + + ydot_nuc(He4) = 0.0_rt; + + ydot_nuc(C12) = 0.0_rt; + + ydot_nuc(N13) = 0.0_rt; + + ydot_nuc(N14) = 0.0_rt; + + ydot_nuc(O16) = 0.0_rt; + + ydot_nuc(F18) = 0.0_rt; + + ydot_nuc(Ne20) = 0.0_rt; + + ydot_nuc(Ne21) = 0.0_rt; + + ydot_nuc(Na22) = 0.0_rt; + + ydot_nuc(Na23) = 0.0_rt; + + ydot_nuc(Mg24) = 0.0_rt; + + ydot_nuc(Al27) = 0.0_rt; + + ydot_nuc(Si28) = 0.0_rt; + + ydot_nuc(P31) = 0.0_rt; + + ydot_nuc(S32) = 0.0_rt; + + ydot_nuc(Ar36) = 0.0_rt; + + ydot_nuc(Ca40) = 0.0_rt; + + ydot_nuc(Ti44) = 0.0_rt; + + ydot_nuc(Cr48) = 0.0_rt; + + ydot_nuc(Fe52) = 0.0_rt; + + ydot_nuc(Ni56) = 0.0_rt; + + enuc_weak = rate_eval.enuc_weak; +} +#endif + + AMREX_GPU_HOST_DEVICE AMREX_INLINE void rhs_nuc(const burn_t& state, amrex::Array1D& ydot_nuc, diff --git a/nse_solver/Make.package b/nse_solver/Make.package index 66fe590429..f5344d27ff 100644 --- a/nse_solver/Make.package +++ b/nse_solver/Make.package @@ -1,4 +1,5 @@ ifeq ($(USE_NSE_NET), TRUE) CEXE_headers += nse_solver.H CEXE_headers += nse_check.H + CEXE_headers += nse_eos.H endif diff --git a/nse_solver/nse_eos.H b/nse_solver/nse_eos.H new file mode 100644 index 0000000000..69d4fccc44 --- /dev/null +++ b/nse_solver/nse_eos.H @@ -0,0 +1,174 @@ +#ifndef NSE_EOS_H +#define NSE_EOS_H + +#include + +#include + +#include +#include + + +AMREX_GPU_HOST_DEVICE AMREX_INLINE +amrex::Real nse_abar(const amrex::Real T, const amrex::Real rho, + const amrex::Real Ye, amrex::Real &mu_p, + amrex::Real &mu_n) { + /// + /// This function calculates abar from NSE using + /// Temp, rho, and Ye + /// + + burn_t burn_state; + burn_state.rho = rho; + burn_state.y_e = Ye; + burn_state.T = T; + burn_state.mu_p = mu_p; + burn_state.mu_n = mu_n; + + auto nse_state = get_actual_nse_state(burn_state, 1.0e-10_rt, true); + + amrex::Real abar{0.0_rt}; + for (int n = 0; n < NumSpec; ++n) { + abar += nse_state.xn[n] * aion_inv[n]; + } + + // update mu_p and mu_n + + mu_p = burn_state.mu_p; + mu_n = burn_state.mu_n; + + return abar; +} + + +AMREX_GPU_HOST_DEVICE AMREX_INLINE +amrex::Real nse_dabar_dT(const amrex::Real T, const amrex::Real rho, + const amrex::Real Ye, amrex::Real &mu_p, + amrex::Real &mu_n) { + /// + /// This function constructs dabar_dT + /// This should be 2nd order accurate + /// + + // deviation in temperature + + const amrex::Real dT = 1.0e-6_rt * T; + + // Calculate derivative using five-point stencil method + + amrex::Real dabar_dT = (-nse_abar(T + 2.0_rt*dT, rho, Ye, mu_p, mu_n) + + 8.0_rt * nse_abar(T+dT, rho, Ye, mu_p, mu_n) - + 8.0_rt * nse_abar(T-dT, rho, Ye, mu_p, mu_n) + + nse_abar(T - 2.0_rt*dT, rho, Ye, mu_p, mu_n) + ) / (12.0_rt * dT); + + // Calculate derivative using central differencing + + // amrex::Real dabar_dT = 0.5_rt * (nse_abar(T + dT, rho, Ye, mu_p, mu_n) - + // nse_abar(T - dT, rho, Ye, mu_p, mu_n)) / dT; + + return dabar_dT; +} + +/// +/// This function inverts this form of the EOS to find the T +/// that satisfies the EOS and NSE given an input e and rho. +/// +/// if we are in NSE, then the entire thermodynamic state is just +/// a function of rho, T, Ye. We can write the energy as: +/// +/// e = e(rho, T, Y_e, Abar(rho, T, Ye)) +/// +/// where we note that Abar is a function of those same inputs. +/// +/// The basic idea is that Abar and Zbar are both functions of +/// rho, T, Ye through NSE calculations, so we express the energy +/// as: +/// +/// e = e(rho, T, Abar(rho, T, Ye), Zbar(rho, T, Ye) +/// +/// and NR on that. Note that Zbar = Ye Abar, so we can group +/// those derivative terms together. +/// +/// T and abar come in as initial guesses and are updated +/// on output +/// + +AMREX_GPU_HOST_DEVICE AMREX_INLINE +void +nse_T_abar_from_e(const amrex::Real rho, const amrex::Real e_in, + const amrex::Real Ye, + amrex::Real &T, amrex::Real &abar, + amrex::Real &mu_p, amrex::Real &mu_n) { + + constexpr Real ttol{1.e-8_rt}; + constexpr int max_iter{100}; + + bool converged{false}; + + int iter{0}; + + // initialize burn_state + burn_t burn_state; + burn_state.rho = rho; + burn_state.y_e = Ye; + burn_state.mu_p = mu_p; + burn_state.mu_n = mu_n; + + while (not converged && iter < max_iter) { + + // update Temperature + + burn_state.T = T; + + auto nse_state = get_actual_nse_state(burn_state, 1.0e-10_rt, true); + + // call the EOS with the initial guess for T + // as well as density and NSE mass fractions + + eos_extra_t eos_state; + for (int n = 0; n < NumSpec; ++n) { + eos_state.xn[n] = nse_state.xn[n]; + } + + // Call the EOS to get internal energy + // abar, zbar, etc are calculated automatically in EOS + + eos_state.rho = rho; + eos_state.T = T; + eos(eos_input_rt, eos_state); + + // f is the quantity we want to zero + + amrex::Real f = eos_state.e - e_in; + + amrex::Real dabar_dT = nse_dabar_dT(T, rho, Ye, + nse_state.mu_p, nse_state.mu_n); + + // compute the correction to our guess + + Real dT = -f / (eos_state.dedT + eos_state.dedA * dabar_dT + + Ye * eos_state.dedZ * dabar_dT); + + // update the temperature and abar + + T = std::clamp(T + dT, 0.25 * T, 4.0 * T); + abar = eos_state.abar; + + // update mu_p and mu_n + + mu_p = nse_state.mu_p; + mu_n = nse_state.mu_n; + + // check convergence + + if (std::abs(dT) < ttol * T) { + converged = true; + } + iter++; + + } + +} + +#endif From 24fc8547090297bacdefe30f14bc2c37e53069f4 Mon Sep 17 00:00:00 2001 From: Zhi Chen <62574124+zhichen3@users.noreply.github.com> Date: Sun, 14 Jul 2024 20:53:47 -0400 Subject: [PATCH 4/5] include higher taylor terms for fast_atan (#1611) Pointed out by Eric: in screening, atan is used in cases like gamma - atan(gamma), so we should use a higher order term rather than simply letting atan(gamma) ~ gamma when gamma is small. --- .../ci-benchmarks/nse_net_unit_test.out | 52 +- unit_test/test_nse_net/make_table/burn_cell.H | 10 +- .../nse_net_make_table_unit_test.out | 704 +++++++++--------- util/approx_math/approx_math.H | 25 +- 4 files changed, 399 insertions(+), 392 deletions(-) diff --git a/unit_test/test_nse_net/ci-benchmarks/nse_net_unit_test.out b/unit_test/test_nse_net/ci-benchmarks/nse_net_unit_test.out index 18e17a9981..2af0ee8961 100644 --- a/unit_test/test_nse_net/ci-benchmarks/nse_net_unit_test.out +++ b/unit_test/test_nse_net/ci-benchmarks/nse_net_unit_test.out @@ -1,5 +1,5 @@ -Initializing AMReX (23.07-395-ge6c93bf22695)... -AMReX (23.07-395-ge6c93bf22695) initialized +Initializing AMReX (23.07-445-gc09c45c31007)... +AMReX (23.07-445-gc09c45c31007) initialized starting the single zone burn... reading in network electron-capture / beta-decay tables... chemical potential of proton is -3 @@ -8,28 +8,28 @@ State Density (g/cm^3): 10000000 State Temperature (K): 6000000000 electron fraction is 0.5 NSE state: -n : 0.000607569351 -H1 : 0.001992432526 -He4 : 0.5190078752 -C12 : 1.326101988e-05 -N13 : 9.477376579e-11 -N14 : 2.56548683e-09 -O16 : 3.225455214e-05 -F18 : 5.422064122e-11 -Ne20 : 7.427509352e-07 -Ne21 : 4.99891257e-08 -Na22 : 2.608008157e-09 -Na23 : 1.008842941e-06 -Mg24 : 0.0001022222498 -Al27 : 0.0005548826564 -Si28 : 0.03680966808 -P31 : 0.0422922376 -S32 : 0.0386577195 -Ar36 : 0.02464742425 -Ca40 : 0.0288515261 -Ti44 : 0.001661231048 -Cr48 : 0.009967017622 -Fe52 : 0.05968005276 -Ni56 : 0.2351208186 +n : 0.0006075779147 +H1 : 0.001992190328 +He4 : 0.519153618 +C12 : 1.326242057e-05 +N13 : 9.476378074e-11 +N14 : 2.565252696e-09 +O16 : 3.225118799e-05 +F18 : 5.421427915e-11 +Ne20 : 7.426541003e-07 +Ne21 : 4.998331295e-08 +Na22 : 2.607634144e-09 +Na23 : 1.008712481e-06 +Mg24 : 0.0001022062575 +Al27 : 0.0005547964335 +Si28 : 0.03680294966 +P31 : 0.0422845632 +S32 : 0.03864965592 +Ar36 : 0.02464164049 +Ca40 : 0.02884400365 +Ti44 : 0.001660754609 +Cr48 : 0.009963899266 +Fe52 : 0.05965982508 +Ni56 : 0.2350349989 We're in NSE. -AMReX (23.07-395-ge6c93bf22695) finalized +AMReX (23.07-445-gc09c45c31007) finalized diff --git a/unit_test/test_nse_net/make_table/burn_cell.H b/unit_test/test_nse_net/make_table/burn_cell.H index 886338bf59..17f056ef06 100644 --- a/unit_test/test_nse_net/make_table/burn_cell.H +++ b/unit_test/test_nse_net/make_table/burn_cell.H @@ -37,9 +37,13 @@ void burn_cell_c() state.rho = rho; state.y_e = Ye; - if (state.y_e > 0.52_rt){ - state.mu_p = -4.0_rt; - state.mu_n = -14.0_rt; + if (state.y_e >= 0.63_rt) { + state.mu_p = -3.0_rt; + state.mu_n = -14.9_rt; + } + else if (state.y_e > 0.52_rt){ + state.mu_p = -6.0_rt; + state.mu_n = -11.0_rt; } else if (state.y_e > 0.48_rt){ state.mu_p = -8.0_rt; diff --git a/unit_test/test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out b/unit_test/test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out index 7ec0293319..b6c5f387b3 100644 --- a/unit_test/test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out +++ b/unit_test/test_nse_net/make_table/ci-benchmarks/nse_net_make_table_unit_test.out @@ -1,355 +1,355 @@ -Initializing AMReX (23.07-395-ge6c93bf22695)... -AMReX (23.07-395-ge6c93bf22695) initialized +Initializing AMReX (23.07-445-gc09c45c31007)... +AMReX (23.07-445-gc09c45c31007) initialized starting the single zone burn... reading in network electron-capture / beta-decay tables... - 1.0000000000e+06 4.0000000000e+09 0.4000000000 -12.9555025204 -4.7372062842 0 - 1.0000000000e+06 4.4286927156e+09 0.4000000000 -12.4096892935 -5.3017297034 1 - 1.0000000000e+06 4.9033297922e+09 0.4000000000 -11.8124863350 -5.9247157509 1 - 1.0000000000e+06 5.4288352332e+09 0.4000000000 -11.2186398049 -6.5894751704 1 - 1.0000000000e+06 6.0106607628e+09 0.4000000000 -10.8221627699 -7.3603511111 1 - 1.0000000000e+06 6.6548423840e+09 0.4000000000 -10.4363763722 -8.2252373664 1 - 1.0000000000e+06 7.3680629973e+09 0.4000000000 -10.1445056150 -9.0849662772 1 - 1.0000000000e+06 8.1577217310e+09 0.4000000000 -10.2582950656 -9.7983309243 1 - 1.0000000000e+06 9.0320107014e+09 0.4000000000 -11.0251962339 -10.7011623998 1 - 1.0000000000e+06 1.0000000000e+10 0.4000000000 -12.3061952746 -11.9583139126 1 - 1.0000000000e+07 4.0000000000e+09 0.4000000000 -13.7539410035 -3.9436985166 0 - 1.0000000000e+07 4.4286927156e+09 0.4000000000 -13.2910063988 -4.4245064973 1 - 1.0000000000e+07 4.9033297922e+09 0.4000000000 -12.7738571048 -4.9622430425 1 - 1.0000000000e+07 5.4288352332e+09 0.4000000000 -12.2022519657 -5.5597835137 1 - 1.0000000000e+07 6.0106607628e+09 0.4000000000 -11.5968451667 -6.2091582807 1 - 1.0000000000e+07 6.6548423840e+09 0.4000000000 -11.0838653206 -6.9172956472 1 - 1.0000000000e+07 7.3680629973e+09 0.4000000000 -10.7117195040 -7.7493017684 1 - 1.0000000000e+07 8.1577217310e+09 0.4000000000 -10.3608224687 -8.6311809609 1 - 1.0000000000e+07 9.0320107014e+09 0.4000000000 -10.2512028838 -9.4182840030 1 - 1.0000000000e+07 1.0000000000e+10 0.4000000000 -10.6008153728 -10.1495273727 1 - 1.0000000000e+08 4.0000000000e+09 0.4000000000 -14.5607973487 -3.1500296846 0 - 1.0000000000e+08 4.4286927156e+09 0.4000000000 -14.1824332333 -3.5459083093 1 - 1.0000000000e+08 4.9033297922e+09 0.4000000000 -13.7577662180 -3.9903922193 1 - 1.0000000000e+08 5.4288352332e+09 0.4000000000 -13.2819583494 -4.4889815352 1 - 1.0000000000e+08 6.0106607628e+09 0.4000000000 -12.7514980414 -5.0464198084 1 - 1.0000000000e+08 6.6548423840e+09 0.4000000000 -12.1712468695 -5.6630196541 1 - 1.0000000000e+08 7.3680629973e+09 0.4000000000 -11.5722427319 -6.3289690263 1 - 1.0000000000e+08 8.1577217310e+09 0.4000000000 -11.0954919155 -7.0732199085 1 - 1.0000000000e+08 9.0320107014e+09 0.4000000000 -10.7313790718 -7.9307043937 1 - 1.0000000000e+08 1.0000000000e+10 0.4000000000 -10.4459234786 -8.7985272846 1 + 1.0000000000e+06 4.0000000000e+09 0.4000000000 -12.9555025985 -4.7372062322 0 + 1.0000000000e+06 4.4286927156e+09 0.4000000000 -12.4096886951 -5.3017300984 1 + 1.0000000000e+06 4.9033297922e+09 0.4000000000 -11.8124827277 -5.9247179903 1 + 1.0000000000e+06 5.4288352332e+09 0.4000000000 -11.2186195287 -6.5894809028 1 + 1.0000000000e+06 6.0106607628e+09 0.4000000000 -10.8221238417 -7.3603510805 1 + 1.0000000000e+06 6.6548423840e+09 0.4000000000 -10.4363399360 -8.2252368745 1 + 1.0000000000e+06 7.3680629973e+09 0.4000000000 -10.1444749829 -9.0849628652 1 + 1.0000000000e+06 8.1577217310e+09 0.4000000000 -10.2582729101 -9.7983255913 1 + 1.0000000000e+06 9.0320107014e+09 0.4000000000 -11.0251835300 -10.7011614674 1 + 1.0000000000e+06 1.0000000000e+10 0.4000000000 -12.3061845155 -11.9583138919 1 + 1.0000000000e+07 4.0000000000e+09 0.4000000000 -13.7539410130 -3.9436985103 0 + 1.0000000000e+07 4.4286927156e+09 0.4000000000 -13.2910064791 -4.4245064438 1 + 1.0000000000e+07 4.9033297922e+09 0.4000000000 -12.7738576654 -4.9622426711 1 + 1.0000000000e+07 5.4288352332e+09 0.4000000000 -12.2022552335 -5.5597814123 1 + 1.0000000000e+07 6.0106607628e+09 0.4000000000 -11.5968613147 -6.2091497038 1 + 1.0000000000e+07 6.6548423840e+09 0.4000000000 -11.0839445473 -6.9172943556 1 + 1.0000000000e+07 7.3680629973e+09 0.4000000000 -10.7118071972 -7.7493029627 1 + 1.0000000000e+07 8.1577217310e+09 0.4000000000 -10.3609046402 -8.6311911606 1 + 1.0000000000e+07 9.0320107014e+09 0.4000000000 -10.2511246531 -9.4182694259 1 + 1.0000000000e+07 1.0000000000e+10 0.4000000000 -10.6007602543 -10.1495148201 1 + 1.0000000000e+08 4.0000000000e+09 0.4000000000 -14.5607973496 -3.1500296840 0 + 1.0000000000e+08 4.4286927156e+09 0.4000000000 -14.1824332412 -3.5459083040 1 + 1.0000000000e+08 4.9033297922e+09 0.4000000000 -13.7577662756 -3.9903921809 1 + 1.0000000000e+08 5.4288352332e+09 0.4000000000 -13.2819587028 -4.4889813003 1 + 1.0000000000e+08 6.0106607628e+09 0.4000000000 -12.7514998816 -5.0464186008 1 + 1.0000000000e+08 6.6548423840e+09 0.4000000000 -12.1712550423 -5.6630145848 1 + 1.0000000000e+08 7.3680629973e+09 0.4000000000 -11.5722750148 -6.3289545463 1 + 1.0000000000e+08 8.1577217310e+09 0.4000000000 -11.0956062224 -7.0732197715 1 + 1.0000000000e+08 9.0320107014e+09 0.4000000000 -10.7314991991 -7.9307107472 1 + 1.0000000000e+08 1.0000000000e+10 0.4000000000 -10.4460221082 -8.7985576148 1 1.0000000000e+09 4.0000000000e+09 0.4000000000 -15.3871516246 -2.3563444183 0 - 1.0000000000e+09 4.4286927156e+09 0.4000000000 -15.0933118585 -2.6671723129 1 - 1.0000000000e+09 4.9033297922e+09 0.4000000000 -14.7621490253 -3.0175641915 1 - 1.0000000000e+09 5.4288352332e+09 0.4000000000 -14.3892566411 -3.4124048451 1 - 1.0000000000e+09 6.0106607628e+09 0.4000000000 -13.9696328335 -3.8570453479 1 - 1.0000000000e+09 6.6548423840e+09 0.4000000000 -13.4990892461 -4.3569663734 1 - 1.0000000000e+09 7.3680629973e+09 0.4000000000 -12.9751070035 -4.9165801598 1 - 1.0000000000e+09 8.1577217310e+09 0.4000000000 -12.4029615922 -5.5362921583 1 - 1.0000000000e+09 9.0320107014e+09 0.4000000000 -11.8075492548 -6.2097916688 1 - 1.0000000000e+09 1.0000000000e+10 0.4000000000 -11.2831770772 -6.9503077459 1 - 1.0000000000e+10 4.0000000000e+09 0.4000000000 -16.2582598176 -1.5626567041 0 - 1.0000000000e+10 4.4286927156e+09 0.4000000000 -16.0486738826 -1.7884208125 1 - 1.0000000000e+10 4.9033297922e+09 0.4000000000 -15.8108280448 -2.0446360334 1 - 1.0000000000e+10 5.4288352332e+09 0.4000000000 -15.5412395186 -2.3352407233 1 - 1.0000000000e+10 6.0106607628e+09 0.4000000000 -15.2357440992 -2.6646600961 1 - 1.0000000000e+10 6.6548423840e+09 0.4000000000 -14.8903553277 -3.0378300141 1 - 1.0000000000e+10 7.3680629973e+09 0.4000000000 -14.5003230778 -3.4601446829 1 - 1.0000000000e+10 8.1577217310e+09 0.4000000000 -14.0614570153 -3.9371817304 1 - 1.0000000000e+10 9.0320107014e+09 0.4000000000 -13.5711724217 -4.4739432697 1 - 1.0000000000e+10 1.0000000000e+10 0.4000000000 -13.0312175162 -5.0733655763 1 - 1.0000000000e+06 4.0000000000e+09 0.4500000000 -12.6199618053 -5.0493825987 0 - 1.0000000000e+06 4.4286927156e+09 0.4500000000 -12.0411581961 -5.6445822151 1 - 1.0000000000e+06 4.9033297922e+09 0.4500000000 -11.4232312199 -6.2868479094 1 - 1.0000000000e+06 5.4288352332e+09 0.4500000000 -10.8518428398 -6.9298657584 1 - 1.0000000000e+06 6.0106607628e+09 0.4500000000 -10.4348708730 -7.7172921982 1 - 1.0000000000e+06 6.6548423840e+09 0.4500000000 -10.0424621871 -8.5877955647 1 - 1.0000000000e+06 7.3680629973e+09 0.4500000000 -9.8849946892 -9.3214605955 1 - 1.0000000000e+06 8.1577217310e+09 0.4500000000 -10.1350747780 -9.9053651295 1 - 1.0000000000e+06 9.0320107014e+09 0.4500000000 -10.9317730375 -10.7720016856 1 - 1.0000000000e+06 1.0000000000e+10 0.4500000000 -12.2046629248 -12.0333675970 1 - 1.0000000000e+07 4.0000000000e+09 0.4500000000 -13.4182964409 -4.2562148091 0 - 1.0000000000e+07 4.4286927156e+09 0.4500000000 -12.9196141769 -4.7702516152 1 - 1.0000000000e+07 4.9033297922e+09 0.4500000000 -12.3646700494 -5.3431239685 1 - 1.0000000000e+07 5.4288352332e+09 0.4500000000 -11.7607919864 -5.9706982948 1 - 1.0000000000e+07 6.0106607628e+09 0.4500000000 -11.1529755129 -6.6222947524 1 - 1.0000000000e+07 6.6548423840e+09 0.4500000000 -10.6533363165 -7.3144541390 1 - 1.0000000000e+07 7.3680629973e+09 0.4500000000 -10.2538115130 -8.1710937198 1 - 1.0000000000e+07 8.1577217310e+09 0.4500000000 -9.9753897162 -8.9845986143 1 - 1.0000000000e+07 9.0320107014e+09 0.4500000000 -10.0361542054 -9.6115697958 1 - 1.0000000000e+07 1.0000000000e+10 0.4500000000 -10.4770937073 -10.2529776199 1 - 1.0000000000e+08 4.0000000000e+09 0.4500000000 -14.2257599760 -3.4625787945 0 - 1.0000000000e+08 4.4286927156e+09 0.4500000000 -13.8113628847 -3.8919431610 1 - 1.0000000000e+08 4.9033297922e+09 0.4500000000 -13.3469605437 -4.3733521393 1 - 1.0000000000e+08 5.4288352332e+09 0.4500000000 -12.8281601561 -4.9118963225 1 - 1.0000000000e+08 6.0106607628e+09 0.4500000000 -12.2549827771 -5.5090675693 1 - 1.0000000000e+08 6.6548423840e+09 0.4500000000 -11.6452553291 -6.1532183846 1 - 1.0000000000e+08 7.3680629973e+09 0.4500000000 -11.0519169747 -6.8138046108 1 - 1.0000000000e+08 8.1577217310e+09 0.4500000000 -10.5744119695 -7.5537001211 1 - 1.0000000000e+08 9.0320107014e+09 0.4500000000 -10.2097691952 -8.4106502242 1 - 1.0000000000e+08 1.0000000000e+10 0.4500000000 -10.0527955898 -9.1577337848 1 - 1.0000000000e+09 4.0000000000e+09 0.4500000000 -15.0536228396 -2.6688948472 0 - 1.0000000000e+09 4.4286927156e+09 0.4500000000 -14.7237134072 -3.0132309283 1 - 1.0000000000e+09 4.9033297922e+09 0.4500000000 -14.3526160867 -3.4007197563 1 - 1.0000000000e+09 5.4288352332e+09 0.4500000000 -13.9355877638 -3.8365577467 1 - 1.0000000000e+09 6.0106607628e+09 0.4500000000 -13.4675881563 -4.3261360138 1 - 1.0000000000e+09 6.6548423840e+09 0.4500000000 -12.9456617325 -4.8738185708 1 - 1.0000000000e+09 7.3680629973e+09 0.4500000000 -12.3724764785 -5.4792958152 1 - 1.0000000000e+09 8.1577217310e+09 0.4500000000 -11.7661588560 -6.1311892901 1 - 1.0000000000e+09 9.0320107014e+09 0.4500000000 -11.1668612132 -6.8085657276 1 - 1.0000000000e+09 1.0000000000e+10 0.4500000000 -10.6512855625 -7.5351949184 1 - 1.0000000000e+10 4.0000000000e+09 0.4500000000 -15.9281485305 -1.8752022445 0 - 1.0000000000e+10 4.4286927156e+09 0.4500000000 -15.6824856145 -2.1344670415 1 - 1.0000000000e+10 4.9033297922e+09 0.4500000000 -15.4047001813 -2.4277705489 1 - 1.0000000000e+10 5.4288352332e+09 0.4500000000 -15.0909188396 -2.7594143274 1 - 1.0000000000e+10 6.0106607628e+09 0.4500000000 -14.7366091103 -3.1341694060 1 - 1.0000000000e+10 6.6548423840e+09 0.4500000000 -14.3375648682 -3.5571890304 1 - 1.0000000000e+10 7.3680629973e+09 0.4500000000 -13.8893261475 -4.0336347952 1 - 1.0000000000e+10 8.1577217310e+09 0.4500000000 -13.3893312974 -4.5676563215 1 - 1.0000000000e+10 9.0320107014e+09 0.4500000000 -12.8393727587 -5.1603436508 1 - 1.0000000000e+10 1.0000000000e+10 0.4500000000 -12.2492192359 -5.8073345674 1 - 1.0000000000e+06 4.0000000000e+09 0.5000000000 -8.0664287063 -9.4890300885 0 - 1.0000000000e+06 4.4286927156e+09 0.5000000000 -8.2111123981 -9.3757297125 1 - 1.0000000000e+06 4.9033297922e+09 0.5000000000 -8.3918399895 -9.2344583793 1 - 1.0000000000e+06 5.4288352332e+09 0.5000000000 -8.7773708532 -8.9447667662 1 - 1.0000000000e+06 6.0106607628e+09 0.5000000000 -9.0637265807 -9.0647441208 1 - 1.0000000000e+06 6.6548423840e+09 0.5000000000 -9.3058731093 -9.3070059375 1 - 1.0000000000e+06 7.3680629973e+09 0.5000000000 -9.5983915599 -9.5996554213 1 - 1.0000000000e+06 8.1577217310e+09 0.5000000000 -10.0168541391 -10.0182626012 1 - 1.0000000000e+06 9.0320107014e+09 0.5000000000 -10.8474197344 -10.8489878418 1 - 1.0000000000e+06 1.0000000000e+10 0.5000000000 -12.1138199092 -12.1155643414 1 - 1.0000000000e+07 4.0000000000e+09 0.5000000000 -7.8839123533 -9.6536792188 0 - 1.0000000000e+07 4.4286927156e+09 0.5000000000 -8.0019180002 -9.5630334034 1 - 1.0000000000e+07 4.9033297922e+09 0.5000000000 -8.1354121681 -9.4609849615 1 - 1.0000000000e+07 5.4288352332e+09 0.5000000000 -8.2961986080 -9.3384887105 1 - 1.0000000000e+07 6.0106607628e+09 0.5000000000 -8.5443326081 -9.1484326209 1 - 1.0000000000e+07 6.6548423840e+09 0.5000000000 -8.9669523913 -8.9714790676 1 - 1.0000000000e+07 7.3680629973e+09 0.5000000000 -9.2004322394 -9.2015703723 1 - 1.0000000000e+07 8.1577217310e+09 0.5000000000 -9.4726257222 -9.4739144064 1 - 1.0000000000e+07 9.0320107014e+09 0.5000000000 -9.8194771269 -9.8209313702 1 - 1.0000000000e+07 1.0000000000e+10 0.5000000000 -10.3609081088 -10.3625443016 1 - 1.0000000000e+08 4.0000000000e+09 0.5000000000 -7.7564252605 -9.7786841281 0 - 1.0000000000e+08 4.4286927156e+09 0.5000000000 -7.8492019112 -9.7096717498 1 - 1.0000000000e+08 4.9033297922e+09 0.5000000000 -7.9527961073 -9.6330305105 1 - 1.0000000000e+08 5.4288352332e+09 0.5000000000 -8.0692464391 -9.5476093419 1 - 1.0000000000e+08 6.0106607628e+09 0.5000000000 -8.2063900317 -9.4473474639 1 - 1.0000000000e+08 6.6548423840e+09 0.5000000000 -8.3860213927 -9.3159021644 1 - 1.0000000000e+08 7.3680629973e+09 0.5000000000 -8.6892088477 -9.0970021780 1 - 1.0000000000e+08 8.1577217310e+09 0.5000000000 -9.0475207044 -9.0501554594 1 - 1.0000000000e+08 9.0320107014e+09 0.5000000000 -9.2986740372 -9.2997699354 1 - 1.0000000000e+08 1.0000000000e+10 0.5000000000 -9.5981877903 -9.5994826015 1 - 1.0000000000e+09 4.0000000000e+09 0.5000000000 -7.7470638015 -9.8201507833 0 - 1.0000000000e+09 4.4286927156e+09 0.5000000000 -7.8147769492 -9.7725305169 1 - 1.0000000000e+09 4.9033297922e+09 0.5000000000 -7.8899734067 -9.7200977838 1 - 1.0000000000e+09 5.4288352332e+09 0.5000000000 -7.9729788613 -9.6630648524 1 - 1.0000000000e+09 6.0106607628e+09 0.5000000000 -8.0659663127 -9.5999364916 1 - 1.0000000000e+09 6.6548423840e+09 0.5000000000 -8.1710732556 -9.5300539621 1 - 1.0000000000e+09 7.3680629973e+09 0.5000000000 -8.2980125161 -9.4462720608 1 - 1.0000000000e+09 8.1577217310e+09 0.5000000000 -8.4687660372 -9.3330282318 1 - 1.0000000000e+09 9.0320107014e+09 0.5000000000 -8.7320500861 -9.1630033338 1 - 1.0000000000e+09 1.0000000000e+10 0.5000000000 -9.0632943617 -9.0828318763 1 - 1.0000000000e+10 4.0000000000e+09 0.5000000000 -7.9879953396 -9.6879468798 1 - 1.0000000000e+10 4.4286927156e+09 0.5000000000 -8.0302314490 -9.6619337577 1 - 1.0000000000e+10 4.9033297922e+09 0.5000000000 -8.0763284304 -9.6342740546 1 - 1.0000000000e+10 5.4288352332e+09 0.5000000000 -8.1258459897 -9.6058175046 1 - 1.0000000000e+10 6.0106607628e+09 0.5000000000 -8.1795088599 -9.5762491105 1 - 1.0000000000e+10 6.6548423840e+09 0.5000000000 -8.2365078557 -9.5472461611 1 - 1.0000000000e+10 7.3680629973e+09 0.5000000000 -8.2988571338 -9.5177368349 1 - 1.0000000000e+10 8.1577217310e+09 0.5000000000 -8.3695493416 -9.4863947153 1 - 1.0000000000e+10 9.0320107014e+09 0.5000000000 -8.4559431120 -9.4487939586 1 - 1.0000000000e+10 1.0000000000e+10 0.5000000000 -8.5723779233 -9.3960835543 1 - 1.0000000000e+06 4.0000000000e+09 0.5500000000 -4.9263567583 -12.6306847771 0 - 1.0000000000e+06 4.4286927156e+09 0.5500000000 -5.5125982821 -12.0759967435 1 - 1.0000000000e+06 4.9033297922e+09 0.5500000000 -6.1679059616 -11.4605293767 1 - 1.0000000000e+06 5.4288352332e+09 0.5500000000 -6.9002510955 -10.8361899606 1 - 1.0000000000e+06 6.0106607628e+09 0.5500000000 -7.7162795920 -10.4358987387 1 - 1.0000000000e+06 6.6548423840e+09 0.5500000000 -8.5866665410 -10.0436084115 1 - 1.0000000000e+06 7.3680629973e+09 0.5500000000 -9.3202024829 -9.8862686067 1 - 1.0000000000e+06 8.1577217310e+09 0.5500000000 -9.9039624201 -10.1364903176 1 - 1.0000000000e+06 9.0320107014e+09 0.5500000000 -10.7704366546 -10.9333444696 1 - 1.0000000000e+06 1.0000000000e+10 0.5500000000 -12.0316256213 -12.2064099419 1 - 1.0000000000e+07 4.0000000000e+09 0.5500000000 -4.1328499283 -13.4067273031 0 - 1.0000000000e+07 4.4286927156e+09 0.5500000000 -4.6340211015 -12.9330486665 1 - 1.0000000000e+07 4.9033297922e+09 0.5500000000 -5.1951522570 -12.4035178568 1 - 1.0000000000e+07 5.4288352332e+09 0.5500000000 -5.8233286811 -11.8138168049 1 - 1.0000000000e+07 6.0106607628e+09 0.5500000000 -6.5263096328 -11.1696374790 1 - 1.0000000000e+07 6.6548423840e+09 0.5500000000 -7.3118205884 -10.6533164642 1 - 1.0000000000e+07 7.3680629973e+09 0.5500000000 -8.1699653523 -10.2549768923 1 - 1.0000000000e+07 8.1577217310e+09 0.5500000000 -8.9833218276 -9.9767011458 1 - 1.0000000000e+07 9.0320107014e+09 0.5500000000 -9.6101289590 -10.0376261258 1 - 1.0000000000e+07 1.0000000000e+10 0.5500000000 -10.2513560243 -10.4787469547 1 - 1.0000000000e+08 4.0000000000e+09 0.5500000000 -3.3394222565 -14.1986157657 0 - 1.0000000000e+08 4.4286927156e+09 0.5500000000 -3.7555058127 -13.8064185762 1 - 1.0000000000e+08 4.9033297922e+09 0.5500000000 -4.2224338348 -13.3665796932 1 - 1.0000000000e+08 5.4288352332e+09 0.5500000000 -4.7463244462 -12.8738625707 1 - 1.0000000000e+08 6.0106607628e+09 0.5500000000 -5.3340170370 -12.3231580847 1 - 1.0000000000e+08 6.6548423840e+09 0.5500000000 -5.9930630460 -11.7121630303 1 - 1.0000000000e+08 7.3680629973e+09 0.5500000000 -6.7314189359 -11.0586117888 1 - 1.0000000000e+08 8.1577217310e+09 0.5500000000 -7.5518998378 -10.5748765826 1 - 1.0000000000e+08 9.0320107014e+09 0.5500000000 -8.4095851780 -10.2109427835 1 - 1.0000000000e+08 1.0000000000e+10 0.5500000000 -9.1564767533 -10.0541540547 1 - 1.0000000000e+09 4.0000000000e+09 0.5500000000 -2.5470603946 -15.0251960133 0 - 1.0000000000e+09 4.4286927156e+09 0.5500000000 -2.8780137697 -14.7144476551 1 - 1.0000000000e+09 4.9033297922e+09 0.5500000000 -3.2506978413 -14.3646520162 1 - 1.0000000000e+09 5.4288352332e+09 0.5500000000 -3.6702538871 -13.9712053620 1 - 1.0000000000e+09 6.0106607628e+09 0.5500000000 -4.1424458459 -13.5290133725 1 - 1.0000000000e+09 6.6548423840e+09 0.5500000000 -4.6737326992 -13.0330537776 1 - 1.0000000000e+09 7.3680629973e+09 0.5500000000 -5.2713212052 -12.4785205588 1 - 1.0000000000e+09 8.1577217310e+09 0.5500000000 -5.9429385089 -11.8635582114 1 - 1.0000000000e+09 9.0320107014e+09 0.5500000000 -6.6949519374 -11.2026555889 1 - 1.0000000000e+09 1.0000000000e+10 0.5500000000 -7.5247904700 -10.6481506638 1 - 1.0000000000e+10 4.0000000000e+09 0.5500000000 -1.7575285290 -15.9281119586 1 - 1.0000000000e+10 4.4286927156e+09 0.5500000000 -2.0032273072 -15.6987367857 1 - 1.0000000000e+10 4.9033297922e+09 0.5500000000 -2.2815466531 -15.4389665570 1 - 1.0000000000e+10 5.4288352332e+09 0.5500000000 -2.5966508213 -15.1450469031 1 - 1.0000000000e+10 6.0106607628e+09 0.5500000000 -2.9532220861 -14.8127045193 1 - 1.0000000000e+10 6.6548423840e+09 0.5500000000 -3.3565239474 -14.4375393899 1 - 1.0000000000e+10 7.3680629973e+09 0.5500000000 -3.8124716671 -14.0145582613 1 - 1.0000000000e+10 8.1577217310e+09 0.5500000000 -4.3277071808 -13.5387148150 1 - 1.0000000000e+10 9.0320107014e+09 0.5500000000 -4.9096377922 -13.0053181043 1 - 1.0000000000e+10 1.0000000000e+10 0.5500000000 -5.5661097741 -12.4114337141 1 - 1.0000000000e+06 4.0000000000e+09 0.6000000000 -4.6874372498 -12.8713326831 0 - 1.0000000000e+06 4.4286927156e+09 0.6000000000 -5.2480726287 -12.3424870400 1 - 1.0000000000e+06 4.9033297922e+09 0.6000000000 -5.8750376557 -11.7561606661 1 - 1.0000000000e+06 5.4288352332e+09 0.6000000000 -6.5760898652 -11.1834318133 1 - 1.0000000000e+06 6.0106607628e+09 0.6000000000 -7.3593425290 -10.8232051193 1 - 1.0000000000e+06 6.6548423840e+09 0.6000000000 -8.2241108560 -10.4375374976 1 - 1.0000000000e+06 7.3680629973e+09 0.6000000000 -9.0837111804 -10.1457926302 1 - 1.0000000000e+06 8.1577217310e+09 0.6000000000 -9.7969327525 -10.2597190165 1 - 1.0000000000e+06 9.0320107014e+09 0.6000000000 -10.6996001970 -11.0267712343 1 - 1.0000000000e+06 1.0000000000e+10 0.6000000000 -11.9565742814 -12.3079450256 1 - 1.0000000000e+07 4.0000000000e+09 0.6000000000 -3.8939383372 -13.6477417081 0 - 1.0000000000e+07 4.4286927156e+09 0.6000000000 -4.3695027901 -13.1998278106 1 - 1.0000000000e+07 4.9033297922e+09 0.6000000000 -4.9022832206 -12.6988733383 1 - 1.0000000000e+07 5.4288352332e+09 0.6000000000 -5.4990786383 -12.1410754239 1 - 1.0000000000e+07 6.0106607628e+09 0.6000000000 -6.1674541533 -11.5347440371 1 - 1.0000000000e+07 6.6548423840e+09 0.6000000000 -6.9156298388 -11.0830753295 1 - 1.0000000000e+07 7.3680629973e+09 0.6000000000 -7.7481817497 -10.7129141315 1 - 1.0000000000e+07 8.1577217310e+09 0.6000000000 -8.6299116983 -10.3621615417 1 - 1.0000000000e+07 9.0320107014e+09 0.6000000000 -9.4168075327 -10.2525633805 1 - 1.0000000000e+07 1.0000000000e+10 0.6000000000 -10.1478709457 -10.6024085069 1 - 1.0000000000e+08 4.0000000000e+09 0.6000000000 -3.1005259273 -14.4405028547 0 - 1.0000000000e+08 4.4286927156e+09 0.6000000000 -3.4910018429 -14.0740545856 1 - 1.0000000000e+08 4.9033297922e+09 0.6000000000 -3.9295779936 -13.6627315310 1 - 1.0000000000e+08 5.4288352332e+09 0.6000000000 -4.4220785916 -13.2016151666 1 - 1.0000000000e+08 6.0106607628e+09 0.6000000000 -4.9750223848 -12.6860143020 1 - 1.0000000000e+08 6.6548423840e+09 0.6000000000 -5.5956833129 -12.1143859321 1 - 1.0000000000e+08 7.3680629973e+09 0.6000000000 -6.2923470974 -11.5100241165 1 - 1.0000000000e+08 8.1577217310e+09 0.6000000000 -7.0719798685 -11.0955902794 1 - 1.0000000000e+08 9.0320107014e+09 0.6000000000 -7.9296619962 -10.7326405101 1 - 1.0000000000e+08 1.0000000000e+10 0.6000000000 -8.7972934615 -10.4473618809 1 - 1.0000000000e+09 4.0000000000e+09 0.6000000000 -2.3082235668 -15.2690142179 0 - 1.0000000000e+09 4.4286927156e+09 0.6000000000 -2.6135665219 -14.9840071714 1 - 1.0000000000e+09 4.9033297922e+09 0.6000000000 -2.9578960346 -14.6627124682 1 - 1.0000000000e+09 5.4288352332e+09 0.6000000000 -3.3460592682 -14.3008241528 1 - 1.0000000000e+09 6.0106607628e+09 0.6000000000 -3.7834946071 -13.8935836540 1 - 1.0000000000e+09 6.6548423840e+09 0.6000000000 -4.2763022554 -13.4363612027 1 - 1.0000000000e+09 7.3680629973e+09 0.6000000000 -4.8313168013 -12.9248402164 1 - 1.0000000000e+09 8.1577217310e+09 0.6000000000 -5.4561228019 -12.3577854824 1 - 1.0000000000e+09 9.0320107014e+09 0.6000000000 -6.1587039833 -11.7516690135 1 - 1.0000000000e+09 1.0000000000e+10 0.6000000000 -6.9454524603 -11.2745753097 1 - 1.0000000000e+10 4.0000000000e+09 0.6000000000 -1.5188716438 -16.1761351108 0 - 1.0000000000e+10 4.4286927156e+09 0.6000000000 -1.7389525113 -15.9724978714 1 - 1.0000000000e+10 4.9033297922e+09 0.6000000000 -1.9889099833 -15.7412229679 1 - 1.0000000000e+10 5.4288352332e+09 0.6000000000 -2.2726142074 -15.4788512292 1 - 1.0000000000e+10 6.0106607628e+09 0.6000000000 -2.5944218278 -15.1814362147 1 - 1.0000000000e+10 6.6548423840e+09 0.6000000000 -2.9592358252 -14.8449412444 1 - 1.0000000000e+10 7.3680629973e+09 0.6000000000 -3.3725725553 -14.4647794265 1 - 1.0000000000e+10 8.1577217310e+09 0.6000000000 -3.8406360372 -14.0363659318 1 - 1.0000000000e+10 9.0320107014e+09 0.6000000000 -4.3703910764 -13.5555230831 1 - 1.0000000000e+10 1.0000000000e+10 0.6000000000 -4.9695596343 -13.0196507108 1 - 1.0000000000e+06 4.0000000000e+09 0.6500000000 -4.5476796742 -13.0130039947 0 - 1.0000000000e+06 4.4286927156e+09 0.6500000000 -5.0933362712 -12.4994174293 1 - 1.0000000000e+06 4.9033297922e+09 0.6500000000 -5.7037185995 -11.9307011924 1 - 1.0000000000e+06 5.4288352332e+09 0.6500000000 -6.3864250027 -11.4022341176 1 - 1.0000000000e+06 6.0106607628e+09 0.6500000000 -7.1497090923 -11.0673819184 1 - 1.0000000000e+06 6.6548423840e+09 0.6500000000 -7.9986396847 -10.7004385555 1 - 1.0000000000e+06 7.3680629973e+09 0.6500000000 -8.8932281298 -10.3705502418 1 - 1.0000000000e+06 8.1577217310e+09 0.6500000000 -9.6962470719 -10.3881605331 1 - 1.0000000000e+06 9.0320107014e+09 0.6500000000 -10.6340099963 -11.1317007071 1 - 1.0000000000e+06 1.0000000000e+10 0.6500000000 -11.8875232376 -12.4230320085 1 - 1.0000000000e+07 4.0000000000e+09 0.6500000000 -3.7541883557 -13.7897569631 0 - 1.0000000000e+07 4.4286927156e+09 0.6500000000 -4.2147736231 -13.3570050888 1 - 1.0000000000e+07 4.9033297922e+09 0.6500000000 -4.7309696318 -12.8729015431 1 - 1.0000000000e+07 5.4288352332e+09 0.6500000000 -5.3094047388 -12.3341136606 1 - 1.0000000000e+07 6.0106607628e+09 0.6500000000 -5.9574765589 -11.7528091119 1 - 1.0000000000e+07 6.6548423840e+09 0.6500000000 -6.6833511744 -11.3535641629 1 - 1.0000000000e+07 7.3680629973e+09 0.6500000000 -7.4942271332 -11.0088492097 1 - 1.0000000000e+07 8.1577217310e+09 0.6500000000 -8.3793161375 -10.6554902076 1 - 1.0000000000e+07 9.0320107014e+09 0.6500000000 -9.2440396955 -10.4605434295 1 - 1.0000000000e+07 1.0000000000e+10 0.6500000000 -10.0502569485 -10.7352236557 1 - 1.0000000000e+08 4.0000000000e+09 0.6500000000 -2.9607906232 -14.5833440181 0 - 1.0000000000e+08 4.4286927156e+09 0.6500000000 -3.3362864793 -14.2320397904 1 - 1.0000000000e+08 4.9033297922e+09 0.6500000000 -3.7582773168 -13.8374930298 1 - 1.0000000000e+08 5.4288352332e+09 0.6500000000 -4.2324153362 -13.3949878626 1 - 1.0000000000e+08 6.0106607628e+09 0.6500000000 -4.7650296834 -12.9001483764 1 - 1.0000000000e+08 6.6548423840e+09 0.6500000000 -5.3631981018 -12.3522858919 1 - 1.0000000000e+08 7.3680629973e+09 0.6500000000 -6.0351054434 -11.7844816604 1 - 1.0000000000e+08 8.1577217310e+09 0.6500000000 -6.7883154852 -11.4259937475 1 - 1.0000000000e+08 9.0320107014e+09 0.6500000000 -7.6264687221 -11.0863271605 1 - 1.0000000000e+08 1.0000000000e+10 0.6500000000 -8.5199483965 -10.7735064153 1 - 1.0000000000e+09 4.0000000000e+09 0.6500000000 -2.1685452785 -15.4136828685 0 - 1.0000000000e+09 4.4286927156e+09 0.6500000000 -2.4589055197 -15.1438122255 1 - 1.0000000000e+09 4.9033297922e+09 0.6500000000 -2.7866471569 -14.8392773913 1 - 1.0000000000e+09 5.4288352332e+09 0.6500000000 -3.1564453048 -14.4959492528 1 - 1.0000000000e+09 6.0106607628e+09 0.6500000000 -3.5735478746 -14.1092731402 1 - 1.0000000000e+09 6.6548423840e+09 0.6500000000 -4.0438444754 -13.6748739586 1 - 1.0000000000e+09 7.3680629973e+09 0.6500000000 -4.5739409491 -13.1888195133 1 - 1.0000000000e+09 8.1577217310e+09 0.6500000000 -5.1712163628 -12.6507561248 1 - 1.0000000000e+09 9.0320107014e+09 0.6500000000 -5.8437444564 -12.0829813034 1 - 1.0000000000e+09 1.0000000000e+10 0.6500000000 -6.5994405224 -11.6765113113 1 - 1.0000000000e+10 4.0000000000e+09 0.6500000000 -1.3793655235 -16.3247814756 0 - 1.0000000000e+10 4.4286927156e+09 0.6500000000 -1.5844565484 -16.1362774695 1 - 1.0000000000e+10 4.9033297922e+09 0.6500000000 -1.8178191791 -15.9217568381 1 - 1.0000000000e+10 5.4288352332e+09 0.6500000000 -2.0831515252 -15.6779340065 1 - 1.0000000000e+10 6.0106607628e+09 0.6500000000 -2.3846197513 -15.4010550224 1 - 1.0000000000e+10 6.6548423840e+09 0.6500000000 -2.7269159436 -15.0872983686 1 - 1.0000000000e+10 7.3680629973e+09 0.6500000000 -3.1153229999 -14.7323215524 1 - 1.0000000000e+10 8.1577217310e+09 0.6500000000 -3.5557871014 -14.3318302438 1 - 1.0000000000e+10 9.0320107014e+09 0.6500000000 -4.0549949963 -13.8820225244 1 - 1.0000000000e+10 1.0000000000e+10 0.6500000000 -4.6204256792 -13.3808512952 1 - 1.0000000000e+06 4.0000000000e+09 0.7000000000 -4.4485209891 -13.1143251183 0 - 1.0000000000e+06 4.4286927156e+09 0.7000000000 -4.9835499142 -12.6117081124 1 - 1.0000000000e+06 4.9033297922e+09 0.7000000000 -5.5821660601 -12.0561394505 1 - 1.0000000000e+06 5.4288352332e+09 0.7000000000 -6.2518496554 -11.5720368075 1 - 1.0000000000e+06 6.0106607628e+09 0.7000000000 -7.0008334984 -11.2561651451 1 - 1.0000000000e+06 6.6548423840e+09 0.7000000000 -7.8361333603 -10.9067344717 1 - 1.0000000000e+06 7.3680629973e+09 0.7000000000 -8.7392243787 -10.5687792383 1 - 1.0000000000e+06 8.1577217310e+09 0.7000000000 -9.6021916703 -10.5229827098 1 - 1.0000000000e+06 9.0320107014e+09 0.7000000000 -10.5729938190 -11.2517560750 1 - 1.0000000000e+06 1.0000000000e+10 0.7000000000 -11.8235850796 -12.5558647574 1 - 1.0000000000e+07 4.0000000000e+09 0.7000000000 -3.6550369761 -13.8914003282 0 - 1.0000000000e+07 4.4286927156e+09 0.7000000000 -4.1049942046 -13.4694919130 1 - 1.0000000000e+07 4.9033297922e+09 0.7000000000 -4.6094232304 -12.9974776854 1 - 1.0000000000e+07 5.4288352332e+09 0.7000000000 -5.1748308102 -12.4725550651 1 - 1.0000000000e+07 6.0106607628e+09 0.7000000000 -5.8084870128 -11.9131933897 1 - 1.0000000000e+07 6.6548423840e+09 0.7000000000 -6.5184625267 -11.5626253210 1 - 1.0000000000e+07 7.3680629973e+09 0.7000000000 -7.3127774909 -11.2390485970 1 - 1.0000000000e+07 8.1577217310e+09 0.7000000000 -8.1901397748 -10.8966710259 1 - 1.0000000000e+07 9.0320107014e+09 0.7000000000 -9.0920672122 -10.6613615921 1 - 1.0000000000e+07 1.0000000000e+10 0.7000000000 -9.9582875767 -10.8797496537 1 - 1.0000000000e+08 4.0000000000e+09 0.7000000000 -2.8616533998 -14.6857719715 0 - 1.0000000000e+08 4.4286927156e+09 0.7000000000 -3.2265203806 -14.3452898486 1 - 1.0000000000e+08 4.9033297922e+09 0.7000000000 -3.6367434104 -13.9627369387 1 - 1.0000000000e+08 5.4288352332e+09 0.7000000000 -4.0978526310 -13.5335594957 1 - 1.0000000000e+08 6.0106607628e+09 0.7000000000 -4.6160423318 -13.0536842248 1 - 1.0000000000e+08 6.6548423840e+09 0.7000000000 -5.1982455454 -12.5235029605 1 - 1.0000000000e+08 7.3680629973e+09 0.7000000000 -5.8525289991 -11.9932157006 1 - 1.0000000000e+08 8.1577217310e+09 0.7000000000 -6.5865650606 -11.6818564067 1 - 1.0000000000e+08 9.0320107014e+09 0.7000000000 -7.4069463043 -11.3651455549 1 - 1.0000000000e+08 1.0000000000e+10 0.7000000000 -8.3014810042 -11.0535468230 1 - 1.0000000000e+09 4.0000000000e+09 0.7000000000 -2.0694628674 -15.5178474652 0 - 1.0000000000e+09 4.4286927156e+09 0.7000000000 -2.3491916894 -15.2587909978 1 - 1.0000000000e+09 4.9033297922e+09 0.7000000000 -2.6651630630 -14.9662311049 1 - 1.0000000000e+09 5.4288352332e+09 0.7000000000 -3.0219300320 -14.6361670520 1 - 1.0000000000e+09 6.0106607628e+09 0.7000000000 -3.4246053577 -14.2641972416 1 - 1.0000000000e+09 6.6548423840e+09 0.7000000000 -3.8789291126 -13.8461550124 1 - 1.0000000000e+09 7.3680629973e+09 0.7000000000 -4.3913430475 -13.3784836554 1 - 1.0000000000e+09 8.1577217310e+09 0.7000000000 -4.9690605512 -12.8619701183 1 - 1.0000000000e+09 9.0320107014e+09 0.7000000000 -5.6200778085 -12.3291976596 1 - 1.0000000000e+09 1.0000000000e+10 0.7000000000 -6.3527292745 -11.9889184561 1 - 1.0000000000e+10 4.0000000000e+09 0.7000000000 -1.2804483941 -16.4327248844 0 - 1.0000000000e+10 4.4286927156e+09 0.7000000000 -1.4749011906 -16.2550319759 1 - 1.0000000000e+10 4.9033297922e+09 0.7000000000 -1.6964869006 -16.0524804719 1 - 1.0000000000e+10 5.4288352332e+09 0.7000000000 -1.9487815730 -15.8219089869 1 - 1.0000000000e+10 6.0106607628e+09 0.7000000000 -2.2358162273 -15.5597016802 1 - 1.0000000000e+10 6.6548423840e+09 0.7000000000 -2.5621333174 -15.2621928113 1 - 1.0000000000e+10 7.3680629973e+09 0.7000000000 -2.9328500820 -14.9252211738 1 - 1.0000000000e+10 8.1577217310e+09 0.7000000000 -3.3537294917 -14.5447208945 1 - 1.0000000000e+10 9.0320107014e+09 0.7000000000 -3.8312578887 -14.1172251767 1 - 1.0000000000e+10 1.0000000000e+10 0.7000000000 -4.3727140958 -13.6413288908 1 -AMReX (23.07-395-ge6c93bf22695) finalized + 1.0000000000e+09 4.4286927156e+09 0.4000000000 -15.0933118588 -2.6671723127 1 + 1.0000000000e+09 4.9033297922e+09 0.4000000000 -14.7621490285 -3.0175641894 1 + 1.0000000000e+09 5.4288352332e+09 0.4000000000 -14.3892566633 -3.4124048303 1 + 1.0000000000e+09 6.0106607628e+09 0.4000000000 -13.9696329645 -3.8570452607 1 + 1.0000000000e+09 6.6548423840e+09 0.4000000000 -13.4990899024 -4.3569659385 1 + 1.0000000000e+09 7.3680629973e+09 0.4000000000 -12.9751097903 -4.9165783433 1 + 1.0000000000e+09 8.1577217310e+09 0.4000000000 -12.4029718543 -5.5362858494 1 + 1.0000000000e+09 9.0320107014e+09 0.4000000000 -11.8075831340 -6.2097753499 1 + 1.0000000000e+09 1.0000000000e+10 0.4000000000 -11.2833005757 -6.9502991890 1 + 1.0000000000e+10 4.0000000000e+09 0.4000000000 -16.2582552424 -1.5626567049 0 + 1.0000000000e+10 4.4286927156e+09 0.4000000000 -16.0486695196 -1.7884208146 1 + 1.0000000000e+10 4.9033297922e+09 0.4000000000 -15.8108280439 -2.0446360339 1 + 1.0000000000e+10 5.4288352332e+09 0.4000000000 -15.5412395159 -2.3352407251 1 + 1.0000000000e+10 6.0106607628e+09 0.4000000000 -15.2357440912 -2.6646601015 1 + 1.0000000000e+10 6.6548423840e+09 0.4000000000 -14.8903553070 -3.0378300279 1 + 1.0000000000e+10 7.3680629973e+09 0.4000000000 -14.5003230284 -3.4601447158 1 + 1.0000000000e+10 8.1577217310e+09 0.4000000000 -14.0614573455 -3.9371815118 1 + 1.0000000000e+10 9.0320107014e+09 0.4000000000 -13.5711738956 -4.4739423061 1 + 1.0000000000e+10 1.0000000000e+10 0.4000000000 -13.0312228621 -5.0733622012 1 + 1.0000000000e+06 4.0000000000e+09 0.4500000000 -12.6199619882 -5.0493824492 0 + 1.0000000000e+06 4.4286927156e+09 0.4500000000 -12.0411567544 -5.6445833881 1 + 1.0000000000e+06 4.9033297922e+09 0.4500000000 -11.4232232546 -6.2868541638 1 + 1.0000000000e+06 5.4288352332e+09 0.4500000000 -10.8518135789 -6.9298795402 1 + 1.0000000000e+06 6.0106607628e+09 0.4500000000 -10.4348297049 -7.7172920610 1 + 1.0000000000e+06 6.6548423840e+09 0.4500000000 -10.0424249574 -8.5877937253 1 + 1.0000000000e+06 7.3680629973e+09 0.4500000000 -9.8849653906 -9.3214540023 1 + 1.0000000000e+06 8.1577217310e+09 0.4500000000 -10.1350525691 -9.9053582092 1 + 1.0000000000e+06 9.0320107014e+09 0.4500000000 -10.9317596586 -10.7720005416 1 + 1.0000000000e+06 1.0000000000e+10 0.4500000000 -12.2046515153 -12.0333675716 1 + 1.0000000000e+07 4.0000000000e+09 0.4500000000 -13.4182964630 -4.2562147910 0 + 1.0000000000e+07 4.4286927156e+09 0.4500000000 -12.9196143636 -4.7702514625 1 + 1.0000000000e+07 4.9033297922e+09 0.4500000000 -12.3646713416 -5.3431229151 1 + 1.0000000000e+07 5.4288352332e+09 0.4500000000 -11.7607991622 -5.9706925457 1 + 1.0000000000e+07 6.0106607628e+09 0.4500000000 -11.1530047519 -6.6222738557 1 + 1.0000000000e+07 6.6548423840e+09 0.4500000000 -10.6534176006 -7.3144511945 1 + 1.0000000000e+07 7.3680629973e+09 0.4500000000 -10.2538965760 -8.1710984842 1 + 1.0000000000e+07 8.1577217310e+09 0.4500000000 -9.9754563098 -8.9846245478 1 + 1.0000000000e+07 9.0320107014e+09 0.4500000000 -10.0361948008 -9.6116151464 1 + 1.0000000000e+07 1.0000000000e+10 0.4500000000 -10.4770374370 -10.2529619738 1 + 1.0000000000e+08 4.0000000000e+09 0.4500000000 -14.2257599781 -3.4625787929 0 + 1.0000000000e+08 4.4286927156e+09 0.4500000000 -13.8113629028 -3.8919431463 1 + 1.0000000000e+08 4.9033297922e+09 0.4500000000 -13.3469606757 -4.3733520313 1 + 1.0000000000e+08 5.4288352332e+09 0.4500000000 -12.8281609635 -4.9118956632 1 + 1.0000000000e+08 6.0106607628e+09 0.4500000000 -12.2549868876 -5.5090642377 1 + 1.0000000000e+08 6.6548423840e+09 0.4500000000 -11.6452720540 -6.1532052811 1 + 1.0000000000e+08 7.3680629973e+09 0.4500000000 -11.0519682942 -6.8137717934 1 + 1.0000000000e+08 8.1577217310e+09 0.4500000000 -10.5745250527 -7.5537012214 1 + 1.0000000000e+08 9.0320107014e+09 0.4500000000 -10.2098734359 -8.4106723339 1 + 1.0000000000e+08 1.0000000000e+10 0.4500000000 -10.0528603240 -9.1577967123 1 + 1.0000000000e+09 4.0000000000e+09 0.4500000000 -15.0536228397 -2.6688948470 0 + 1.0000000000e+09 4.4286927156e+09 0.4500000000 -14.7237134075 -3.0132309281 1 + 1.0000000000e+09 4.9033297922e+09 0.4500000000 -14.3526160907 -3.4007197531 1 + 1.0000000000e+09 5.4288352332e+09 0.4500000000 -13.9355878021 -3.8365577153 1 + 1.0000000000e+09 6.0106607628e+09 0.4500000000 -13.4675884154 -4.3261358021 1 + 1.0000000000e+09 6.6548423840e+09 0.4500000000 -12.9456631769 -4.8738173937 1 + 1.0000000000e+09 7.3680629973e+09 0.4500000000 -12.3724824027 -5.4792910380 1 + 1.0000000000e+09 8.1577217310e+09 0.4500000000 -11.7661786775 -6.1311738817 1 + 1.0000000000e+09 9.0320107014e+09 0.4500000000 -11.1669148875 -6.8085299247 1 + 1.0000000000e+09 1.0000000000e+10 0.4500000000 -10.6514166226 -7.5351746507 1 + 1.0000000000e+10 4.0000000000e+09 0.4500000000 -15.9281437071 -1.8752022471 0 + 1.0000000000e+10 4.4286927156e+09 0.4500000000 -15.6824809806 -2.1344670468 1 + 1.0000000000e+10 4.9033297922e+09 0.4500000000 -15.4047001663 -2.4277705612 1 + 1.0000000000e+10 5.4288352332e+09 0.4500000000 -15.0909187931 -2.7594143654 1 + 1.0000000000e+10 6.0106607628e+09 0.4500000000 -14.7366089802 -3.1341695123 1 + 1.0000000000e+10 6.6548423840e+09 0.4500000000 -14.3375645420 -3.5571892970 1 + 1.0000000000e+10 7.3680629973e+09 0.4500000000 -13.8893254086 -4.0336353982 1 + 1.0000000000e+10 8.1577217310e+09 0.4500000000 -13.3893298655 -4.5676574853 1 + 1.0000000000e+10 9.0320107014e+09 0.4500000000 -12.8393736234 -5.1603429531 1 + 1.0000000000e+10 1.0000000000e+10 0.4500000000 -12.2492265137 -5.8073288375 1 + 1.0000000000e+06 4.0000000000e+09 0.5000000000 -8.0664194354 -9.4890393815 0 + 1.0000000000e+06 4.4286927156e+09 0.5000000000 -8.2111046670 -9.3757376097 1 + 1.0000000000e+06 4.9033297922e+09 0.5000000000 -8.3918234083 -9.2344738095 1 + 1.0000000000e+06 5.4288352332e+09 0.5000000000 -8.7772909062 -8.9448205348 1 + 1.0000000000e+06 6.0106607628e+09 0.5000000000 -9.0636971271 -9.0647301302 1 + 1.0000000000e+06 6.6548423840e+09 0.5000000000 -9.3058452726 -9.3069927949 1 + 1.0000000000e+06 7.3680629973e+09 0.5000000000 -9.5983657200 -9.5996435483 1 + 1.0000000000e+06 8.1577217310e+09 0.5000000000 -10.0168321560 -10.0182538936 1 + 1.0000000000e+06 9.0320107014e+09 0.5000000000 -10.8474057609 -10.8489864867 1 + 1.0000000000e+06 1.0000000000e+10 0.5000000000 -12.1138078856 -12.1155643114 1 + 1.0000000000e+07 4.0000000000e+09 0.5000000000 -7.8838825215 -9.6537090536 0 + 1.0000000000e+07 4.4286927156e+09 0.5000000000 -8.0018897889 -9.5630616368 1 + 1.0000000000e+07 4.9033297922e+09 0.5000000000 -8.1353862644 -9.4610110058 1 + 1.0000000000e+07 5.4288352332e+09 0.5000000000 -8.2961790243 -9.3385091392 1 + 1.0000000000e+07 6.0106607628e+09 0.5000000000 -8.5443481247 -9.1484241687 1 + 1.0000000000e+07 6.6548423840e+09 0.5000000000 -8.9669746383 -8.9715401000 1 + 1.0000000000e+07 7.3680629973e+09 0.5000000000 -9.2004551976 -9.2016374220 1 + 1.0000000000e+07 8.1577217310e+09 0.5000000000 -9.4726511197 -9.4739817205 1 + 1.0000000000e+07 9.0320107014e+09 0.5000000000 -9.8195004568 -9.8209945474 1 + 1.0000000000e+07 1.0000000000e+10 0.5000000000 -10.3608512486 -10.3625253205 1 + 1.0000000000e+08 4.0000000000e+09 0.5000000000 -7.7564783284 -9.7786310606 0 + 1.0000000000e+08 4.4286927156e+09 0.5000000000 -7.8492582957 -9.7096153678 1 + 1.0000000000e+08 4.9033297922e+09 0.5000000000 -7.9528574806 -9.6329691536 1 + 1.0000000000e+08 5.4288352332e+09 0.5000000000 -8.0693126926 -9.5475431816 1 + 1.0000000000e+08 6.0106607628e+09 0.5000000000 -8.2064635264 -9.4472744455 1 + 1.0000000000e+08 6.6548423840e+09 0.5000000000 -8.3859624449 -9.3159632798 1 + 1.0000000000e+08 7.3680629973e+09 0.5000000000 -8.6892137197 -9.0970157579 1 + 1.0000000000e+08 8.1577217310e+09 0.5000000000 -9.0475146736 -9.0502773705 1 + 1.0000000000e+08 9.0320107014e+09 0.5000000000 -9.2986738307 -9.2998953589 1 + 1.0000000000e+08 1.0000000000e+10 0.5000000000 -9.5981915843 -9.5996058564 1 + 1.0000000000e+09 4.0000000000e+09 0.5000000000 -7.7470944343 -9.8201166535 0 + 1.0000000000e+09 4.4286927156e+09 0.5000000000 -7.8148151276 -9.7724887691 1 + 1.0000000000e+09 4.9033297922e+09 0.5000000000 -7.8900133105 -9.7200543142 1 + 1.0000000000e+09 5.4288352332e+09 0.5000000000 -7.9730314841 -9.6630087712 1 + 1.0000000000e+09 6.0106607628e+09 0.5000000000 -8.0660308619 -9.5998690277 1 + 1.0000000000e+09 6.6548423840e+09 0.5000000000 -8.1711502365 -9.5299746078 1 + 1.0000000000e+09 7.3680629973e+09 0.5000000000 -8.2981193470 -9.4461659460 1 + 1.0000000000e+09 8.1577217310e+09 0.5000000000 -8.4688942226 -9.3329031919 1 + 1.0000000000e+09 9.0320107014e+09 0.5000000000 -8.7322319942 -9.1628388393 1 + 1.0000000000e+09 1.0000000000e+10 0.5000000000 -9.0635018458 -9.0827468205 1 + 1.0000000000e+10 4.0000000000e+09 0.5000000000 -7.9880361426 -9.6879018121 1 + 1.0000000000e+10 4.4286927156e+09 0.5000000000 -8.0302837106 -9.6618767633 1 + 1.0000000000e+10 4.9033297922e+09 0.5000000000 -8.0763904321 -9.6342068332 1 + 1.0000000000e+10 5.4288352332e+09 0.5000000000 -8.1258757774 -9.6057820019 1 + 1.0000000000e+10 6.0106607628e+09 0.5000000000 -8.1795463024 -9.5762054669 1 + 1.0000000000e+10 6.6548423840e+09 0.5000000000 -8.2365564026 -9.5471909647 1 + 1.0000000000e+10 7.3680629973e+09 0.5000000000 -8.2989177606 -9.5176692137 1 + 1.0000000000e+10 8.1577217310e+09 0.5000000000 -8.3696249312 -9.4863119797 1 + 1.0000000000e+10 9.0320107014e+09 0.5000000000 -8.4560400376 -9.4486905350 1 + 1.0000000000e+10 1.0000000000e+10 0.5000000000 -8.5725042842 -9.3959527205 1 + 1.0000000000e+06 4.0000000000e+09 0.5500000000 -4.9263369035 -12.6307046567 0 + 1.0000000000e+06 4.4286927156e+09 0.5500000000 -5.5125794078 -12.0760158044 1 + 1.0000000000e+06 4.9033297922e+09 0.5500000000 -6.1678880195 -11.4605459733 1 + 1.0000000000e+06 5.4288352332e+09 0.5500000000 -6.9002340024 -10.8361713231 1 + 1.0000000000e+06 6.0106607628e+09 0.5500000000 -7.7162632259 -10.4358694538 1 + 1.0000000000e+06 6.6548423840e+09 0.5500000000 -8.5866490983 -10.0435826768 1 + 1.0000000000e+06 7.3680629973e+09 0.5500000000 -9.3201805492 -9.8862508729 1 + 1.0000000000e+06 8.1577217310e+09 0.5500000000 -9.9039408495 -10.1364796930 1 + 1.0000000000e+06 9.0320107014e+09 0.5500000000 -10.7704221565 -10.9333429140 1 + 1.0000000000e+06 1.0000000000e+10 0.5500000000 -12.0316130148 -12.2064099076 1 + 1.0000000000e+07 4.0000000000e+09 0.5500000000 -4.1327873451 -13.4067898896 0 + 1.0000000000e+07 4.4286927156e+09 0.5500000000 -4.6339615895 -12.9331082034 1 + 1.0000000000e+07 4.9033297922e+09 0.5500000000 -5.1950956687 -12.4035746041 1 + 1.0000000000e+07 5.4288352332e+09 0.5500000000 -5.8232748773 -11.8138715673 1 + 1.0000000000e+07 6.0106607628e+09 0.5500000000 -6.5262585311 -11.1696968936 1 + 1.0000000000e+07 6.6548423840e+09 0.5500000000 -7.3117723591 -10.6534502323 1 + 1.0000000000e+07 7.3680629973e+09 0.5500000000 -8.1699241087 -10.2551099463 1 + 1.0000000000e+07 8.1577217310e+09 0.5500000000 -8.9833050294 -9.9768125421 1 + 1.0000000000e+07 9.0320107014e+09 0.5500000000 -9.6101346160 -10.0377081495 1 + 1.0000000000e+07 1.0000000000e+10 0.5500000000 -10.2512990136 -10.4787244918 1 + 1.0000000000e+08 4.0000000000e+09 0.5500000000 -3.3395295007 -14.1985085219 0 + 1.0000000000e+08 4.4286927156e+09 0.5500000000 -3.7556211435 -13.8063032478 1 + 1.0000000000e+08 4.9033297922e+09 0.5500000000 -4.2225575794 -13.3664559653 1 + 1.0000000000e+08 5.4288352332e+09 0.5500000000 -4.7464569494 -12.8737301627 1 + 1.0000000000e+08 6.0106607628e+09 0.5500000000 -5.3341586613 -12.3230169469 1 + 1.0000000000e+08 6.6548423840e+09 0.5500000000 -5.9932141535 -11.7120144662 1 + 1.0000000000e+08 7.3680629973e+09 0.5500000000 -6.7312739156 -11.0587794004 1 + 1.0000000000e+08 8.1577217310e+09 0.5500000000 -7.5517647875 -10.5751303759 1 + 1.0000000000e+08 9.0320107014e+09 0.5500000000 -8.4094768040 -10.2111789283 1 + 1.0000000000e+08 1.0000000000e+10 0.5500000000 -9.1564178754 -10.0543419536 1 + 1.0000000000e+09 4.0000000000e+09 0.5500000000 -2.5471797995 -15.0250730353 0 + 1.0000000000e+09 4.4286927156e+09 0.5500000000 -2.8781461534 -14.7143116038 1 + 1.0000000000e+09 4.9033297922e+09 0.5500000000 -3.2508437389 -14.3645024308 1 + 1.0000000000e+09 5.4288352332e+09 0.5500000000 -3.6704138417 -13.9710418050 1 + 1.0000000000e+09 6.0106607628e+09 0.5500000000 -4.1426204134 -13.5288354140 1 + 1.0000000000e+09 6.6548423840e+09 0.5500000000 -4.6739224526 -13.0328615008 1 + 1.0000000000e+09 7.3680629973e+09 0.5500000000 -5.2715267295 -12.4783156933 1 + 1.0000000000e+09 8.1577217310e+09 0.5500000000 -5.9431603152 -11.8633393484 1 + 1.0000000000e+09 9.0320107014e+09 0.5500000000 -6.6951903623 -11.2024350566 1 + 1.0000000000e+09 1.0000000000e+10 0.5500000000 -7.5250446792 -10.6480296449 1 + 1.0000000000e+10 4.0000000000e+09 0.5500000000 -1.7576144396 -15.9280217909 1 + 1.0000000000e+10 4.4286927156e+09 0.5500000000 -2.0033314082 -15.6986279517 1 + 1.0000000000e+10 4.9033297922e+09 0.5500000000 -2.2816705223 -15.4388374567 1 + 1.0000000000e+10 5.4288352332e+09 0.5500000000 -2.5967959027 -15.1448960799 1 + 1.0000000000e+10 6.0106607628e+09 0.5500000000 -2.9533897196 -14.8125306367 1 + 1.0000000000e+10 6.6548423840e+09 0.5500000000 -3.3567153951 -14.4373412142 1 + 1.0000000000e+10 7.3680629973e+09 0.5500000000 -3.8126881350 -14.0143346733 1 + 1.0000000000e+10 8.1577217310e+09 0.5500000000 -4.3279498369 -13.5384648003 1 + 1.0000000000e+10 9.0320107014e+09 0.5500000000 -4.9099077591 -13.0050412453 1 + 1.0000000000e+10 1.0000000000e+10 0.5500000000 -5.5664079579 -12.4111302237 1 + 1.0000000000e+06 4.0000000000e+09 0.6000000000 -4.6874165138 -12.8713534471 0 + 1.0000000000e+06 4.4286927156e+09 0.6000000000 -5.2480529163 -12.3425069636 1 + 1.0000000000e+06 4.9033297922e+09 0.6000000000 -5.8750189170 -11.7561777822 1 + 1.0000000000e+06 5.4288352332e+09 0.6000000000 -6.5760720464 -11.1834054668 1 + 1.0000000000e+06 6.0106607628e+09 0.6000000000 -7.3593255585 -10.8231744084 1 + 1.0000000000e+06 6.6548423840e+09 0.6000000000 -8.2240941596 -10.4375089980 1 + 1.0000000000e+06 7.3680629973e+09 0.6000000000 -9.0836917067 -10.1457704344 1 + 1.0000000000e+06 8.1577217310e+09 0.6000000000 -9.7969116846 -10.2597064355 1 + 1.0000000000e+06 9.0320107014e+09 0.6000000000 -10.6995852342 -11.0267695004 1 + 1.0000000000e+06 1.0000000000e+10 0.6000000000 -11.9565611188 -12.3079449877 1 + 1.0000000000e+07 4.0000000000e+09 0.6000000000 -3.8938729829 -13.6478070662 0 + 1.0000000000e+07 4.4286927156e+09 0.6000000000 -4.3694406419 -13.1998899866 1 + 1.0000000000e+07 4.9033297922e+09 0.6000000000 -4.9022241247 -12.6989326129 1 + 1.0000000000e+07 5.4288352332e+09 0.6000000000 -5.4990224481 -12.1411327078 1 + 1.0000000000e+07 6.0106607628e+09 0.6000000000 -6.1674007417 -11.5348079253 1 + 1.0000000000e+07 6.6548423840e+09 0.6000000000 -6.9155791417 -11.0832126364 1 + 1.0000000000e+07 7.3680629973e+09 0.6000000000 -7.7481347730 -10.7130540280 1 + 1.0000000000e+07 8.1577217310e+09 0.6000000000 -8.6298769859 -10.3622929168 1 + 1.0000000000e+07 9.0320107014e+09 0.6000000000 -9.4167978029 -10.2526624488 1 + 1.0000000000e+07 1.0000000000e+10 0.6000000000 -10.1478611795 -10.6024621642 1 + 1.0000000000e+08 4.0000000000e+09 0.6000000000 -3.1006340152 -14.4403947671 0 + 1.0000000000e+08 4.4286927156e+09 0.6000000000 -3.4911181626 -14.0739382686 1 + 1.0000000000e+08 4.9033297922e+09 0.6000000000 -3.9297028770 -13.6626066662 1 + 1.0000000000e+08 5.4288352332e+09 0.6000000000 -4.4222123883 -13.2014814753 1 + 1.0000000000e+08 6.0106607628e+09 0.6000000000 -4.9751654635 -12.6858717662 1 + 1.0000000000e+08 6.6548423840e+09 0.6000000000 -5.5958360589 -12.1142361002 1 + 1.0000000000e+08 7.3680629973e+09 0.6000000000 -6.2921951408 -11.5102059086 1 + 1.0000000000e+08 8.1577217310e+09 0.6000000000 -7.0718360681 -11.0958538127 1 + 1.0000000000e+08 9.0320107014e+09 0.6000000000 -7.9295315244 -10.7329004795 1 + 1.0000000000e+08 1.0000000000e+10 0.6000000000 -8.7971964340 -10.4475908920 1 + 1.0000000000e+09 4.0000000000e+09 0.6000000000 -2.3083427797 -15.2688913657 0 + 1.0000000000e+09 4.4286927156e+09 0.6000000000 -2.6136989393 -14.9838710003 1 + 1.0000000000e+09 4.9033297922e+09 0.6000000000 -2.9580422061 -14.6625625008 1 + 1.0000000000e+09 5.4288352332e+09 0.6000000000 -3.3462197500 -14.3006599290 1 + 1.0000000000e+09 6.0106607628e+09 0.6000000000 -3.7836699672 -13.8934047622 1 + 1.0000000000e+09 6.6548423840e+09 0.6000000000 -4.2764930786 -13.4361677469 1 + 1.0000000000e+09 7.3680629973e+09 0.6000000000 -4.8315236918 -12.9246340459 1 + 1.0000000000e+09 8.1577217310e+09 0.6000000000 -5.4563463669 -12.3575652127 1 + 1.0000000000e+09 9.0320107014e+09 0.6000000000 -6.1589448250 -11.7514500907 1 + 1.0000000000e+09 1.0000000000e+10 0.6000000000 -6.9457107700 -11.2744549299 1 + 1.0000000000e+10 4.0000000000e+09 0.6000000000 -1.5189550687 -16.1760474376 0 + 1.0000000000e+10 4.4286927156e+09 0.6000000000 -1.7390541688 -15.9723914831 1 + 1.0000000000e+10 4.9033297922e+09 0.6000000000 -1.9890315435 -15.7410961693 1 + 1.0000000000e+10 5.4288352332e+09 0.6000000000 -2.2727571929 -15.4787024815 1 + 1.0000000000e+10 6.0106607628e+09 0.6000000000 -2.5945876455 -15.1812641094 1 + 1.0000000000e+10 6.6548423840e+09 0.6000000000 -2.9594257944 -14.8447444844 1 + 1.0000000000e+10 7.3680629973e+09 0.6000000000 -3.3727879321 -14.4645568357 1 + 1.0000000000e+10 8.1577217310e+09 0.6000000000 -3.8408780345 -14.0361164439 1 + 1.0000000000e+10 9.0320107014e+09 0.6000000000 -4.3706608771 -13.5552459151 1 + 1.0000000000e+10 1.0000000000e+10 0.6000000000 -4.9698583638 -13.0193464360 1 + 1.0000000000e+06 4.0000000000e+09 0.6500000000 -4.5476580932 -13.0130256076 0 + 1.0000000000e+06 4.4286927156e+09 0.6500000000 -5.0933157554 -12.4994381879 1 + 1.0000000000e+06 4.9033297922e+09 0.6500000000 -5.7036990969 -11.9307186936 1 + 1.0000000000e+06 5.4288352332e+09 0.6500000000 -6.3864064622 -11.4022027916 1 + 1.0000000000e+06 6.0106607628e+09 0.6500000000 -7.1496914530 -11.0673499344 1 + 1.0000000000e+06 6.6548423840e+09 0.6500000000 -7.9986226626 -10.7004084995 1 + 1.0000000000e+06 7.3680629973e+09 0.6500000000 -8.8932097425 -10.3705250415 1 + 1.0000000000e+06 8.1577217310e+09 0.6500000000 -9.6962265139 -10.3881460503 1 + 1.0000000000e+06 9.0320107014e+09 0.6500000000 -10.6339946188 -11.1316988308 1 + 1.0000000000e+06 1.0000000000e+10 0.6500000000 -11.8875095425 -12.4230319679 1 + 1.0000000000e+07 4.0000000000e+09 0.6500000000 -3.7541203442 -13.7898249789 0 + 1.0000000000e+07 4.4286927156e+09 0.6500000000 -4.2147089469 -13.3570697965 1 + 1.0000000000e+07 4.9033297922e+09 0.6500000000 -4.7309081311 -12.8729632474 1 + 1.0000000000e+07 5.4288352332e+09 0.6500000000 -5.3093462612 -12.3341734101 1 + 1.0000000000e+07 6.0106607628e+09 0.6500000000 -5.9574209643 -11.7528786941 1 + 1.0000000000e+07 6.6548423840e+09 0.6500000000 -6.6832983524 -11.3537044311 1 + 1.0000000000e+07 7.3680629973e+09 0.6500000000 -7.4941774420 -11.0089926944 1 + 1.0000000000e+07 8.1577217310e+09 0.6500000000 -8.3792738800 -10.6556306767 1 + 1.0000000000e+07 9.0320107014e+09 0.6500000000 -9.2440179618 -10.4606562951 1 + 1.0000000000e+07 1.0000000000e+10 0.6500000000 -10.0502406180 -10.7352829937 1 + 1.0000000000e+08 4.0000000000e+09 0.6500000000 -2.9608994723 -14.5832351693 0 + 1.0000000000e+08 4.4286927156e+09 0.6500000000 -3.3364036964 -14.2319225764 1 + 1.0000000000e+08 4.9033297922e+09 0.6500000000 -3.7584032382 -13.8373671291 1 + 1.0000000000e+08 5.4288352332e+09 0.6500000000 -4.2325503160 -13.3948530013 1 + 1.0000000000e+08 6.0106607628e+09 0.6500000000 -4.7651740951 -12.9000045811 1 + 1.0000000000e+08 6.6548423840e+09 0.6500000000 -5.3633523384 -12.3521350793 1 + 1.0000000000e+08 7.3680629973e+09 0.6500000000 -6.0349472150 -11.7846818528 1 + 1.0000000000e+08 8.1577217310e+09 0.6500000000 -6.7881653086 -11.4262642633 1 + 1.0000000000e+08 9.0320107014e+09 0.6500000000 -7.6263287309 -11.0865976478 1 + 1.0000000000e+08 1.0000000000e+10 0.6500000000 -8.5198307182 -10.7737584967 1 + 1.0000000000e+09 4.0000000000e+09 0.6500000000 -2.1686642542 -15.4135601959 0 + 1.0000000000e+09 4.4286927156e+09 0.6500000000 -2.4590379093 -15.1436760063 1 + 1.0000000000e+09 4.9033297922e+09 0.6500000000 -2.7867935242 -14.8391271318 1 + 1.0000000000e+09 5.4288352332e+09 0.6500000000 -3.1566062187 -14.4957844807 1 + 1.0000000000e+09 6.0106607628e+09 0.6500000000 -3.5737239146 -14.1090934401 1 + 1.0000000000e+09 6.6548423840e+09 0.6500000000 -4.0440362374 -13.6746794679 1 + 1.0000000000e+09 7.3680629973e+09 0.6500000000 -4.5741490489 -13.1886122154 1 + 1.0000000000e+09 8.1577217310e+09 0.6500000000 -5.1714414325 -12.6505348485 1 + 1.0000000000e+09 9.0320107014e+09 0.6500000000 -5.8439871511 -12.0827672096 1 + 1.0000000000e+09 1.0000000000e+10 0.6500000000 -6.5997012945 -11.6763905107 1 + 1.0000000000e+10 4.0000000000e+09 0.6500000000 -1.3794466033 -16.3246961572 0 + 1.0000000000e+10 4.4286927156e+09 0.6500000000 -1.5845558746 -16.1361734166 1 + 1.0000000000e+10 4.9033297922e+09 0.6500000000 -1.8179385128 -15.9216322620 1 + 1.0000000000e+10 5.4288352332e+09 0.6500000000 -2.0832924664 -15.6777872875 1 + 1.0000000000e+10 6.0106607628e+09 0.6500000000 -2.3847837732 -15.4008846806 1 + 1.0000000000e+10 6.6548423840e+09 0.6500000000 -2.7271044222 -15.0871030464 1 + 1.0000000000e+10 7.3680629973e+09 0.6500000000 -3.1155372402 -14.7321000031 1 + 1.0000000000e+10 8.1577217310e+09 0.6500000000 -3.5560283588 -14.3315813824 1 + 1.0000000000e+10 9.0320107014e+09 0.6500000000 -4.0552644947 -13.8817455175 1 + 1.0000000000e+10 1.0000000000e+10 0.6500000000 -4.6207246098 -13.3805467847 1 + 1.0000000000e+06 4.0000000000e+09 0.7000000000 -4.4484985951 -13.1143475495 0 + 1.0000000000e+06 4.4286927156e+09 0.7000000000 -4.9835286255 -12.6117296860 1 + 1.0000000000e+06 4.9033297922e+09 0.7000000000 -5.5821458226 -12.0561621904 1 + 1.0000000000e+06 5.4288352332e+09 0.7000000000 -6.2518304171 -11.5720026737 1 + 1.0000000000e+06 6.0106607628e+09 0.7000000000 -7.0008152024 -11.2561319512 1 + 1.0000000000e+06 6.6548423840e+09 0.7000000000 -7.8361158291 -10.9067031322 1 + 1.0000000000e+06 7.3680629973e+09 0.7000000000 -8.7392063208 -10.5687519073 1 + 1.0000000000e+06 8.1577217310e+09 0.7000000000 -9.6021715722 -10.5229664649 1 + 1.0000000000e+06 9.0320107014e+09 0.7000000000 -10.5729780670 -11.2517541081 1 + 1.0000000000e+06 1.0000000000e+10 0.7000000000 -11.8235708727 -12.5558647155 1 + 1.0000000000e+07 4.0000000000e+09 0.7000000000 -3.6549664086 -13.8914709006 0 + 1.0000000000e+07 4.4286927156e+09 0.7000000000 -4.1049270968 -13.4695590574 1 + 1.0000000000e+07 4.9033297922e+09 0.7000000000 -4.6093594165 -12.9975417362 1 + 1.0000000000e+07 5.4288352332e+09 0.7000000000 -5.1747701322 -12.4726172619 1 + 1.0000000000e+07 6.0106607628e+09 0.7000000000 -5.8084293227 -11.9132714976 1 + 1.0000000000e+07 6.6548423840e+09 0.7000000000 -6.5184076951 -11.5627682932 1 + 1.0000000000e+07 7.3680629973e+09 0.7000000000 -7.3127256426 -11.2391950016 1 + 1.0000000000e+07 8.1577217310e+09 0.7000000000 -8.1900933044 -10.8968169016 1 + 1.0000000000e+07 9.0320107014e+09 0.7000000000 -9.0920365426 -10.6614850215 1 + 1.0000000000e+07 1.0000000000e+10 0.7000000000 -9.9582648682 -10.8798139159 1 + 1.0000000000e+08 4.0000000000e+09 0.7000000000 -2.8617629403 -14.6856624315 0 + 1.0000000000e+08 4.4286927156e+09 0.7000000000 -3.2266384172 -14.3451718155 1 + 1.0000000000e+08 4.9033297922e+09 0.7000000000 -3.6368702837 -13.9626100890 1 + 1.0000000000e+08 5.4288352332e+09 0.7000000000 -4.0979886993 -13.5334235631 1 + 1.0000000000e+08 6.0106607628e+09 0.7000000000 -4.6161879732 -13.0535392989 1 + 1.0000000000e+08 6.6548423840e+09 0.7000000000 -5.1984011582 -12.5233515125 1 + 1.0000000000e+08 7.3680629973e+09 0.7000000000 -5.8523647926 -11.9934421813 1 + 1.0000000000e+08 8.1577217310e+09 0.7000000000 -6.5864090547 -11.6821331659 1 + 1.0000000000e+08 9.0320107014e+09 0.7000000000 -7.4067994755 -11.3654235744 1 + 1.0000000000e+08 1.0000000000e+10 0.7000000000 -8.3013506853 -11.0538132814 1 + 1.0000000000e+09 4.0000000000e+09 0.7000000000 -2.0695815696 -15.5177250154 0 + 1.0000000000e+09 4.4286927156e+09 0.7000000000 -2.3493240008 -15.2586547889 1 + 1.0000000000e+09 4.9033297922e+09 0.7000000000 -2.6653095613 -14.9660806276 1 + 1.0000000000e+09 5.4288352332e+09 0.7000000000 -3.0220912983 -14.6360018221 1 + 1.0000000000e+09 6.0106607628e+09 0.7000000000 -3.4247819832 -14.2640168398 1 + 1.0000000000e+09 6.6548423840e+09 0.7000000000 -3.8791217034 -13.8459596123 1 + 1.0000000000e+09 7.3680629973e+09 0.7000000000 -4.3915522296 -13.3782729765 1 + 1.0000000000e+09 8.1577217310e+09 0.7000000000 -4.9692869712 -12.8617482228 1 + 1.0000000000e+09 9.0320107014e+09 0.7000000000 -5.6203221431 -12.3289938297 1 + 1.0000000000e+09 1.0000000000e+10 0.7000000000 -6.3529920884 -11.9887960558 1 + 1.0000000000e+10 4.0000000000e+09 0.7000000000 -1.2805272580 -16.4326417920 0 + 1.0000000000e+10 4.4286927156e+09 0.7000000000 -1.4749982900 -16.2549301551 1 + 1.0000000000e+10 4.9033297922e+09 0.7000000000 -1.6966040866 -16.0523580418 1 + 1.0000000000e+10 5.4288352332e+09 0.7000000000 -1.9489205221 -15.8217642481 1 + 1.0000000000e+10 6.0106607628e+09 0.7000000000 -2.2359784782 -15.5595330834 1 + 1.0000000000e+10 6.6548423840e+09 0.7000000000 -2.5623203021 -15.2619989382 1 + 1.0000000000e+10 7.3680629973e+09 0.7000000000 -2.9330631536 -14.9250007236 1 + 1.0000000000e+10 8.1577217310e+09 0.7000000000 -3.3539699474 -14.5444727373 1 + 1.0000000000e+10 9.0320107014e+09 0.7000000000 -3.8315269899 -14.1169483246 1 + 1.0000000000e+10 1.0000000000e+10 0.7000000000 -4.3730130767 -13.6410243937 1 +AMReX (23.07-445-gc09c45c31007) finalized diff --git a/util/approx_math/approx_math.H b/util/approx_math/approx_math.H index 9f4a6a4835..97fa6b43b3 100644 --- a/util/approx_math/approx_math.H +++ b/util/approx_math/approx_math.H @@ -3,6 +3,7 @@ #include #include +#include #include #include #include @@ -38,6 +39,19 @@ amrex::Real fast_atan_1(const amrex::Real x) { constexpr amrex::Real A = 0.0776509570923569_rt; constexpr amrex::Real B = -0.287434475393028_rt; constexpr amrex::Real C = GCEM_PI / 4.0_rt - A - B; + constexpr amrex::Real third = 1.0_rt / 3.0_rt; + + /// + /// If x < 0.378, then using arctan(x) ~ x - x^3/3 + x^5/5 + /// gives better answer than the approximation below. + /// And accuracy increase as x << 0.378. + /// + + if (std::abs(x) < 0.378_rt) { + // return -0.3333333333_rt*amrex::Math::powi<3>(x) + x; + amrex::Real x2 = x*x; + return ((0.2_rt*x2 - third)*x2 + 1.0_rt)*x; + } amrex::Real x2 = x*x; return ((A*x2 + B)*x2 + C)*x; @@ -52,17 +66,6 @@ amrex::Real fast_atan(const amrex::Real x) { constexpr amrex::Real PI_2 = 0.5_rt * GCEM_PI; - /// - /// If x < 0.113, then using arctan(x) ~ x - /// gives better answer than the approximation below. - /// And accuracy increase as x << 0.113. - /// Also significantly faster. - /// - - if (std::abs(x) < 0.113_rt) { - return x; - } - // Check for large number, close to infinity. // Error is ~1e-8 rad by not checking actual inf From 8fd0d4aed4b7934a3a704fe938614685df988305 Mon Sep 17 00:00:00 2001 From: Michael Zingale Date: Mon, 15 Jul 2024 17:13:08 -0400 Subject: [PATCH 5/5] rework the unit test docs (#1570) We split the current page into 3 separate pages and add descriptions for all of the unit tests --- sphinx_docs/source/comprehensive_tests.rst | 119 +++++++++ sphinx_docs/source/index.rst | 11 +- sphinx_docs/source/one_zone_tests.rst | 103 ++++++++ sphinx_docs/source/unit_tests.rst | 280 ++++++++------------- 4 files changed, 332 insertions(+), 181 deletions(-) create mode 100644 sphinx_docs/source/comprehensive_tests.rst create mode 100644 sphinx_docs/source/one_zone_tests.rst diff --git a/sphinx_docs/source/comprehensive_tests.rst b/sphinx_docs/source/comprehensive_tests.rst new file mode 100644 index 0000000000..5ea20ab663 --- /dev/null +++ b/sphinx_docs/source/comprehensive_tests.rst @@ -0,0 +1,119 @@ +************************ +Comprehensive Unit Tests +************************ + +Generally, for each test, you simply type ``make`` in the test +directory. There are a number of runtime parameters that can +control the behavior. These are specified (along with defaults) +in ``_parameters`` files in each test directory and can be +overridden in an inputs file or on the commandline. + +Some additional details on a few of the comprehensive unit tests +are given below. + +EOS test (``test_eos``) +======================= + +``Microphysics/unit_test/test_eos/`` is a unit test for the equations +of state in Microphysics. It sets up a cube of data, with +:math:`\rho`, :math:`T`, and :math:`X_k` varying along the three +dimensions, and then calls the EOS in each zone. Calls are done to +exercise all modes of calling the EOS, in order: + +- ``eos_input_rt``: We call the EOS using :math:`\rho, T`. this is the + reference call, and we save the :math:`h`, :math:`e`, :math:`p`, and + :math:`s` from here to use in subsequent calls. + +- ``eos_input_rh``: We call the EOS using :math:`\rho, h`, to recover + the original :math:`T`. To give the root finder some work to do, we + perturb the initial temperature. + + We store the relative error in :math:`T` in the output file. + +- ``eos_input_tp``: We call the EOS using :math:`T, p`, to recover the + original :math:`\rho`. To give the root finder some work to do, we + perturb the initial density. + + We store the relative error in :math:`\rho` in the output file. + +- ``eos_input_rp``: We call the EOS using :math:`\rho, p`, to recover + the original :math:`T`. To give the root finder some work to do, we + perturb the initial temperature. + + We store the relative error in :math:`T` in the output file. + +- ``eos_input_re``: We call the EOS using :math:`\rho, e`, to recover + the original :math:`T`. To give the root finder some work to do, we + perturb the initial temperature. + + We store the relative error in :math:`T` in the output file. + +- ``eos_input_ps``: We call the EOS using :math:`p, s`, to recover the + original :math:`\rho, T`. To give the root finder some work to do, + we perturb the initial density and temperature. + + Note: entropy is not well-defined for some EOSs, so we only attempt + the root find if :math:`s > 0`. + + We store the relative error in :math:`\rho, T` in the output file. + +- ``eos_input_ph``: We call the EOS using :math:`p, h`, to recover the + original :math:`\rho, T`. To give the root finder some work to do, + we perturb the initial density and temperature. + + We store the relative error in :math:`\rho, T` in the output file. + +- ``eos_input_th``: We call the EOS using :math:`T, h`, to recover the + original :math:`\rho`. To give the root finder some work to do, we + perturb the initial density. + + Note: for some EOSs, :math:`h = h(\rho)` (e.g., an ideal gas), so there + is no temperature dependence, and we do not do this test. + + We store the relative error in :math:`\rho` in the output file. + +This unit test is marked up with OpenMP directives and therefore also +tests whether the EOS is threadsafe. + +To compile for a specific EOS, e.g., helmholtz, do:: + + make EOS_DIR=helmholtz -j 4 + +Examining the output (an AMReX plotfile) will show you how big the +errors are. You can use the ``amrex/Tools/Plotfile/`` tool +``fextrema`` to display the maximum error for each variable. + + +Network test (``test_react``) +============================= + +``Microphysics/unit_test/test_react/`` is a unit test for the nuclear +reaction networks in Microphysics. It sets up a cube of data, with +:math:`\rho`, :math:`T`, and :math:`X_k` varying along the three +dimensions (as a :math:`16^3` domain), and then calls the EOS in each +zone. This test does the entire ODE integration of the network for +each zone. + +The state in each zone of our data cube is determined by the runtime +parameters ``dens_min``, ``dens_max``, ``temp_min``, and ``temp_max`` +for :math:`(\rho, T)`. Because each network carries different +compositions, we specify the composition through runtime parameters in +the ``&extern`` namelist: ``primary_species_1``, +``primary_species_2``, ``primary_species_3``. These primary species +will vary from X = 0.2 to X = 0.7 to 0.9 (depending on the number). +Only one primary species varies at a time. The non-primary species +will be set equally to share whatever fraction of 1 is not accounted +for by the primary species mass fractions. + +This test calls the network on each zone, running for a time +``tmax``. The full state, including new mass fractions and energy +release is output to a AMReX plotfile. + +You can compile for a specific integrator (e.g., ``VODE``) or +network (e.g., ``aprox13``) as:: + + make NETWORK_DIR=aprox13 INTEGRATOR_DIR=VODE -j 4 + +The loop over the burner is marked up for OpenMP and CUDA and +therefore this test can be used to assess threadsafety of the burners +as well as to optimize the GPU performance of the burners. diff --git a/sphinx_docs/source/index.rst b/sphinx_docs/source/index.rst index 5c547fcd43..7b48af3de1 100644 --- a/sphinx_docs/source/index.rst +++ b/sphinx_docs/source/index.rst @@ -24,7 +24,6 @@ for astrophysical simulation codes. data_structures autodiff rp_intro - unit_tests .. toctree:: :maxdepth: 1 @@ -53,7 +52,15 @@ for astrophysical simulation codes. .. toctree:: :maxdepth: 1 - :caption: references + :caption: Unit tests + + unit_tests + comprehensive_tests + one_zone_tests + +.. toctree:: + :maxdepth: 1 + :caption: References zreferences diff --git a/sphinx_docs/source/one_zone_tests.rst b/sphinx_docs/source/one_zone_tests.rst new file mode 100644 index 0000000000..ec9706a05e --- /dev/null +++ b/sphinx_docs/source/one_zone_tests.rst @@ -0,0 +1,103 @@ +************** +One Zone Tests +************** + +There are several tests that let you call the EOS or reaction network +on a single zone to inspect the output directly. + + +``burn_cell`` +============= + +``burn_cell`` is a simple one-zone burn that will evolve a state with +a network for a specified amount of time. This can be used to +understand the timescales involved in a reaction sequence or to +determine the needed ODE tolerances. + + +Getting Started +--------------- + +The ``burn_cell`` code are located in +``Microphysics/unit_test/burn_cell``. To run a simulation, ensure that +both an input file and an initial conditions file have been created +and are in the same directory as the executable. + +Input File +---------- + +These files are typically named as ``inputs_burn_network`` where network +is the network you wish to use for your testing. + +The structure of this file is is fairly self-explanatory. The run +prefix defined should be unique to the tests that will be run as they +will be used to identify all of the output files. Typically, the run +prefix involves the name of the network being tested. The ``atol`` +variables define absolute tolerances of the ordinary differential +equations and the ``rtol`` variables define the relative tolerances. The +second section of the input file collects the inputs that ``main.f90`` +asks for so that the user does not have to input all 5+ +parameters that are required every time the test is run. Each input +required is defined and initialized on the lines following +``&cellparams``. The use of the parameters is show below: + +.. table:: The definition of parameters used in the burn_cell unit tests and specified in the second half of each inputs file. + + +-----------------------+----------------------------------------+ + | ``tmax`` | Maximum Time (s) | + +-----------------------+----------------------------------------+ + | ``nsteps`` | Number of time subdivisions | + +-----------------------+----------------------------------------+ + | ``density`` | State Density (:math:`\frac{g}{cm^3}`) | + +-----------------------+----------------------------------------+ + | ``temperature`` | State Temperature (K) | + +-----------------------+----------------------------------------+ + | ``massfractions(i)`` | Mass Fraction for element i | + +-----------------------+----------------------------------------+ + +Running the Code +---------------- + +To run the code, enter the burn_cell directory and run:: + + ./main3d.gnu.ex inputs + +where ``inputs`` is the name of your inputs file. + +For each of the ``numsteps`` steps defined in the inputs +file, the code will output a files into a new directory titled +``run_prefix_output`` where ``run_prefix`` is the run prefix defined in the +inputs file. Each output file will be named using the run prefix +defined in the inputs file and the corresponding timestep. + +Next, run ``burn_cell.py`` using python 3.x, giving the defined run prefix as an argument. +For example:: + + python3 burn_cell.py react_aprox13 + +The ``burn_cell.py`` code will gather information from all of the +output files and compile them into three graphs explained below. + +Graphs Output by ``burn_cell.py`` +--------------------------------- + +The file ``run-prefix_logX.png`` and ``run-prefix_logX.eps`` will display a +graph of the chemical abundances as a function of the time, both on +logarithmic scales, for all species involved in the simulation. An +example of this graph is shown below. + +.. figure:: react_aprox13_logX.png + :alt: An example of a plot output by the burn_cell unit test. This is the logX output corresponding to the network aprox13. + :width: 4.5in + + An example of a plot output by the burn_cell unit test. This is the + logX output corresponding to the network aprox13. + + + +The file ``run-prefix_ydot.png`` and ``run-prefix_ydot.eps`` will display the +molar fraction (mass fraction / atomic weight) as a function of time, +both on logarithmic scales, for all species involved in the code. + +The file ``run-prefix_T-edot.png`` and ``run-prefix_T-edot.eps`` will display +the temperature and the energy generation rate as a function of time. diff --git a/sphinx_docs/source/unit_tests.rst b/sphinx_docs/source/unit_tests.rst index ab4d58e4b4..d9ba9c5d1e 100644 --- a/sphinx_docs/source/unit_tests.rst +++ b/sphinx_docs/source/unit_tests.rst @@ -1,10 +1,6 @@ -********** -Unit Tests -********** - - -Comprehensive Unit Tests -======================== +********************** +Overview of Unit Tests +********************** There are a few unit tests in Microphysics that operate on a generic EOS, reaction network, conductivity, or some smaller component to @@ -23,229 +19,155 @@ script. Most of these tests work with MPI+OpenMP, MPI+CUDA, and MPI+HIP +Tests are divided into three categories: -EOS test --------- - -``Microphysics/unit_test/test_eos/`` is a unit test for the equations -of state in Microphysics. It sets up a cube of data, with -:math:`\rho`, :math:`T`, and :math:`X_k` varying along the three -dimensions, and then calls the EOS in each zone. Calls are done to -exercise all modes of calling the EOS, in order: - -- ``eos_input_rt``: We call the EOS using :math:`\rho, T`. this is the - reference call, and we save the :math:`h`, :math:`e`, :math:`p`, and - :math:`s` from here to use in subsequent calls. - -- ``eos_input_rh``: We call the EOS using :math:`\rho, h`, to recover - the original :math:`T`. To give the root finder some work to do, we - perturb the initial temperature. - - We store the relative error in :math:`T` in the output file. - -- ``eos_input_tp``: We call the EOS using :math:`T, p`, to recover the - original :math:`\rho`. To give the root finder some work to do, we - perturb the initial density. - - We store the relative error in :math:`\rho` in the output file. - -- ``eos_input_rp``: We call the EOS using :math:`\rho, p`, to recover - the original :math:`T`. To give the root finder some work to do, we - perturb the initial temperature. - - We store the relative error in :math:`T` in the output file. - -- ``eos_input_re``: We call the EOS using :math:`\rho, e`, to recover - the original :math:`T`. To give the root finder some work to do, we - perturb the initial temperature. - - We store the relative error in :math:`T` in the output file. +* *comprehensive tests* work on a cube of data (usually + $\rho$, $T$, and composition varying along the three dimensions) and + are meant to exercise a wide range of input conditions. -- ``eos_input_ps``: We call the EOS using :math:`p, s`, to recover the - original :math:`\rho, T`. To give the root finder some work to do, - we perturb the initial density and temperature. + These are mainly used for regression testing. - Note: entropy is not well-defined for some EOSs, so we only attempt - the root find if :math:`s > 0`. +* *one-zone tests* allow you to evaluate the conditions for a + particular thermodynamic state. - We store the relative error in :math:`\rho, T` in the output file. + These are often used for interactive explorations and within the CI. -- ``eos_input_ph``: We call the EOS using :math:`p, h`, to recover the - original :math:`\rho, T`. To give the root finder some work to do, - we perturb the initial density and temperature. +* *infrastructure tests* test small bits of the solver, function + interfaces, or runtime infrastructure. - We store the relative error in :math:`\rho, T` in the output file. + These are not really meant for exploring the actual thermodynamic + state. -- ``eos_input_th``: We call the EOS using :math:`T, h`, to recover the - original :math:`\rho`. To give the root finder some work to do, we - perturb the initial density. - Note: for some EOSs, :math:`h = h(\rho)` (e.g., an ideal gas), so there - is no temperature dependence, and we do not do this test. - We store the relative error in :math:`\rho` in the output file. +Comprehensive tests +=================== -This unit test is marked up with OpenMP directives and therefore also -tests whether the EOS is threadsafe. +Each of these tests sets up a cube of data, $(\rho, T, X_k)$, with the +range of $T$ and $\rho$, and the species to focus on for $X_k$ controlled +by options in the input file. -To compile for a specific EOS, e.g., helmholtz, do:: +* ``test_aprox_rates`` : - make EOS_DIR=helmholtz -j 4 + call each of the hardcoded rate functions in ``Microphysics/rates/`` + on each cell in the data cube and store the output in a plotfile. -Examining the output (an AMReX plotfile) will show you how big the -errors are. You can use the ``amrex/Tools/Plotfile/`` tool -``fextrema`` to display the maximum error for each variable. +* ``test_conductivity`` : + call one of the conductivity routines (set via ``CONDUCTIVITY_DIR``) + on each cell in the data cube and store the output in a plotfile. -Network test ------------- +* ``test_eos`` : -``Microphysics/unit_test/test_react/`` is a unit test for the nuclear -reaction networks in Microphysics. It sets up a cube of data, with -:math:`\rho`, :math:`T`, and :math:`X_k` varying along the three -dimensions (as a :math:`16^3` domain), and then calls the EOS in each -zone. This test does the entire ODE integration of the network for -each zone. + call one of the equations of state (set via ``EOS_DIR``) on each + cell in the data cube. We first call it with $(\rho, T, X_k)$ as + input (``eos_input_rt``), and then test each of the other EOS modes + (``eos_input_rh``, ``eos_input_tp``, ``eos_input_rp``, + ``eos_input_re``, ``eos_input_ps``, ``eos_input_ph``, + ``eos_input_th``) and for each of these modes, we compute the error + in the recovered $T$ and/or $\rho$ (as appropriate). The full + thermodynamic state and errors are stored in a plotfile. -The state in each zone of our data cube is determined by the runtime -parameters ``dens_min``, ``dens_max``, ``temp_min``, and ``temp_max`` -for :math:`(\rho, T)`. Because each network carries different -compositions, we specify the composition through runtime parameters in -the ``&extern`` namelist: ``primary_species_1``, -``primary_species_2``, ``primary_species_3``. These primary species -will vary from X = 0.2 to X = 0.7 to 0.9 (depending on the number). -Only one primary species varies at a time. The non-primary species -will be set equally to share whatever fraction of 1 is not accounted -for by the primary species mass fractions. +* ``test_jac`` : -This test calls the network on each zone, running for a time -``tmax``. The full state, including new mass fractions and energy -release is output to a AMReX plotfile. + for each cell in the data cube, and for a specific network (set via + ``NETWORK_DIR``) call the analytic Jacobian provided by the network + and compute the numerical Jacobian (via finite differencing) and + store the relative difference between each approximation for each + Jacobian element in a plotfile. -You can compile for a specific integrator (e.g., ``VODE``) or -network (e.g., ``aprox13``) as:: +* ``test_neutrino_cooling`` : - make NETWORK_DIR=aprox13 INTEGRATOR_DIR=VODE -j 4 + for each cell in the data cube, call the neutrino cooling function + and store the output in a plotfile. -The loop over the burner is marked up for OpenMP and CUDA and -therefore this test can be used to assess threadsafety of the burners -as well as to optimize the GPU performance of the burners. +* ``test_react`` : + for each cell in the data cube, call the reaction network and + integrate to a specified time. Statistics about the number of RHS + calls are reported at the end. A lot of options can be set via the + inputs file to control the integration. -Aprox Rates Test ----------------- +* ``test_rhs`` : -``Microphysics/unit_test/test_aprox_rates`` just evaluates the -instantaneous reaction rates in ``Microphysics/rates/`` used by the -``iso7``, ``aprox13``, ``aprox19``, and ``aprox21`` reaction networks. -This uses the same basic ideas as the tests above---a cube of data is -setup and the rates are evaluated using each zone's thermodynamic -conditions. This test is not really network specific---it tests all -of the available rates. + for each cell in the data cube, call the network's righthand side and + Jacobian functions and store the output in a plotfile. The network + is controlled by the ``NETWORK_DIR`` variable. +* ``test_screening_templated`` : -Screening Test --------------- + for any of the templated-C++ networks, this test will loop over all of + the rates and calculate the screening factor (the screening routine can + be set via ``SCREEN_METHOD``). The factors for each cell in the data + cube are written to a plotfile. -``Microphysics/unit_test/test_screening`` just evaluates the screening -routine, using the ``aprox21`` reaction network. -This uses the same basic ideas as the tests above---a cube of data is -setup and the rates are evaluated using each zone's thermodynamic -conditions. +* ``test_sdc`` : + similar to ``test_react``, except now we exercise the SDC + integration infrastructure. The conserved state that is input to + the burner is chosen to have a Mach number of $0.1$ (and otherwise + use the thermodynamics from the data cube). No advective terms are + modeled. -``burn_cell`` -============= -``burn_cell`` is a simple one-zone burn that will evolve a state with -a network for a specified amount of time. This can be used to -understand the timescales involved in a reaction sequence or to -determine the needed ODE tolerances. +One-zone tests +============== +* ``burn_cell`` : -Getting Started ---------------- + given a $\rho$, $T$, and $X_k$, integrate a reaction network through a specified time + and output the new state. -The ``burn_cell`` code are located in -``Microphysics/unit_test/burn_cell``. To run a simulation, ensure that -both an input file and an initial conditions file have been created -and are in the same directory as the executable. +* ``burn_cell_primordial_chem`` : -Input File ----------- + similar to ``burn_cell`` except specific for the primordial chemistry network. -These files are typically named as ``inputs_burn_network`` where network -is the network you wish to use for your testing. +* ``burn_cell_sdc`` : -The structure of this file is is fairly self-explanatory. The run -prefix defined should be unique to the tests that will be run as they -will be used to identify all of the output files. Typically, the run -prefix involves the name of the network being tested. The ``atol`` -variables define absolute tolerances of the ordinary differential -equations and the ``rtol`` variables define the relative tolerances. The -second section of the input file collects the inputs that ``main.f90`` -asks for so that the user does not have to input all 5+ -parameters that are required every time the test is run. Each input -required is defined and initialized on the lines following -``&cellparams``. The use of the parameters is show below: + similar to ``burn_cell`` except this uses the SDC integration code paths. -.. table:: The definition of parameters used in the burn_cell unit tests and specified in the second half of each inputs file. +* ``eos_cell`` : - +-----------------------+----------------------------------------+ - | ``tmax`` | Maximum Time (s) | - +-----------------------+----------------------------------------+ - | ``numsteps`` | Number of time subdivisions | - +-----------------------+----------------------------------------+ - | ``density`` | State Density (:math:`\frac{g}{cm^3}`) | - +-----------------------+----------------------------------------+ - | ``temperature`` | State Temperature (K) | - +-----------------------+----------------------------------------+ - | ``massfractions(i)`` | Mass Fraction for element i | - +-----------------------+----------------------------------------+ + given a $\rho$, $T$, and $X_k$, call the equation of state and print out + the thermodynamic information. -Running the Code ----------------- +* ``nse_table_cell`` : -To run the code, enter the burn_cell directory and run:: + given a $\rho$, $T$, and $Y_e$, evaluate the NSE state via table interpolation + and print it out. - ./main3d.gnu.exe with inputs +* ``test_ase`` : -where ``inputs`` is the name of your inputs file. + for the self-consistent NSE, take a $\rho$, $T$, and $Y_e$, and solve for the NSE + state. Then check the NSE condition to see if we are actually satisfying the NSE + criteria for the network. -For each of the ``numsteps`` steps defined in the inputs -file, the code will output a files into a new directory titled -``run_prefix_output`` where ``run_prefix`` is the run prefix defined in the -inputs file. Each output file will be named using the run prefix -defined in the inputs file and the corresponding timestep. +* ``test_part_func`` -Next, run ``burn_cell.py`` using python 3.x, giving the defined run prefix as an argument. -For example:: + exercise the partition function interpolation for a few select nuclei. - python3 burn_cell.py react_aprox13 -The ``burn_cell.py`` code will gather information from all of the -output files and compile them into three graphs explained below. +Infrastructure tests +==================== -Graphs Output by ``burn_cell.py`` ---------------------------------- +* ``test_linear_algebra`` : -The file ``run-prefix_logX.png`` and ``run-prefix_logX.eps`` will display a -graph of the chemical abundances as a function of the time, both on -logarithmic scales, for all species involved in the simulation. An -example of this graph is shown below. + create a diagonally dominant matrix, multiply it by a test vector, $x$, + to get $b = Ax$, and then call the linear algebra routines to see if we + we recover $x$ from $b$. -.. figure:: react_aprox13_logX.png - :alt: An example of a plot output by the burn_cell unit test. This is the logX output corresponding to the network aprox13. - :width: 4.5in +* ``test_nse_interp`` : - An example of a plot output by the burn_cell unit test. This is the - logX output corresponding to the network aprox13. + run various tests of the NSE interpolation routines. +* ``test_parameters`` : + a simple setup that initializes the runtime parameters and can be + used to test if we can override them at runtime via inputs or the + commandline. This uses both the global data and the struct form + of the runtime parameters. -The file ``run-prefix_ydot.png`` and ``run-prefix_ydot.eps`` will display the -molar fraction (mass fraction / atomic weight) as a function of time, -both on logarithmic scales, for all species involved in the code. +* ``test_sdc_vode_rhs`` : -The file ``run-prefix_T-edot.png`` and ``run-prefix_T-edot.eps`` will display -the temperature and the energy generation rate as a function of time. + a simple driver for the SDC RHS routines. Given a thermodynamic + state, it outputs the RHS that the integrator will see.