Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
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
1 change: 1 addition & 0 deletions networks/CNO_extras/Make.package
Original file line number Diff line number Diff line change
Expand Up @@ -4,6 +4,7 @@ ifeq ($(USE_REACT),TRUE)
CEXE_sources += actual_network_data.cpp
CEXE_headers += actual_network.H
CEXE_headers += tfactors.H
CEXE_headers += interp_tools.H
CEXE_headers += partition_functions.H
CEXE_sources += partition_functions_data.cpp
CEXE_headers += actual_rhs.H
Expand Down
10 changes: 1 addition & 9 deletions networks/CNO_extras/actual_network.H
Original file line number Diff line number Diff line change
Expand Up @@ -289,14 +289,6 @@ namespace Rates
NumRates = k_He4_Mg24_to_C12_O16_derived
};

// number of reaclib rates

const int NrateReaclib = 63;

// number of tabular rates

const int NrateTabular = 8;

// rate names -- note: the rates are 1-based, not zero-based, so we pad
// this vector with rate_names[0] = "" so the indices line up with the
// NetworkRates enum
Expand Down Expand Up @@ -389,7 +381,7 @@ namespace NSE_INDEX
// First 3 row indices for reactants, followed by 3 product indices
// last index is the corresponding reverse rate index.

extern AMREX_GPU_MANAGED amrex::Array2D<int, 1, Rates::NumRates, 1, 7, amrex::Order::C> rate_indices;
extern AMREX_GPU_MANAGED amrex::Array2D<std::int8_t, 1, Rates::NumRates, 1, 7, amrex::Order::C> rate_indices;
}
#endif

