Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Use amrex::Math::powi for RHS #709

Merged
merged 2 commits into from
Dec 31, 2023
Merged

Use amrex::Math::powi for RHS #709

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
5005cce
use powi for amrex rhs
zhichen3 Dec 21, 2023
1916630
update benchmark
zhichen3 Dec 21, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
100 changes: 50 additions & 50 deletions .github/workflows/microphysics-benchmarks/subch_approx_unit_test.out
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
Initializing AMReX (23.11-4-ga7afcba3cffd)...
AMReX (23.11-4-ga7afcba3cffd) initialized
Initializing AMReX (23.05-278-g75571e2dcbf2)...
AMReX (23.05-278-g75571e2dcbf2) initialized
starting the single zone burn...
reading in network electron-capture / beta-decay tables...
Maximum Time (s): 1e-05
Expand Down Expand Up @@ -52,59 +52,59 @@ RHS at t = 0
Ni56 2.574225196e-29
------------------------------------
successful? 1
- Hnuc = 7.242504668e+22
- added e = 7.242504668e+17
- final T = 3332738253
- Hnuc = 7.242979451e+22
- added e = 7.242979451e+17
- final T = 3332756660
------------------------------------
e initial = 1.396711859e+18
e final = 2.120962326e+18
e final = 2.121009804e+18
------------------------------------
new mass fractions:
H1 0.1076931392
He4 0.08271604755
C12 1.784659301e-05
N13 3.392571632e-08
N14 3.299221207e-07
O16 0.108645045
F18 1.417109763e-07
Ne20 0.007004350626
Ne21 5.73319595e-06
Na22 0.1571361633
Na23 4.823680511e-07
Mg24 0.004455549374
Al27 6.197826902e-06
Si28 0.06147651501
P31 9.545310692e-06
S32 0.1151668114
Ar36 0.09896332458
Ca40 0.2507603092
Ti44 0.005262916908
Cr48 0.0006665524196
Fe52 1.291730016e-05
Ni56 4.728298219e-08
H1 0.107693139
He4 0.08268984836
C12 1.785169978e-05
N13 3.393388377e-08
N14 3.301659593e-07
O16 0.1086277596
F18 1.41780265e-07
Ne20 0.007003004431
Ne21 5.733983976e-06
Na22 0.157136162
Na23 4.821438747e-07
Mg24 0.004454634907
Al27 6.194771841e-06
Si28 0.06146859376
P31 9.541393692e-06
S32 0.1151759569
Ar36 0.09898610611
Ca40 0.2507842611
Ti44 0.00526137162
Cr48 0.0006659118953
Fe52 1.289330863e-05
Ni56 4.714295333e-08
------------------------------------
species creation rates:
omegadot(H1): 769.3139202
omegadot(He4): -41728.39524
omegadot(C12): -9998.215341
omegadot(H1): 769.3138989
omegadot(He4): -41731.01516
omegadot(C12): -9998.21483
omegadot(N13): -9999.996607
omegadot(N14): -9999.967008
omegadot(O16): 864.5045001
omegadot(F18): 0.01417109763
omegadot(Ne20): 700.4350626
omegadot(Ne21): 0.573319595
omegadot(Na22): 15713.61633
omegadot(Na23): 0.04823680511
omegadot(Mg24): 445.5549374
omegadot(Al27): 0.6197826902
omegadot(Si28): 6147.651501
omegadot(P31): 0.9545310692
omegadot(S32): 11516.68114
omegadot(Ar36): 9896.332458
omegadot(Ca40): 25076.03092
omegadot(Ti44): 526.2916908
omegadot(Cr48): 66.65524196
omegadot(Fe52): 1.291730016
omegadot(Ni56): 0.004728298219
omegadot(N14): -9999.966983
omegadot(O16): 862.7759613
omegadot(F18): 0.0141780265
omegadot(Ne20): 700.3004431
omegadot(Ne21): 0.5733983976
omegadot(Na22): 15713.6162
omegadot(Na23): 0.04821438747
omegadot(Mg24): 445.4634907
omegadot(Al27): 0.6194771841
omegadot(Si28): 6146.859376
omegadot(P31): 0.9541393692
omegadot(S32): 11517.59569
omegadot(Ar36): 9898.610611
omegadot(Ca40): 25078.42611
omegadot(Ti44): 526.137162
omegadot(Cr48): 66.59118953
omegadot(Fe52): 1.289330863
omegadot(Ni56): 0.004714295333
number of steps taken: 594
AMReX (23.11-4-ga7afcba3cffd) finalized
AMReX (23.05-278-g75571e2dcbf2) finalized
9 changes: 9 additions & 0 deletions pynucastro/networks/amrexastro_cxx_network.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@

import glob
import os
import re

from pynucastro.networks.base_cxx_network import BaseCxxNetwork
from pynucastro.nucdata import Nucleus
Expand Down Expand Up @@ -58,6 +59,14 @@ def _rate_param_tests(self, n_indent, of):
of.write(f"{self.indent*n_indent} }}\n")
of.write(f"{self.indent*n_indent}}}\n\n")

def _cxxify(self, s):
# Replace std::pow(x, n) with amrex::Math::powi<n>(x) for amrexastro_cxx_network

cxx_code = super()._cxxify(s)
std_pow_pattern = r"std::pow\(([^,]+),\s*(\d+)\)"
amrex_powi_replacement = r"amrex::Math::powi<\2>(\1)"
return re.sub(std_pow_pattern, amrex_powi_replacement, cxx_code)

