From cbedde4a4625f8eccbe41690d95ac1743009d177 Mon Sep 17 00:00:00 2001 From: Luke Shingles Date: Wed, 20 Nov 2024 14:02:50 +0000 Subject: [PATCH] Update nltepop.cc --- nltepop.cc | 28 ++++++++++++++-------------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/nltepop.cc b/nltepop.cc index 21f903aff..70907e200 100644 --- a/nltepop.cc +++ b/nltepop.cc @@ -15,6 +15,7 @@ #include #include #include +#include #include #include @@ -450,10 +451,8 @@ void nltepop_matrix_add_boundbound(const int nonemptymgi, const int element, con // de-excitation const int ndowntrans = get_ndowntrans(element, ion, level); - const auto *const leveldowntranslist = get_downtranslist(element, ion, level); - const auto ndowntransindices = std::ranges::iota_view{0, ndowntrans}; - std::for_each(EXEC_PAR ndowntransindices.begin(), ndowntransindices.end(), [&](const auto i) { - const auto &downtransition = leveldowntranslist[i]; + const auto leveldowntrans = std::span(get_downtranslist(element, ion, level), ndowntrans); + std::for_each(EXEC_PAR leveldowntrans.begin(), leveldowntrans.end(), [&](const auto &downtransition) { const double A_ul = downtransition.einstein_A; const int lower = downtransition.targetlevelindex; @@ -529,7 +528,8 @@ void nltepop_matrix_add_ionisation(const int modelgridindex, const int element, const int nionisinglevels = get_nlevels_ionising(element, ion); const int maxrecombininglevel = get_maxrecombininglevel(element, ion + 1); - for (int level = 0; level < nionisinglevels; level++) { + const auto levels = std::ranges::iota_view{0, nionisinglevels}; + std::for_each(EXEC_PAR levels.begin(), levels.end(), [&](const auto level) { const int lower_index = get_nlte_vector_index(element, ion, level); // thermal collisional ionization, photoionisation and recombination processes @@ -548,10 +548,10 @@ void nltepop_matrix_add_ionisation(const int modelgridindex, const int element, const double C_ionisation = col_ionization_ratecoeff(T_e, nne, element, ion, level, phixstargetindex, epsilon_trans); - *gsl_matrix_ptr(rate_matrix_rad_bf, lower_index, lower_index) -= R_ionisation * s_renorm[level]; - *gsl_matrix_ptr(rate_matrix_rad_bf, upper_index, lower_index) += R_ionisation * s_renorm[level]; - *gsl_matrix_ptr(rate_matrix_coll_bf, lower_index, lower_index) -= C_ionisation * s_renorm[level]; - *gsl_matrix_ptr(rate_matrix_coll_bf, upper_index, lower_index) += C_ionisation * s_renorm[level]; + atomicadd(*gsl_matrix_ptr(rate_matrix_rad_bf, lower_index, lower_index), -R_ionisation * s_renorm[level]); + atomicadd(*gsl_matrix_ptr(rate_matrix_rad_bf, upper_index, lower_index), R_ionisation * s_renorm[level]); + atomicadd(*gsl_matrix_ptr(rate_matrix_coll_bf, lower_index, lower_index), -C_ionisation * s_renorm[level]); + atomicadd(*gsl_matrix_ptr(rate_matrix_coll_bf, upper_index, lower_index), C_ionisation * s_renorm[level]); if ((R_ionisation < 0) || (C_ionisation < 0)) { printout(" WARNING: Negative ionization rate from ionstage %d level %d phixstargetindex %d\n", @@ -564,10 +564,10 @@ void nltepop_matrix_add_ionisation(const int modelgridindex, const int element, const double R_recomb = rad_recombination_ratecoeff(T_e, nne, element, ion + 1, upper, level, nonemptymgi); const double C_recomb = col_recombination_ratecoeff(T_e, nne, element, ion + 1, upper, level, epsilon_trans); - *gsl_matrix_ptr(rate_matrix_rad_bf, upper_index, upper_index) -= R_recomb * s_renorm[upper]; - *gsl_matrix_ptr(rate_matrix_rad_bf, lower_index, upper_index) += R_recomb * s_renorm[upper]; - *gsl_matrix_ptr(rate_matrix_coll_bf, upper_index, upper_index) -= C_recomb * s_renorm[upper]; - *gsl_matrix_ptr(rate_matrix_coll_bf, lower_index, upper_index) += C_recomb * s_renorm[upper]; + atomicadd(*gsl_matrix_ptr(rate_matrix_rad_bf, upper_index, upper_index), -R_recomb * s_renorm[upper]); + atomicadd(*gsl_matrix_ptr(rate_matrix_rad_bf, lower_index, upper_index), R_recomb * s_renorm[upper]); + atomicadd(*gsl_matrix_ptr(rate_matrix_coll_bf, upper_index, upper_index), -C_recomb * s_renorm[upper]); + atomicadd(*gsl_matrix_ptr(rate_matrix_coll_bf, lower_index, upper_index), C_recomb * s_renorm[upper]); if ((R_recomb < 0) || (C_recomb < 0)) { printout(" WARNING: Negative recombination rate to ionstage %d level %d phixstargetindex %d\n", @@ -575,7 +575,7 @@ void nltepop_matrix_add_ionisation(const int modelgridindex, const int element, } } } - } + }); } void nltepop_matrix_add_nt_ionisation(const int nonemptymgi, const int element, const int ion,