Expand Down
144 changes: 72 additions & 72 deletions networks/CNO_extras/actual_network_data.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -4,78 +4,78 @@
#ifdef NSE_NET
namespace NSE_INDEX
{
AMREX_GPU_MANAGED amrex::Array2D<int, 1, Rates::NumRates, 1, 7, amrex::Order::C> rate_indices {
-1, -1, 4, -1, -1, 3, -1,
-1, -1, 7, -1, -1, 5, -1,
-1, -1, 8, -1, -1, 6, -1,
-1, 0, 2, -1, -1, 4, 42,
-1, 1, 2, -1, -1, 9, 47,
-1, 0, 3, -1, -1, 5, 43,
-1, 0, 4, -1, -1, 7, 44,
-1, 0, 5, -1, -1, 8, 45,
-1, 1, 5, -1, -1, 13, 50,
-1, 0, 6, -1, -1, 9, 46,
-1, 1, 6, -1, -1, 14, 52,
-1, 1, 7, -1, -1, 15, 54,
-1, 1, 8, -1, -1, 16, 56,
-1, 0, 9, -1, -1, 12, 48,
-1, 1, 9, -1, -1, 17, 58,
-1, 0, 10, -1, -1, 13, 49,
-1, 0, 11, -1, -1, 14, 51,
-1, 0, 12, -1, -1, 15, 53,
-1, 0, 13, -1, -1, 16, 55,
-1, 0, 14, -1, -1, 17, 57,
-1, 1, 15, -1, -1, 18, 59,
-1, 1, 17, -1, -1, 19, 60,
-1, 2, 2, -1, 1, 17, 70,
-1, 1, 4, -1, 0, 9, 66,
-1, 0, 6, -1, 1, 2, 62,
-1, 1, 7, -1, 0, 12, 68,
-1, 2, 9, -1, 1, 19, 71,
-1, 0, 10, -1, 1, 5, 63,
-1, 0, 11, -1, 1, 6, 64,
-1, 0, 13, -1, 1, 8, 65,
-1, 0, 14, -1, 1, 9, 67,
-1, 0, 17, -1, 1, 12, -1,
1, 1, 1, -1, -1, 2, 61,
-1, -1, 12, -1, -1, 10, 35,
-1, -1, 10, -1, -1, 12, -1,
-1, -1, 13, -1, -1, 15, -1,
-1, -1, 13, -1, -1, 11, 39,
-1, -1, 15, -1, -1, 13, 36,
-1, -1, 11, -1, -1, 13, -1,
-1, -1, 14, -1, -1, 16, -1,
-1, -1, 16, -1, -1, 14, 40,
-1, -1, 4, -1, 0, 2, -1,
-1, -1, 5, -1, 0, 3, -1,
-1, -1, 7, -1, 0, 4, -1,
-1, -1, 8, -1, 0, 5, -1,
-1, -1, 9, -1, 0, 6, -1,
-1, -1, 9, -1, 1, 2, -1,
-1, -1, 12, -1, 0, 9, -1,
-1, -1, 13, -1, 0, 10, -1,
-1, -1, 13, -1, 1, 5, -1,
-1, -1, 14, -1, 0, 11, -1,
-1, -1, 14, -1, 1, 6, -1,
-1, -1, 15, -1, 0, 12, -1,
-1, -1, 15, -1, 1, 7, -1,
-1, -1, 16, -1, 0, 13, -1,
-1, -1, 16, -1, 1, 8, -1,
-1, -1, 17, -1, 0, 14, -1,
-1, -1, 17, -1, 1, 9, -1,
-1, -1, 18, -1, 1, 15, -1,
-1, -1, 19, -1, 1, 17, -1,
-1, -1, 2, 1, 1, 1, -1,
-1, 1, 2, -1, 0, 6, -1,
-1, 1, 5, -1, 0, 10, -1,
-1, 1, 6, -1, 0, 11, -1,
-1, 1, 8, -1, 0, 13, -1,
-1, 0, 9, -1, 1, 4, -1,
-1, 1, 9, -1, 0, 14, -1,
-1, 0, 12, -1, 1, 7, -1,
-1, 1, 12, -1, 0, 17, 32,
-1, 1, 17, -1, 2, 2, -1,
-1, 1, 19, -1, 2, 9, -1
AMREX_GPU_MANAGED amrex::Array2D<std::int8_t, 1, Rates::NumRates, 1, 7, amrex::Order::C> rate_indices {
-1, -1, 4, -1, -1, 3, -1, // N13_to_C13_weak_wc12
-1, -1, 7, -1, -1, 5, -1, // O14_to_N14_weak_wc12
-1, -1, 8, -1, -1, 6, -1, // O15_to_N15_weak_wc12
-1, 0, 2, -1, -1, 4, 42, // p_C12_to_N13
-1, 1, 2, -1, -1, 9, 47, // He4_C12_to_O16
-1, 0, 3, -1, -1, 5, 43, // p_C13_to_N14
-1, 0, 4, -1, -1, 7, 44, // p_N13_to_O14
-1, 0, 5, -1, -1, 8, 45, // p_N14_to_O15
-1, 1, 5, -1, -1, 13, 50, // He4_N14_to_F18
-1, 0, 6, -1, -1, 9, 46, // p_N15_to_O16
-1, 1, 6, -1, -1, 14, 52, // He4_N15_to_F19
-1, 1, 7, -1, -1, 15, 54, // He4_O14_to_Ne18
-1, 1, 8, -1, -1, 16, 56, // He4_O15_to_Ne19
-1, 0, 9, -1, -1, 12, 48, // p_O16_to_F17
-1, 1, 9, -1, -1, 17, 58, // He4_O16_to_Ne20
-1, 0, 10, -1, -1, 13, 49, // p_O17_to_F18
-1, 0, 11, -1, -1, 14, 51, // p_O18_to_F19
-1, 0, 12, -1, -1, 15, 53, // p_F17_to_Ne18
-1, 0, 13, -1, -1, 16, 55, // p_F18_to_Ne19
-1, 0, 14, -1, -1, 17, 57, // p_F19_to_Ne20
-1, 1, 15, -1, -1, 18, 59, // He4_Ne18_to_Mg22
-1, 1, 17, -1, -1, 19, 60, // He4_Ne20_to_Mg24
-1, 2, 2, -1, 1, 17, 70, // C12_C12_to_He4_Ne20
-1, 1, 4, -1, 0, 9, 66, // He4_N13_to_p_O16
-1, 0, 6, -1, 1, 2, 62, // p_N15_to_He4_C12
-1, 1, 7, -1, 0, 12, 68, // He4_O14_to_p_F17
-1, 2, 9, -1, 1, 19, 71, // C12_O16_to_He4_Mg24
-1, 0, 10, -1, 1, 5, 63, // p_O17_to_He4_N14
-1, 0, 11, -1, 1, 6, 64, // p_O18_to_He4_N15
-1, 0, 13, -1, 1, 8, 65, // p_F18_to_He4_O15
-1, 0, 14, -1, 1, 9, 67, // p_F19_to_He4_O16
-1, 0, 17, -1, 1, 12, -1, // p_Ne20_to_He4_F17
1, 1, 1, -1, -1, 2, 61, // He4_He4_He4_to_C12
-1, -1, 12, -1, -1, 10, 35, // F17_to_O17
-1, -1, 10, -1, -1, 12, -1, // O17_to_F17
-1, -1, 13, -1, -1, 15, -1, // F18_to_Ne18
-1, -1, 13, -1, -1, 11, 39, // F18_to_O18
-1, -1, 15, -1, -1, 13, 36, // Ne18_to_F18
-1, -1, 11, -1, -1, 13, -1, // O18_to_F18
-1, -1, 14, -1, -1, 16, -1, // F19_to_Ne19
-1, -1, 16, -1, -1, 14, 40, // Ne19_to_F19
-1, -1, 4, -1, 0, 2, -1, // N13_to_p_C12_derived
-1, -1, 5, -1, 0, 3, -1, // N14_to_p_C13_derived
-1, -1, 7, -1, 0, 4, -1, // O14_to_p_N13_derived
-1, -1, 8, -1, 0, 5, -1, // O15_to_p_N14_derived
-1, -1, 9, -1, 0, 6, -1, // O16_to_p_N15_derived
-1, -1, 9, -1, 1, 2, -1, // O16_to_He4_C12_derived
-1, -1, 12, -1, 0, 9, -1, // F17_to_p_O16_derived
-1, -1, 13, -1, 0, 10, -1, // F18_to_p_O17_derived
-1, -1, 13, -1, 1, 5, -1, // F18_to_He4_N14_derived
-1, -1, 14, -1, 0, 11, -1, // F19_to_p_O18_derived
-1, -1, 14, -1, 1, 6, -1, // F19_to_He4_N15_derived
-1, -1, 15, -1, 0, 12, -1, // Ne18_to_p_F17_derived
-1, -1, 15, -1, 1, 7, -1, // Ne18_to_He4_O14_derived
-1, -1, 16, -1, 0, 13, -1, // Ne19_to_p_F18_derived
-1, -1, 16, -1, 1, 8, -1, // Ne19_to_He4_O15_derived
-1, -1, 17, -1, 0, 14, -1, // Ne20_to_p_F19_derived
-1, -1, 17, -1, 1, 9, -1, // Ne20_to_He4_O16_derived
-1, -1, 18, -1, 1, 15, -1, // Mg22_to_He4_Ne18_derived
-1, -1, 19, -1, 1, 17, -1, // Mg24_to_He4_Ne20_derived
-1, -1, 2, 1, 1, 1, -1, // C12_to_He4_He4_He4_derived
-1, 1, 2, -1, 0, 6, -1, // He4_C12_to_p_N15_derived
-1, 1, 5, -1, 0, 10, -1, // He4_N14_to_p_O17_derived
-1, 1, 6, -1, 0, 11, -1, // He4_N15_to_p_O18_derived
-1, 1, 8, -1, 0, 13, -1, // He4_O15_to_p_F18_derived
-1, 0, 9, -1, 1, 4, -1, // p_O16_to_He4_N13_derived
-1, 1, 9, -1, 0, 14, -1, // He4_O16_to_p_F19_derived
-1, 0, 12, -1, 1, 7, -1, // p_F17_to_He4_O14_derived
-1, 1, 12, -1, 0, 17, 32, // He4_F17_to_p_Ne20_derived
-1, 1, 17, -1, 2, 2, -1, // He4_Ne20_to_C12_C12_derived
-1, 1, 19, -1, 2, 9, -1 // He4_Mg24_to_C12_O16_derived
};
}
#endif
Expand Down
45 changes: 26 additions & 19 deletions networks/CNO_extras/actual_rhs.H
Original file line number Diff line number Diff line change
Expand Up @@ -38,9 +38,10 @@ void ener_gener_rate(T const& dydt, amrex::Real& enuc)

enuc = 0.0_rt;

for (int n = 1; n <= NumSpec; ++n) {
enuc += dydt(n) * network::mion(n);
}
amrex::constexpr_for<1, NumSpec+1>([&] (auto n)
{
enuc += dydt(n) * network::mion<n>();
});

enuc *= C::enuc_conv2;
}
Expand Down Expand Up @@ -563,64 +564,67 @@ void evaluate_rates(const burn_t& state,

rate_eval.enuc_weak = 0.0_rt;

amrex::Real log_temp = std::log10(state.T);
amrex::Real log_rhoy = std::log10(rhoy);

tabular_evaluate(j_F17_O17_meta, j_F17_O17_rhoy, j_F17_O17_temp, j_F17_O17_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F17_to_O17) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_F17_to_O17) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(F17) * (edot_nu + edot_gamma);