def _write_network(self, odir=None):
"""
This writes the RHS, jacobian and ancillary files for the system of ODEs that
Expand Down
10 changes: 7 additions & 3 deletions pynucastro/networks/base_cxx_network.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -351,6 +351,10 @@ def _compute_tabular_rates(self, n_indent, of):

of.write('\n')

def _cxxify(self, s):
# This is a helper function that converts sympy cxxcode to the actual c++ code we use.
return self.symbol_rates.cxxify(s)

def _ydot(self, n_indent, of):
# Write YDOT
for n in self.unique_nuclei:
Expand All @@ -375,7 +379,7 @@ def _ydot(self, n_indent, of):
of.write("(")

if pair[0] is not None:
sol_value = self.symbol_rates.cxxify(sympy.cxxcode(pair[0], precision=15,
sol_value = self._cxxify(sympy.cxxcode(pair[0], precision=15,
standard="c++11"))

of.write(f"{sol_value}")
Expand All @@ -384,7 +388,7 @@ def _ydot(self, n_indent, of):
of.write(" + ")

if pair[1] is not None:
sol_value = self.symbol_rates.cxxify(sympy.cxxcode(pair[1], precision=15,
sol_value = self._cxxify(sympy.cxxcode(pair[1], precision=15,
standard="c++11"))

of.write(f"{sol_value}")
Expand Down Expand Up @@ -417,7 +421,7 @@ def _jacnuc(self, n_indent, of):
for ini, ni in enumerate(self.unique_nuclei):
jac_idx = n_unique_nuclei*jnj + ini
if not self.jac_null_entries[jac_idx]:
jvalue = self.symbol_rates.cxxify(sympy.cxxcode(self.jac_out_result[jac_idx], precision=15,
jvalue = self._cxxify(sympy.cxxcode(self.jac_out_result[jac_idx], precision=15,
standard="c++11"))
of.write(f"{self.indent*(n_indent)}scratch = {jvalue};\n")
of.write(f"{self.indent*n_indent}jac.set({nj.cindex()}, {ni.cindex()}, scratch);\n\n")
Expand Down
18 changes: 9 additions & 9 deletions pynucastro/networks/tests/_amrexastro_cxx_reference/actual_rhs.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -164,37 +164,37 @@ void rhs_nuc(const burn_t& state,

ydot_nuc(N) =
-screened_rates(k_n_to_p_weak_wc12)*Y(N) +
0.5*screened_rates(k_C12_C12_to_n_Mg23)*std::pow(Y(C12), 2)*state.rho;
0.5*screened_rates(k_C12_C12_to_n_Mg23)*amrex::Math::powi<2>(Y(C12))*state.rho;

ydot_nuc(H1) =
0.5*screened_rates(k_C12_C12_to_p_Na23)*std::pow(Y(C12), 2)*state.rho +
0.5*screened_rates(k_C12_C12_to_p_Na23)*amrex::Math::powi<2>(Y(C12))*state.rho +
screened_rates(k_n_to_p_weak_wc12)*Y(N);

ydot_nuc(He4) =
0.5*screened_rates(k_C12_C12_to_He4_Ne20)*std::pow(Y(C12), 2)*state.rho +
0.5*screened_rates(k_C12_C12_to_He4_Ne20)*amrex::Math::powi<2>(Y(C12))*state.rho +
-screened_rates(k_He4_C12_to_O16)*Y(C12)*Y(He4)*state.rho;

ydot_nuc(C12) =
-screened_rates(k_C12_C12_to_He4_Ne20)*std::pow(Y(C12), 2)*state.rho +
-screened_rates(k_C12_C12_to_p_Na23)*std::pow(Y(C12), 2)*state.rho +
-screened_rates(k_C12_C12_to_He4_Ne20)*amrex::Math::powi<2>(Y(C12))*state.rho +
-screened_rates(k_C12_C12_to_p_Na23)*amrex::Math::powi<2>(Y(C12))*state.rho +
-screened_rates(k_He4_C12_to_O16)*Y(C12)*Y(He4)*state.rho +
-screened_rates(k_C12_C12_to_n_Mg23)*std::pow(Y(C12), 2)*state.rho;
-screened_rates(k_C12_C12_to_n_Mg23)*amrex::Math::powi<2>(Y(C12))*state.rho;

ydot_nuc(O16) =
screened_rates(k_He4_C12_to_O16)*Y(C12)*Y(He4)*state.rho;

ydot_nuc(Ne20) =
0.5*screened_rates(k_C12_C12_to_He4_Ne20)*std::pow(Y(C12), 2)*state.rho;
0.5*screened_rates(k_C12_C12_to_He4_Ne20)*amrex::Math::powi<2>(Y(C12))*state.rho;

ydot_nuc(Ne23) =
(-screened_rates(k_Ne23_to_Na23)*Y(Ne23) + screened_rates(k_Na23_to_Ne23)*Y(Na23));

ydot_nuc(Na23) =
0.5*screened_rates(k_C12_C12_to_p_Na23)*std::pow(Y(C12), 2)*state.rho +
0.5*screened_rates(k_C12_C12_to_p_Na23)*amrex::Math::powi<2>(Y(C12))*state.rho +
(screened_rates(k_Ne23_to_Na23)*Y(Ne23) + -screened_rates(k_Na23_to_Ne23)*Y(Na23));

ydot_nuc(Mg23) =
0.5*screened_rates(k_C12_C12_to_n_Mg23)*std::pow(Y(C12), 2)*state.rho;
0.5*screened_rates(k_C12_C12_to_n_Mg23)*amrex::Math::powi<2>(Y(C12))*state.rho;

}

Expand Down