Skip to content

Fix RHS/Jac coefficients when a species is on both sides #2054

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

Draft
maximumcats wants to merge 2 commits into
base: development
Choose a base branch
from
Draft

Fix RHS/Jac coefficients when a species is on both sides #2054

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
7ca9367
Fix RHS/Jac coefficients when a species is on both sides
maximumcats Sep 5, 2026
fa034fd
Update test_rhs
maximumcats Sep 5, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
53 changes: 24 additions & 29 deletions networks/rhs.H
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -829,35 +829,30 @@ constexpr std::pair<amrex::Real, amrex::Real> rhs_term (const burn_t& state, con

// Now compute the total contribution to this species.

amrex::Real coefficient = 0.0_rt;

if constexpr (data.species_A == species) {
forward_term *= -data.number_A;
reverse_term *= data.number_A;
}
else if constexpr (data.species_B == species) {
forward_term *= -data.number_B;
reverse_term *= data.number_B;
coefficient -= data.number_A;
}
else if constexpr (data.species_C == species) {
forward_term *= -data.number_C;
reverse_term *= data.number_C;
if constexpr (data.species_B == species) {
coefficient -= data.number_B;
}
else if constexpr (data.species_D == species) {
forward_term *= data.number_D;
reverse_term *= -data.number_D;
if constexpr (data.species_C == species) {
coefficient -= data.number_C;
}
else if constexpr (data.species_E == species) {
forward_term *= data.number_E;
reverse_term *= -data.number_E;
if constexpr (data.species_D == species) {
coefficient += data.number_D;
}
else if constexpr (data.species_F == species) {
forward_term *= data.number_F;
reverse_term *= -data.number_F;
if constexpr (data.species_E == species) {
coefficient += data.number_E;
}
else {
forward_term = 0.0;
reverse_term = 0.0;
if constexpr (data.species_F == species) {
coefficient += data.number_F;
}

forward_term *= coefficient;
reverse_term *= -coefficient;

return {forward_term, reverse_term};
}

Expand Down Expand Up @@ -1100,28 +1095,28 @@ constexpr amrex::Real jac_term (const burn_t& state, const rate_t<autodiff::dual

// Now compute the total contribution to this species.

amrex::Real term = 0.0_rt;
amrex::Real coefficient = 0.0_rt;

if constexpr (data.species_A == spec1) {
term = data.number_A * (reverse_term - forward_term);
coefficient -= data.number_A;
}
if constexpr (data.species_B == spec1) {
term = data.number_B * (reverse_term - forward_term);
coefficient -= data.number_B;
}
if constexpr (data.species_C == spec1) {
term = data.number_C * (reverse_term - forward_term);
coefficient -= data.number_C;
}
if constexpr (data.species_D == spec1) {
term = data.number_D * (forward_term - reverse_term);
coefficient += data.number_D;
}
if constexpr (data.species_E == spec1) {
term = data.number_E * (forward_term - reverse_term);
coefficient += data.number_E;
}
if constexpr (data.species_F == spec1) {
term = data.number_F * (forward_term - reverse_term);
coefficient += data.number_F;
}

return term;
return coefficient * (forward_term - reverse_term);
}

template<int rate, typename number_t>
Expand Down
18 changes: 9 additions & 9 deletions unit_test/test_rhs/ci-benchmarks/aprox19.out
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@
temperature 50000000 5000000000
Ydot_hydrogen-1 8.6607464632e-05 1.3073380344e+12
Ydot_helium-3 -1.3076161972e+12 -8.6943749398e-05
Ydot_helium-4 4.3237314804e-05 6.538087849e+11
Ydot_helium-4 4.3550105477e-05 6.538087849e+11
Ydot_carbon-12 -1105871714.6 1868.0861291
Ydot_nitrogen-14 1.1710238604e-08 428182473.08
Ydot_oxygen-16 -3263960.7002 215410.42728
Expand Down Expand Up @@ -41,7 +41,7 @@
Xold_nickel-56 0.00625 0.0375
Xold_neutron 0.00625 0.0375
Xold_proton 0.00625 0.0375
Edot 7.5242949946e+14 8.1592675307e+30
Edot 5.4153457593e+14 8.1592675306e+30
J_hydrogen-1_hydrogen-1 -59732766098 -3.7590547944e-06
J_helium-3_hydrogen-1 2.586350502e-16 0.04958554103
J_helium-4_hydrogen-1 5.2975868011e-09 0.07375836734
Expand All @@ -64,7 +64,7 @@
J_E_hydrogen-1 2.1237371823e+13 3.3757896591e+29
J_hydrogen-1_helium-3 0.083240790137 2.0921859156e+14
J_helium-3_helium-3 -2.0921859156e+14 -0.08339092966
J_helium-4_helium-3 0.041808069473 1.0460929578e+14
J_helium-4_helium-3 0.041732999711 1.0460929578e+14
J_carbon-12_helium-3 0 0
J_nitrogen-14_helium-3 0 0
J_oxygen-16_helium-3 0 0
Expand All @@ -81,10 +81,10 @@
J_nickel-56_helium-3 0 0
J_neutron_helium-3 0 0
J_proton_helium-3 0 0
J_E_helium-3 2.4825899081e+17 1.2979583538e+33
J_E_helium-3 5.183119541e+17 1.2979583538e+33
J_hydrogen-1_helium-4 -0.448 -8.9368763863e-07
J_helium-3_helium-4 -0.448 -8.9368763863e-07
J_helium-4_helium-4 -573505840.55 0.896
J_helium-4_helium-4 -573505841 0.448
J_carbon-12_helium-4 -28066508.652 151095.37704
J_nitrogen-14_helium-4 -107832532.98 -3.5509584119e-31
J_oxygen-16_helium-4 -52056348.745 20152253.695
Expand All @@ -101,7 +101,7 @@
J_nickel-56_helium-4 2.0673829201e-80 3548434.8158
J_neutron_helium-4 0 3.1564509005e-05
J_proton_helium-4 0 474196.3415
J_E_helium-4 -1.6030606197e+21 4.9109161152e+27
J_E_helium-4 1.7069221056e+13 4.9109177269e+27
J_hydrogen-1_carbon-12 -39962772073 -2.8113415685e-06
J_helium-3_carbon-12 0 0
J_helium-4_carbon-12 -543156.94794 16180581376
Expand Down Expand Up @@ -424,13 +424,13 @@
J_E_proton -7.593501065e+15 4.237794952e+30
J_hydrogen-1_E 2.056391585e-19 5.7718078776e-06
J_helium-3_E -5.7716864292e-06 -2.0648367761e-19
J_helium-4_E 1.0265571946e-19 2.885937043e-06
J_helium-4_E 1.0343734791e-19 2.885937043e-06
J_carbon-12_E -1.7761025826e-08 3.1733884044e-10
J_nitrogen-14_E -3.3419345218e-10 9.9166452811e-09
J_oxygen-16_E -4.8278139961e-11 2.1225767505e-12
J_neon-20_E -7.3202938181e-14 8.8914523939e-09
J_magnesium-24_E -3.3222436154e-11 1.7090590032e-11
J_silicon-28_E -7.495337146e-48 6.2608476694e-11
J_silicon-28_E -7.4953371464e-48 6.2608476694e-11
J_sulfur-32_E -4.6637773178e-12 1.1717119105e-13
J_argon-36_E -1.1893314146e-13 5.8730727453e-12
J_calcium-40_E 8.5459256525e-77 8.7687490767e-12
Expand All @@ -441,5 +441,5 @@
J_nickel-56_E -1.41049225e-12 3.8177441789e-11
J_neutron_E 0 1.1761800714e-07
J_proton_E -5.6936091804e-14 1.1296574697e-07
J_E_E 2.0925380937 3.5786987301e+13
J_E_E 1.2863865678 3.5786987301e+13

16 changes: 8 additions & 8 deletions unit_test/test_rhs/ci-benchmarks/aprox21.out
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@
temperature 50000000 5000000000
Ydot_hydrogen-1 -4219179576.2 1.0316760608e+12
Ydot_helium-3 -1.0331595558e+12 -6.7832225696e-05
Ydot_helium-4 3.3710523725e-05 5.1667697447e+11
Ydot_helium-4 3.3985000991e-05 5.1667697447e+11
Ydot_carbon-12 -3487533699.9 -1.8356613712e-07
Ydot_nitrogen-14 1.8357874205e-07 2974640611.6
Ydot_oxygen-16 -36005.322458 632549.00147
Expand Down Expand Up @@ -45,7 +45,7 @@
Xold_nickel-56 0.0055555555556 0.033333333333
Xold_neutron 0.0055555555556 0.033333333333
Xold_proton 0.0055555555556 0.033333333333
Edot 5.8765130693e+14 6.8107090088e+30
Edot 4.2467149955e+14 6.8107090087e+30
J_hydrogen-1_hydrogen-1 -59494588148 -3.6818237055e-06
J_helium-3_hydrogen-1 4.1248545523e-15 0.10494648684
J_helium-4_hydrogen-1 2.6882623706e-09 0.00030487244657
Expand All @@ -70,7 +70,7 @@
J_E_hydrogen-1 2.0815071597e+13 3.3623290626e+29
J_hydrogen-1_helium-3 0.073036477167 1.8596872004e+14
J_helium-3_helium-3 -1.8596872004e+14 -0.07318469489
J_helium-4_helium-3 0.036703510738 9.298436002e+13
J_helium-4_helium-3 0.036629401876 9.298436002e+13
J_carbon-12_helium-3 0 0
J_nitrogen-14_helium-3 0 0
J_oxygen-16_helium-3 0 0
Expand All @@ -89,10 +89,10 @@
J_nickel-56_helium-3 0 0
J_neutron_helium-3 0 0
J_proton_helium-3 0 0
J_E_helium-3 1.8838546949e+17 1.1537199057e+33
J_E_helium-3 4.5498173017e+17 1.1537199057e+33
J_hydrogen-1_helium-4 -0.39822222222 -7.842207602e-07
J_helium-3_helium-4 -0.39822222222 -7.842207602e-07
J_helium-4_helium-4 -502897744.49 0.79644444444
J_helium-4_helium-4 -502897744.88 0.39822222222
J_carbon-12_helium-4 -27977366.316 146014.40756
J_nitrogen-14_helium-4 -95845076.964 -3.1070673892e-31
J_oxygen-16_helium-4 -705164.28103 27539882.999
Expand All @@ -111,7 +111,7 @@
J_nickel-56_helium-4 1.9427222337e-80 3128320.9725
J_neutron_helium-4 0 3.1570367276e-05
J_proton_helium-4 0 1435569.2733
J_E_helium-4 -1.4249427731e+21 4.3383246754e+27
J_E_helium-4 1.4978429749e+13 4.3383261079e+27
J_hydrogen-1_carbon-12 -6.9959804335e+11 -4.4117322214e-05
J_helium-3_carbon-12 0 0
J_helium-4_carbon-12 -536900.29748 16026654461
Expand Down Expand Up @@ -510,7 +510,7 @@
J_E_proton -2.6108991385e+16 3.6817330507e+30
J_hydrogen-1_E -1.0006504976e-07 4.2872891104e-06
J_helium-3_E -4.2866803458e-06 -7.5769167747e-20
J_helium-4_E 3.7852613992e-20 2.1441607162e-06
J_helium-4_E 3.7898785602e-20 2.1441607162e-06
J_carbon-12_E -5.3556320773e-08 1.4041391721e-09
J_nitrogen-14_E -1.9393729314e-09 5.3556258717e-08
J_oxygen-16_E -3.4898941611e-13 6.2789315765e-12
Expand All @@ -529,5 +529,5 @@
J_nickel-56_E -1.3650886568e-12 2.7798205721e-11
J_neutron_E -4.0660804174e-09 5.9256825143e-08
J_proton_E -3.9082377477e-10 8.5866665712e-08
J_E_E 0.65504112892 2.6535626445e+13
J_E_E 0.48894522114 2.6535626445e+13

Loading