Skip to content

Commit

Permalink
Changed MPI routines and noramlisation of BF estimators
Browse files Browse the repository at this point in the history
  • Loading branch information
@jpollin98
jpollin98 committed May 3, 2024
1 parent d9e7d90 commit 43d1bca
Show file tree
Hide file tree
Showing 3 changed files with 60 additions and 18 deletions.
18 changes: 14 additions & 4 deletions nltepop.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -126,8 +126,8 @@ static void filter_nlte_matrix(const int element, gsl_matrix *rate_matrix, gsl_v
}

static auto get_total_rate(const int index_selected, const gsl_matrix *rate_matrix, const gsl_vector *popvec,
const bool into_level, const bool only_levels_below,
const bool only_levels_above) -> double {
const bool into_level, const bool only_levels_below, const bool only_levels_above)
-> double {
double total_rate = 0.;
assert_always(!only_levels_below || !only_levels_above);

Expand Down Expand Up @@ -181,8 +181,8 @@ static auto get_total_rate_in(const int index_to, const gsl_matrix *rate_matrix,
return get_total_rate(index_to, rate_matrix, popvec, true, false, false);
}

static auto get_total_rate_out(const int index_from, const gsl_matrix *rate_matrix,
const gsl_vector *popvec) -> double {
static auto get_total_rate_out(const int index_from, const gsl_matrix *rate_matrix, const gsl_vector *popvec)
-> double {
return get_total_rate(index_from, rate_matrix, popvec, false, false, false);
}

Expand Down Expand Up @@ -808,6 +808,16 @@ void solve_nlte_pops_element(const int element, const int modelgridindex, const
return;
}

double cell_Te = grid::get_Te(modelgridindex);

if (cell_Te < 1001.) {
printout(
"Not solving for NLTE populations in cell %d at timestep %d for element Z=%d due to low temperature Te=%g\n",
modelgridindex, timestep, atomic_number, cell_Te);
set_element_pops_lte(modelgridindex, element);
return;
}

const auto sys_time_start_nltesolver = std::time(nullptr);

const double t_mid = globals::timesteps[timestep].mid;
Expand Down
37 changes: 26 additions & 11 deletions radfield.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -754,8 +754,8 @@ static auto planck_integral(double T_R, double nu_lower, double nu_upper, const
return integral;
}

auto planck_integral_analytic(const double T_R, const double nu_lower, const double nu_upper,
const bool times_nu) -> double {
auto planck_integral_analytic(const double T_R, const double nu_lower, const double nu_upper, const bool times_nu)
-> double {
// return the integral of nu^3 / (exp(h nu / k T) - 1) from nu_lower to nu_upper
// or if times_nu is true, the integral of nu^4 / (exp(h nu / k T) - 1) from nu_lower to nu_upper
double integral = 0.;
Expand Down Expand Up @@ -1043,7 +1043,7 @@ void normalise_J(const int modelgridindex, const double estimator_normfactor_ove
void normalise_bf_estimators(const int modelgridindex, const ptrdiff_t nonemptymgi,
const double estimator_normfactor_over_H) {
if constexpr (DETAILED_BF_ESTIMATORS_ON) {
printout("normalise_bf_estimators for cell %d with factor %g\n", modelgridindex, estimator_normfactor_over_H);
// printout("normalise_bf_estimators for cell %d with factor %g\n", modelgridindex, estimator_normfactor_over_H);
assert_always(nonemptymgi >= 0);
for (int i = 0; i < globals::bfestimcount; i++) {
const auto mgibfindex = nonemptymgi * globals::bfestimcount + i;
Expand Down Expand Up @@ -1138,11 +1138,26 @@ void reduce_estimators()
MPI_Allreduce(MPI_IN_PLACE, J.data(), nonempty_npts_model, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);
MPI_Allreduce(MPI_IN_PLACE, nuJ.data(), nonempty_npts_model, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD);

// if constexpr (DETAILED_BF_ESTIMATORS_ON) {
// for (ptrdiff_t nonemptymgi = 0; nonemptymgi < nonempty_npts_model; nonemptymgi++) {
// MPI_Allreduce(MPI_IN_PLACE, &bfrate_raw[nonemptymgi * globals::bfestimcount], globals::bfestimcount,
// MPI_DOUBLE,
// MPI_SUM, MPI_COMM_WORLD);
// }
// }

if constexpr (DETAILED_BF_ESTIMATORS_ON) {
for (ptrdiff_t nonemptymgi = 0; nonemptymgi < nonempty_npts_model; nonemptymgi++) {
MPI_Allreduce(MPI_IN_PLACE, &bfrate_raw[nonemptymgi * globals::bfestimcount], globals::bfestimcount, MPI_DOUBLE,
MPI_SUM, MPI_COMM_WORLD);
MPI_Allreduce(MPI_IN_PLACE,
&bfrate_raw[static_cast<size_t>(nonemptymgi) * static_cast<size_t>(globals::bfestimcount)],
globals::bfestimcount, MPI_DOUBLE, MPI_SUM, globals::mpi_comm_node);
}

if (globals::rank_in_node == 0) {
MPI_Allreduce(MPI_IN_PLACE, bfrate_raw, nonempty_npts_model * globals::bfestimcount, MPI_DOUBLE, MPI_SUM,
globals::mpi_comm_internode);
}
MPI_Bcast(bfrate_raw, nonempty_npts_model * globals::bfestimcount, MPI_DOUBLE, 0, globals::mpi_comm_node);
}

if constexpr (MULTIBIN_RADFIELD_MODEL_ON) {
Expand Down Expand Up @@ -1203,12 +1218,12 @@ void do_MPI_Bcast(const int modelgridindex, const int root, int root_node_id)
}
}

if constexpr (DETAILED_BF_ESTIMATORS_ON) {
if (globals::rank_in_node == 0) {
MPI_Bcast(&prev_bfrate_normed[nonemptymgi * globals::bfestimcount], globals::bfestimcount, MPI_FLOAT,
root_node_id, globals::mpi_comm_internode);
}
}
// if constexpr (DETAILED_BF_ESTIMATORS_ON) {
// if (globals::rank_in_node == 0) {
// MPI_Bcast(&prev_bfrate_normed[nonemptymgi * globals::bfestimcount], globals::bfestimcount, MPI_FLOAT,
// root_node_id, globals::mpi_comm_internode);
// }
// }

if constexpr (DETAILED_LINE_ESTIMATORS_ON) {
for (int jblueindex = 0; jblueindex < detailed_linecount; jblueindex++) {
Expand Down
23 changes: 20 additions & 3 deletions update_grid.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -1012,9 +1012,9 @@ void update_grid_cell(const int mgi, const int nts, const int nts_prev, const in
// full-spectrum and binned J and nuJ estimators
radfield::fit_parameters(mgi, nts);

if constexpr (DETAILED_BF_ESTIMATORS_ON) {
radfield::normalise_bf_estimators(mgi, nonemptymgi, estimator_normfactor / H);
}
// if constexpr (DETAILED_BF_ESTIMATORS_ON) {
// radfield::normalise_bf_estimators(mgi, nonemptymgi, estimator_normfactor / H);
// }

solve_Te_nltepops(mgi, nonemptymgi, nts, titer, heatingcoolingrates);
}
Expand Down Expand Up @@ -1147,6 +1147,23 @@ void update_grid(FILE *estimators_file, const int nts, const int nts_prev, const
#pragma omp for schedule(dynamic)
#endif

if constexpr (DETAILED_BF_ESTIMATORS_ON) {
if (!(nts == globals::timestep_initial && titer == 0)) {
for (int mgi = 0; mgi < grid::get_npts_model(); mgi++) {
const auto nonemptymgi = grid::get_modelcell_nonemptymgi(mgi);
if (!(globals::lte_iteration || grid::modelgrid[mgi].thick == 1)) {
int assoc_cells = grid::get_numassociatedcells(mgi);
if (assoc_cells > 0) {
double deltaV =
grid::get_modelcell_assocvolume_tmin(mgi) * pow(globals::timesteps[nts_prev].mid / globals::tmin, 3);
double estimator_normfactor = 1 / deltaV / deltat / globals::nprocs;
radfield::normalise_bf_estimators(mgi, nonemptymgi, estimator_normfactor / H);
}
}
}
}
}

for (int mgi = nstart; mgi < nstart + ndo; mgi++) {
/// Check if this task should work on the current model grid cell.
/// If yes, update the cell and write out the estimators
Expand Down

0 comments on commit 43d1bca

Please sign in to comment.