tabular_evaluate(j_O17_F17_meta, j_O17_F17_rhoy, j_O17_F17_temp, j_O17_F17_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_O17_to_F17) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_O17_to_F17) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(O17) * (edot_nu + edot_gamma);

tabular_evaluate(j_F18_Ne18_meta, j_F18_Ne18_rhoy, j_F18_Ne18_temp, j_F18_Ne18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F18_to_Ne18) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_F18_to_Ne18) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(F18) * (edot_nu + edot_gamma);

tabular_evaluate(j_F18_O18_meta, j_F18_O18_rhoy, j_F18_O18_temp, j_F18_O18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F18_to_O18) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_F18_to_O18) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(F18) * (edot_nu + edot_gamma);

tabular_evaluate(j_Ne18_F18_meta, j_Ne18_F18_rhoy, j_Ne18_F18_temp, j_Ne18_F18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_Ne18_to_F18) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_Ne18_to_F18) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(Ne18) * (edot_nu + edot_gamma);

tabular_evaluate(j_O18_F18_meta, j_O18_F18_rhoy, j_O18_F18_temp, j_O18_F18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_O18_to_F18) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_O18_to_F18) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(O18) * (edot_nu + edot_gamma);

tabular_evaluate(j_F19_Ne19_meta, j_F19_Ne19_rhoy, j_F19_Ne19_temp, j_F19_Ne19_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F19_to_Ne19) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_F19_to_Ne19) = drate_dt;
}
rate_eval.enuc_weak += C::n_A * Y(F19) * (edot_nu + edot_gamma);

tabular_evaluate(j_Ne19_F19_meta, j_Ne19_F19_rhoy, j_Ne19_F19_temp, j_Ne19_F19_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_Ne19_to_F19) = rate;
if constexpr (std::is_same_v<T, rate_derivs_t>) {
rate_eval.dscreened_rates_dT(k_Ne19_to_F19) = drate_dt;
Expand Down Expand Up @@ -675,43 +679,46 @@ void get_ydot_weak(const burn_t& state,

// Calculate tabular rates and get ydot_weak

amrex::Real log_temp = std::log10(state.T);
amrex::Real log_rhoy = std::log10(rhoy);

tabular_evaluate(j_F17_O17_meta, j_F17_O17_rhoy, j_F17_O17_temp, j_F17_O17_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F17_to_O17) = rate;
rate_eval.enuc_weak += C::n_A * Y(F17) * (edot_nu + edot_gamma);

tabular_evaluate(j_O17_F17_meta, j_O17_F17_rhoy, j_O17_F17_temp, j_O17_F17_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_O17_to_F17) = rate;
rate_eval.enuc_weak += C::n_A * Y(O17) * (edot_nu + edot_gamma);

tabular_evaluate(j_F18_Ne18_meta, j_F18_Ne18_rhoy, j_F18_Ne18_temp, j_F18_Ne18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F18_to_Ne18) = rate;
rate_eval.enuc_weak += C::n_A * Y(F18) * (edot_nu + edot_gamma);

tabular_evaluate(j_F18_O18_meta, j_F18_O18_rhoy, j_F18_O18_temp, j_F18_O18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F18_to_O18) = rate;
rate_eval.enuc_weak += C::n_A * Y(F18) * (edot_nu + edot_gamma);

tabular_evaluate(j_Ne18_F18_meta, j_Ne18_F18_rhoy, j_Ne18_F18_temp, j_Ne18_F18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_Ne18_to_F18) = rate;
rate_eval.enuc_weak += C::n_A * Y(Ne18) * (edot_nu + edot_gamma);

tabular_evaluate(j_O18_F18_meta, j_O18_F18_rhoy, j_O18_F18_temp, j_O18_F18_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_O18_to_F18) = rate;
rate_eval.enuc_weak += C::n_A * Y(O18) * (edot_nu + edot_gamma);

tabular_evaluate(j_F19_Ne19_meta, j_F19_Ne19_rhoy, j_F19_Ne19_temp, j_F19_Ne19_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_F19_to_Ne19) = rate;
rate_eval.enuc_weak += C::n_A * Y(F19) * (edot_nu + edot_gamma);

tabular_evaluate(j_Ne19_F19_meta, j_Ne19_F19_rhoy, j_Ne19_F19_temp, j_Ne19_F19_data,
rhoy, state.T, rate, drate_dt, edot_nu, edot_gamma);
log_rhoy, log_temp, state.T, rate, drate_dt, edot_nu, edot_gamma);
rate_eval.screened_rates(k_Ne19_to_F19) = rate;
rate_eval.enuc_weak += C::n_A * Y(Ne19) * (edot_nu + edot_gamma);

Expand Down
66 changes: 66 additions & 0 deletions networks/CNO_extras/interp_tools.H
Original file line number Diff line number Diff line change
@@ -0,0 +1,66 @@
#ifndef INTERP_TOOLS_H
#define INTERP_TOOLS_H

#include <AMReX_REAL.H>

namespace interp_net {

// index locator
// return the index i such that x_array[i] <= x0 <= x_array[i+1]

template <typename T>
AMREX_GPU_HOST_DEVICE AMREX_INLINE
int
find_index(const amrex::Real x0, const T& x_array) {

int left = x_array.lo();
int right = x_array.hi();

int idx = -1;

if (x0 >= x_array(left) && x0 < x_array(right)) {

// find the largest x element <= x0 using a
// binary search

while (left < right) {
int mid = (left + right) / 2;
if (x_array(mid) > x0) {
right = mid;
} else {
left = mid + 1;
}
}

idx = right - 1;
}

return idx;
}

// an index locator that extrapolates at the boundaries instead
// of returning -1

template <typename T>
AMREX_GPU_HOST_DEVICE AMREX_INLINE
int
find_index_extrap(const amrex::Real x0, const T& x_array) {

int left = x_array.lo();
int right = x_array.hi();

int idx;

if (x0 < x_array(left)) {
idx = left;
} else if (x0 > x_array(right)) {
idx = right - 1;
} else {
idx = find_index(x0, x_array);
}

return idx;
}

}
#endif
Loading
Loading