diff --git a/.github/workflows/ci-checks.yml b/.github/workflows/ci-checks.yml index 37bb5f256..c7058e23c 100644 --- a/.github/workflows/ci-checks.yml +++ b/.github/workflows/ci-checks.yml @@ -40,7 +40,7 @@ jobs: - name: Run clang-format style check uses: jidicula/clang-format-action@v4.11.0 with: - clang-format-version: '16' + clang-format-version: '18' check-path: . clang-tidy: diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 47cec82bd..d2fb353ff 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -31,6 +31,7 @@ jobs: kilonova_2d_2dgrid, kilonova_2d_3dgrid, nebularonezone_1d_3dgrid, + nebularonezone_1d_3dgrid_limitbfest, ] exclude: - os: self-hosted diff --git a/.idea/.gitignore b/.idea/.gitignore new file mode 100644 index 000000000..13566b81b --- /dev/null +++ b/.idea/.gitignore @@ -0,0 +1,8 @@ +# Default ignored files +/shelf/ +/workspace.xml +# Editor-based HTTP Client requests +/httpRequests/ +# Datasource local storage ignored files +/dataSources/ +/dataSources.local.xml diff --git a/.idea/codeStyles/Project.xml b/.idea/codeStyles/Project.xml new file mode 100644 index 000000000..f60388162 --- /dev/null +++ b/.idea/codeStyles/Project.xml @@ -0,0 +1,7 @@ + + + + + + \ No newline at end of file diff --git a/.idea/codeStyles/codeStyleConfig.xml b/.idea/codeStyles/codeStyleConfig.xml new file mode 100644 index 000000000..79ee123c2 --- /dev/null +++ b/.idea/codeStyles/codeStyleConfig.xml @@ -0,0 +1,5 @@ + + + + \ No newline at end of file diff --git a/.idea/misc.xml b/.idea/misc.xml new file mode 100644 index 000000000..53624c9e1 --- /dev/null +++ b/.idea/misc.xml @@ -0,0 +1,18 @@ + + + + + + + + \ No newline at end of file diff --git a/.idea/vcs.xml b/.idea/vcs.xml new file mode 100644 index 000000000..35eb1ddfb --- /dev/null +++ b/.idea/vcs.xml @@ -0,0 +1,6 @@ + + + + + + \ No newline at end of file diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 99fe3f88b..afd2ac7f4 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -27,7 +27,7 @@ repos: rev: v1.3.5 hooks: - id: clang-format - args: [-i] + args: [--version=18, -i] # - id: clang-tidy # args: [-extra-arg=-std=c++20] # - id: oclint diff --git a/artisoptions_christinenonthermal.h b/artisoptions_christinenonthermal.h index 9c5cf6ead..53d0b6f92 100644 --- a/artisoptions_christinenonthermal.h +++ b/artisoptions_christinenonthermal.h @@ -81,6 +81,8 @@ constexpr bool DETAILED_BF_ESTIMATORS_ON = true; constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP = 13; +constexpr bool LEVEL_HAS_BFEST(int element_z, int ionstage, int level) { return true; } + constexpr bool USE_LUT_PHOTOION = false; constexpr bool USE_LUT_BFHEATING = false; diff --git a/artisoptions_classic.h b/artisoptions_classic.h index 08eb7f280..4e395310f 100644 --- a/artisoptions_classic.h +++ b/artisoptions_classic.h @@ -75,6 +75,8 @@ constexpr bool DETAILED_LINE_ESTIMATORS_ON = false; constexpr bool DETAILED_BF_ESTIMATORS_ON = false; +constexpr bool LEVEL_HAS_BFEST(int element_z, int ionstage, int level) { return false; } + constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP = 13; constexpr bool USE_LUT_PHOTOION = true; diff --git a/artisoptions_kilonova_lte.h b/artisoptions_kilonova_lte.h index 02b8000fc..250824bad 100644 --- a/artisoptions_kilonova_lte.h +++ b/artisoptions_kilonova_lte.h @@ -74,6 +74,8 @@ constexpr bool DETAILED_LINE_ESTIMATORS_ON = false; constexpr bool DETAILED_BF_ESTIMATORS_ON = false; +constexpr bool LEVEL_HAS_BFEST(int element_z, int ionstage, int level) { return false; } + constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP = 13; constexpr bool USE_LUT_PHOTOION = true; diff --git a/artisoptions_nltenebular.h b/artisoptions_nltenebular.h index 4de2507d9..ca2704692 100644 --- a/artisoptions_nltenebular.h +++ b/artisoptions_nltenebular.h @@ -82,6 +82,12 @@ constexpr bool DETAILED_LINE_ESTIMATORS_ON = false; constexpr bool DETAILED_BF_ESTIMATORS_ON = true; +constexpr bool LEVEL_HAS_BFEST(int element_z, int ionstage, int level) { + // To only BF estimators for NLTE levels: + // return LEVEL_IS_NLTE(element_z, ionstage, level); + return true; +} + constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP = 13; constexpr bool USE_LUT_PHOTOION = false; @@ -146,4 +152,4 @@ constexpr bool KEEP_ALL_RESTART_FILES = false; constexpr bool BFCOOLING_USELEVELPOPNOTIONPOP = false; // NOLINTEND(modernize*,misc-unused-parameters) -#endif // ARTISOPTIONS_H \ No newline at end of file +#endif // ARTISOPTIONS_H diff --git a/artisoptions_nltewithoutnonthermal.h b/artisoptions_nltewithoutnonthermal.h index d42344dcb..f295413f4 100644 --- a/artisoptions_nltewithoutnonthermal.h +++ b/artisoptions_nltewithoutnonthermal.h @@ -78,6 +78,8 @@ constexpr bool DETAILED_LINE_ESTIMATORS_ON = false; constexpr bool DETAILED_BF_ESTIMATORS_ON = true; +constexpr bool LEVEL_HAS_BFEST(int element_z, int ionstage, int level) { return true; } + constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP = 13; constexpr bool USE_LUT_PHOTOION = false; diff --git a/atomic.cc b/atomic.cc index d08dd9378..a6f01bebb 100644 --- a/atomic.cc +++ b/atomic.cc @@ -21,8 +21,8 @@ int includedions = 0; int includedions_excludehighest = 0; std::array phixs_file_version_exists; -auto get_continuumindex_phixstargetindex(const int element, const int ion, const int level, const int phixstargetindex) - -> int +auto get_continuumindex_phixstargetindex(const int element, const int ion, const int level, + const int phixstargetindex) -> int /// Returns the index of the continuum associated to the given level. { return globals::elements[element].ions[ion].levels[level].cont_index - phixstargetindex; @@ -90,8 +90,8 @@ auto level_isinsuperlevel(const int element, const int ion, const int level) -> return (!is_nlte(element, ion, level) && level != 0 && (get_nlevels_nlte(element, ion) > 0)); } -auto photoionization_crosssection_fromtable(const float *const photoion_xs, const double nu_edge, const double nu) - -> double +auto photoionization_crosssection_fromtable(const float *const photoion_xs, const double nu_edge, + const double nu) -> double /// Calculates the photoionisation cross-section at frequency nu out of the atomic data. /// Input: - edge frequency nu_edge of the desired bf-continuum /// - nu diff --git a/decay.cc b/decay.cc index 25adb348d..b2f4e3a46 100644 --- a/decay.cc +++ b/decay.cc @@ -751,8 +751,8 @@ static auto sample_decaytime(const int decaypathindex, const double tdecaymin, c } static constexpr auto calculate_decaychain(const double firstinitabund, const std::vector &lambdas, - const int num_nuclides, const double timediff, const bool useexpansionfactor) - -> double { + const int num_nuclides, const double timediff, + const bool useexpansionfactor) -> double { // calculate final number abundance from multiple decays, e.g., Ni56 -> Co56 -> Fe56 (nuc[0] -> nuc[1] -> nuc[2]) // the top nuclide initial abundance is set and the chain-end abundance is returned (all intermediates nuclides // are assumed to start with zero abundance) @@ -912,8 +912,8 @@ static auto get_endecay_to_tinf_per_ejectamass_at_time(const int modelgridindex, } auto get_endecay_per_ejectamass_t0_to_time_withexpansion_chain_numerical(const int modelgridindex, - const int decaypathindex, const double tstart) - -> double + const int decaypathindex, + const double tstart) -> double // just here as as check on the analytic result from get_endecay_per_ejectamass_t0_to_time_withexpansion() // this version does an Euler integration { @@ -979,8 +979,8 @@ auto get_endecay_per_ejectamass_t0_to_time_withexpansion(const int modelgridinde return tot_endecay; } -static auto get_endecay_per_ejectamass_between_times(const int mgi, const int decaypathindex, double tlow, double thigh) - -> double +static auto get_endecay_per_ejectamass_between_times(const int mgi, const int decaypathindex, double tlow, + double thigh) -> double // get decay energy per mass [erg/g] released by a decaypath between times tlow [s] and thigh [s] { assert_always(tlow <= thigh); @@ -1016,8 +1016,8 @@ static auto get_simtime_endecay_per_ejectamass(const int mgi, const int decaypat return chainendecay; } -static auto get_decaypath_power_per_ejectamass(const int decaypathindex, const int modelgridindex, const double time) - -> double +static auto get_decaypath_power_per_ejectamass(const int decaypathindex, const int modelgridindex, + const double time) -> double // total decay power per mass [erg/s/g] for a given decaypath { // only decays at the end of the chain contributed from the initial abundance of the top of the chain are counted diff --git a/globals.cc b/globals.cc index 369193f5e..82664375a 100644 --- a/globals.cc +++ b/globals.cc @@ -12,6 +12,8 @@ std::unique_ptr timesteps = nullptr; double *rpkt_emiss = nullptr; /// Volume estimator for the rpkt emissivity +int bfestimcount = 0; + // for USE_LUT_PHOTOION = true double *corrphotoionrenorm = nullptr; double *gammaestimator = nullptr; diff --git a/globals.h b/globals.h index 81d2914fd..0352f5b76 100644 --- a/globals.h +++ b/globals.h @@ -50,6 +50,7 @@ struct fullphixslist { float *photoion_xs; double probability; int index_in_groundphixslist; + int bfestimindex; }; struct groundphixslist { @@ -209,6 +210,8 @@ extern std::unique_ptr timesteps; extern double *rpkt_emiss; +extern int bfestimcount; + // for USE_LUT_PHOTOION = true extern double *corrphotoionrenorm; extern double *gammaestimator; diff --git a/grid.cc b/grid.cc index d0db59a01..7f5d20627 100644 --- a/grid.cc +++ b/grid.cc @@ -2425,8 +2425,8 @@ static auto get_coordboundary_distances_cylindrical2d(std::span } [[nodiscard]] auto boundary_distance(std::span dir, std::span pos, - const double tstart, int cellindex, int *snext, enum cell_boundary *pkt_last_cross) - -> double + const double tstart, int cellindex, int *snext, + enum cell_boundary *pkt_last_cross) -> double /// Basic routine to compute distance to a cell boundary. { if constexpr (FORCE_SPHERICAL_ESCAPE_SURFACE) { diff --git a/input.cc b/input.cc index f4f19ac15..e886db7a6 100644 --- a/input.cc +++ b/input.cc @@ -1118,8 +1118,8 @@ static void read_atomicdata_files() { printout("cont_index %d\n", cont_index); } -static auto search_groundphixslist(double nu_edge, int *index_in_groundlevelcontestimator, int el, int in, int ll) - -> int +static auto search_groundphixslist(double nu_edge, int *index_in_groundlevelcontestimator, int el, int in, + int ll) -> int /// Return the closest ground level continuum index to the given edge /// frequency. If the given edge frequency is redder than the reddest /// continuum return -1. @@ -1456,6 +1456,7 @@ static void setup_phixs_list() { globals::nbfcontinua * (sizeof(fullphixslist)) / 1024. / 1024.); size_t nbftables = 0; int allcontindex = 0; + globals::bfestimcount = 0; for (int element = 0; element < get_nelements(); element++) { const int nions = get_nions(element); for (int ion = 0; ion < nions - 1; ion++) { @@ -1482,6 +1483,13 @@ static void setup_phixs_list() { nonconstallcont[allcontindex].probability = get_phixsprobability(element, ion, level, phixstargetindex); nonconstallcont[allcontindex].upperlevel = get_phixsupperlevel(element, ion, level, phixstargetindex); + if (LEVEL_HAS_BFEST(get_atomicnumber(element), get_ionstage(element, ion), level)) { + nonconstallcont[allcontindex].bfestimindex = globals::bfestimcount; + globals::bfestimcount++; + } else { + nonconstallcont[allcontindex].bfestimindex = -1; + } + if constexpr (USE_LUT_PHOTOION || USE_LUT_BFHEATING) { int index_in_groundlevelcontestimator = 0; nonconstallcont[allcontindex].index_in_groundphixslist = @@ -1495,6 +1503,7 @@ static void setup_phixs_list() { } } } + printout("[info] BF estimators activated for %d photoionisation transitions\n", globals::bfestimcount); assert_always(allcontindex == globals::nbfcontinua); assert_always(globals::nbfcontinua >= 0); // was initialised as -1 before startup @@ -2189,4 +2198,4 @@ void write_timestep_file() { globals::timesteps[n].width / DAY); } fclose(timestepfile); -} \ No newline at end of file +} diff --git a/ltepop.cc b/ltepop.cc index 376207534..e811e5550 100644 --- a/ltepop.cc +++ b/ltepop.cc @@ -241,8 +241,8 @@ auto calculate_levelpop_lte(int modelgridindex, int element, int ion, int level) exp(-E_aboveground / KB / T_exc)); } -static auto calculate_levelpop_nominpop(int modelgridindex, int element, int ion, int level, bool *skipminpop) - -> double { +static auto calculate_levelpop_nominpop(int modelgridindex, int element, int ion, int level, + bool *skipminpop) -> double { assert_testmodeonly(modelgridindex < grid::get_npts_model()); assert_testmodeonly(element < get_nelements()); assert_testmodeonly(ion < get_nions(element)); @@ -424,8 +424,8 @@ auto get_nnion(int modelgridindex, int element, int ion) -> double grid::modelgrid[modelgridindex].composition[element].partfunct[ion] / stat_weight(element, ion, 0); } -static auto find_uppermost_ion(const int modelgridindex, const int element, const double nne_hi, const bool force_lte) - -> int { +static auto find_uppermost_ion(const int modelgridindex, const int element, const double nne_hi, + const bool force_lte) -> int { const int nions = get_nions(element); if (nions == 0) { return -1; diff --git a/ltepop.h b/ltepop.h index cddec48fa..d48f165d3 100644 --- a/ltepop.h +++ b/ltepop.h @@ -11,13 +11,13 @@ auto get_groundlevelpop(int modelgridindex, int element, int ion) -> double; auto calculate_levelpop(int modelgridindex, int element, int ion, int level) -> double; auto calculate_levelpop_lte(int modelgridindex, int element, int ion, int level) -> double; auto get_levelpop(int modelgridindex, int element, int ion, int level) -> double; -[[nodiscard]] auto calculate_sahafact(int element, int ion, int level, int upperionlevel, double T, double E_threshold) - -> double; +[[nodiscard]] auto calculate_sahafact(int element, int ion, int level, int upperionlevel, double T, + double E_threshold) -> double; [[nodiscard]] auto get_nnion(int modelgridindex, int element, int ion) -> double; void calculate_ion_balance_nne(int modelgridindex); void calculate_cellpartfuncts(int modelgridindex, int element); -[[nodiscard]] auto calculate_ionfractions(int element, int modelgridindex, double nne, bool use_phi_lte) - -> std::vector; +[[nodiscard]] auto calculate_ionfractions(int element, int modelgridindex, double nne, + bool use_phi_lte) -> std::vector; void set_groundlevelpops(int modelgridindex, int element, float nne, bool force_lte); #endif // LTEPOP_H diff --git a/macroatom.cc b/macroatom.cc index b0c52af7a..7cc647cc0 100644 --- a/macroatom.cc +++ b/macroatom.cc @@ -744,8 +744,8 @@ auto rad_deexcitation_ratecoeff(const int modelgridindex, const int element, con } auto rad_excitation_ratecoeff(const int modelgridindex, const int element, const int ion, const int lower, - const int uptransindex, const double epsilon_trans, int lineindex, const double t_current) - -> double + const int uptransindex, const double epsilon_trans, int lineindex, + const double t_current) -> double /// radiative excitation rate: paperII 3.5.2 // multiply by lower level population to get a rate per second { diff --git a/macroatom.h b/macroatom.h index b3df8b56b..f43b788ac 100644 --- a/macroatom.h +++ b/macroatom.h @@ -40,8 +40,8 @@ auto rad_excitation_ratecoeff(int modelgridindex, int element, int ion, int lowe double epsilon_trans, int lineindex, double t_current) -> double; auto rad_recombination_ratecoeff(float T_e, float nne, int element, int upperion, int upperionlevel, int lowerionlevel, int modelgridindex) -> double; -auto stim_recombination_ratecoeff(float nne, int element, int upperion, int upper, int lower, int modelgridindex) - -> double; +auto stim_recombination_ratecoeff(float nne, int element, int upperion, int upper, int lower, + int modelgridindex) -> double; auto col_recombination_ratecoeff(int modelgridindex, int element, int upperion, int upper, int lower, double epsilon_trans) -> double; diff --git a/nltepop.cc b/nltepop.cc index 76f2d038f..a4fd570c4 100644 --- a/nltepop.cc +++ b/nltepop.cc @@ -119,8 +119,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); @@ -174,8 +174,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); } @@ -393,7 +393,7 @@ static void nltepop_reset_element(const int modelgridindex, const int element) { for (int ion = 0; ion < nions; ion++) { const int nlte_start = globals::elements[element].ions[ion].first_nlte; const int nlevels_nlte = get_nlevels_nlte(element, ion); - for (int level = 1; level < nlevels_nlte; level++) { + for (int level = 1; level <= nlevels_nlte; level++) { grid::modelgrid[modelgridindex].nlte_pops[nlte_start + level - 1] = -1.0; // flag to indicate no useful data } diff --git a/nonthermal.cc b/nonthermal.cc index 907f5b424..28b57c94f 100644 --- a/nonthermal.cc +++ b/nonthermal.cc @@ -1551,8 +1551,8 @@ auto nt_ionisation_maxupperion(const int element, const int lowerion) -> int { return maxupper; } -auto nt_random_upperion(const int modelgridindex, const int element, const int lowerion, const bool energyweighted) - -> int { +auto nt_random_upperion(const int modelgridindex, const int element, const int lowerion, + const bool energyweighted) -> int { assert_testmodeonly(lowerion < get_nions(element) - 1); if (NT_SOLVE_SPENCERFANO && NT_MAX_AUGER_ELECTRONS > 0) { while (true) { @@ -1610,8 +1610,8 @@ auto nt_ionization_ratecoeff(const int modelgridindex, const int element, const static auto calculate_nt_excitation_ratecoeff_perdeposition(const int modelgridindex, const int element, const int ion, const int lower, const int uptransindex, - const double statweight_lower, const double epsilon_trans) - -> double + const double statweight_lower, + const double epsilon_trans) -> double // Kozma & Fransson equation 9 divided by level population and epsilon_trans { if (nt_solution[modelgridindex].yfunc == nullptr) { @@ -1962,8 +1962,8 @@ static void analyse_sf_solution(const int modelgridindex, const int timestep, co (excitationindex)++; } } // NT_EXCITATION_ON - } // for t - } // for lower + } // for t + } // for lower printout(" frac_excitation: %g\n", frac_excitation_ion); if (frac_excitation_ion > 1. || !std::isfinite(frac_excitation_ion)) { diff --git a/radfield.cc b/radfield.cc index dfb5749bb..29f89eb90 100644 --- a/radfield.cc +++ b/radfield.cc @@ -7,9 +7,11 @@ #include #include +#include #include #include "atomic.h" +#include "globals.h" #include "grid.h" #include "sn3d.h" #include "vectors.h" @@ -354,7 +356,7 @@ void init(int my_rank, int ndo_nonempty) if (globals::rank_in_node == 0) { my_rank_cells += nonempty_npts_model - (my_rank_cells * globals::node_nprocs); } - auto size = static_cast(my_rank_cells * globals::nbfcontinua * sizeof(float)); + auto size = static_cast(my_rank_cells * globals::bfestimcount * sizeof(float)); int disp_unit = sizeof(float); MPI_Win_allocate_shared(size, disp_unit, MPI_INFO_NULL, globals::mpi_comm_node, &prev_bfrate_normed, &win_prev_bfrate_normed); @@ -362,16 +364,16 @@ void init(int my_rank, int ndo_nonempty) } #else { - prev_bfrate_normed = static_cast(malloc(nonempty_npts_model * globals::nbfcontinua * sizeof(float))); + prev_bfrate_normed = static_cast(malloc(nonempty_npts_model * globals::bfestimcount * sizeof(float))); } #endif printout("[info] mem_usage: detailed bf estimators for non-empty cells occupy %.3f MB (node shared memory)\n", - nonempty_npts_model * globals::nbfcontinua * sizeof(float) / 1024. / 1024.); + nonempty_npts_model * globals::bfestimcount * sizeof(float) / 1024. / 1024.); - bfrate_raw = static_cast(malloc(nonempty_npts_model * globals::nbfcontinua * sizeof(double))); + bfrate_raw = static_cast(malloc(nonempty_npts_model * globals::bfestimcount * sizeof(double))); printout("[info] mem_usage: detailed bf estimator acculumators for non-empty cells occupy %.3f MB\n", - nonempty_npts_model * globals::nbfcontinua * sizeof(double) / 1024. / 1024.); + nonempty_npts_model * globals::bfestimcount * sizeof(double) / 1024. / 1024.); } for (int modelgridindex = 0; modelgridindex < grid::get_npts_model(); modelgridindex++) { @@ -380,9 +382,9 @@ void init(int my_rank, int ndo_nonempty) if constexpr (MULTIBIN_RADFIELD_MODEL_ON) { if (globals::rank_in_node == 0) { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); for (int binindex = 0; binindex < RADFIELDBINCOUNT; binindex++) { - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; radfieldbin_solutions[mgibinindex].W = -1.; radfieldbin_solutions[mgibinindex].T_R = -1.; } @@ -473,23 +475,21 @@ auto get_Jb_lu_contribcount(const int modelgridindex, const int jblueindex) -> i static auto get_bin_J(int modelgridindex, int binindex) -> double // get the normalised J_nu { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); assert_testmodeonly(J_normfactor[nonemptymgi] > 0.0); assert_testmodeonly(modelgridindex < grid::get_npts_model()); assert_testmodeonly(binindex >= 0); assert_testmodeonly(binindex < RADFIELDBINCOUNT); - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; - return radfieldbins[mgibinindex].J_raw * J_normfactor[nonemptymgi]; + return radfieldbins[nonemptymgi * RADFIELDBINCOUNT + binindex].J_raw * J_normfactor[nonemptymgi]; } static auto get_bin_nuJ(int modelgridindex, int binindex) -> double { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); assert_testmodeonly(J_normfactor[nonemptymgi] > 0.0); assert_testmodeonly(modelgridindex < grid::get_npts_model()); assert_testmodeonly(binindex >= 0); assert_testmodeonly(binindex < RADFIELDBINCOUNT); - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; - return radfieldbins[mgibinindex].nuJ_raw * J_normfactor[nonemptymgi]; + return radfieldbins[nonemptymgi * RADFIELDBINCOUNT + binindex].nuJ_raw * J_normfactor[nonemptymgi]; } static inline auto get_bin_nu_bar(int modelgridindex, int binindex) -> double @@ -508,18 +508,18 @@ static inline auto get_bin_nu_lower(int binindex) -> double { } static inline auto get_bin_contribcount(int modelgridindex, int binindex) -> int { - const int mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; - return radfieldbins[mgibinindex].contribcount; + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + return radfieldbins[nonemptymgi * RADFIELDBINCOUNT + binindex].contribcount; } static inline auto get_bin_W(int modelgridindex, int binindex) -> float { - const int mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; - return radfieldbin_solutions[mgibinindex].W; + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + return radfieldbin_solutions[nonemptymgi * RADFIELDBINCOUNT + binindex].W; } static inline auto get_bin_T_R(int modelgridindex, int binindex) -> float { - const int mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; - return radfieldbin_solutions[mgibinindex].T_R; + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + return radfieldbin_solutions[nonemptymgi * RADFIELDBINCOUNT + binindex].T_R; } static inline auto select_bin(double nu) -> int { @@ -540,7 +540,7 @@ static inline auto select_bin(double nu) -> int { void write_to_file(int modelgridindex, int timestep) { assert_always(MULTIBIN_RADFIELD_MODEL_ON); - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); #ifdef _OPENMP #pragma omp critical(out_file) { @@ -650,13 +650,13 @@ void zero_estimators(int modelgridindex) return; } - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); if constexpr (DETAILED_BF_ESTIMATORS_ON) { assert_always(bfrate_raw != nullptr); if (grid::get_numassociatedcells(modelgridindex) > 0) { - for (int i = 0; i < globals::nbfcontinua; i++) { - bfrate_raw[nonemptymgi * globals::nbfcontinua + i] = 0.; + for (int i = 0; i < globals::bfestimcount; i++) { + bfrate_raw[nonemptymgi * globals::bfestimcount + i] = 0.; } } } @@ -678,7 +678,7 @@ void zero_estimators(int modelgridindex) assert_always(radfieldbins != nullptr); for (int binindex = 0; binindex < RADFIELDBINCOUNT; binindex++) { - const int mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; radfieldbins[mgibinindex].J_raw = 0.0; radfieldbins[mgibinindex].nuJ_raw = 0.0; radfieldbins[mgibinindex].contribcount = 0; @@ -708,8 +708,13 @@ static void update_bfestimators(const int nonemptymgi, const double distance_e_c const double nu_max_phixs = nu_edge * last_phixs_nuovernuedge; // nu of the uppermost point in the phixs table if (nu_cmf >= nu_edge && nu_cmf <= nu_max_phixs) { - safeadd(bfrate_raw[nonemptymgi * nbfcontinua + allcontindex], - globals::phixslist[tid].gamma_contr[allcontindex] * distance_e_cmf_over_nu); + int bf_estimator = globals::allcont[allcontindex].bfestimindex; + if (bf_estimator >= 0) { + assert_testmodeonly(bf_estimator < globals::bfestimcount); + + safeadd(bfrate_raw[nonemptymgi * globals::bfestimcount + bf_estimator], + globals::phixslist[tid].gamma_contr[allcontindex] * distance_e_cmf_over_nu); + } } else if (nu_cmf < nu_edge) { // list is sorted by nu_edge, so all remaining will have nu_cmf < nu_edge @@ -731,7 +736,7 @@ void update_estimators(const int nonemptymgi, const double distance_e_cmf, const const int binindex = select_bin(nu_cmf); if (binindex >= 0) { - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; safeadd(radfieldbins[mgibinindex].J_raw, distance_e_cmf); safeadd(radfieldbins[mgibinindex].nuJ_raw, distance_e_cmf * nu_cmf); safeincrement(radfieldbins[mgibinindex].contribcount); @@ -769,7 +774,7 @@ auto radfield(double nu, int modelgridindex) -> double if (globals::timestep >= FIRST_NLTE_RADFIELD_TIMESTEP) { const int binindex = select_bin(nu); if (binindex >= 0) { - const int mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = grid::get_modelcell_nonemptymgi(modelgridindex) * RADFIELDBINCOUNT + binindex; const struct radfieldbin_solution *const bin = &radfieldbin_solutions[mgibinindex]; if (bin->W >= 0.) { const double J_nu = dbb(nu, bin->T_R, bin->W); @@ -1114,7 +1119,7 @@ void fit_parameters(int modelgridindex, int timestep) // T_R_bin = -1; // W_bin = -1; // } - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; radfieldbin_solutions[mgibinindex].T_R = T_R_bin; radfieldbin_solutions[mgibinindex].W = W_bin; } @@ -1159,11 +1164,11 @@ void normalise_J(const int modelgridindex, const double estimator_normfactor_ove void normalise_bf_estimators(const int modelgridindex, 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); - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); assert_always(nonemptymgi >= 0); - for (int i = 0; i < globals::nbfcontinua; i++) { - const int mgibfindex = nonemptymgi * globals::nbfcontinua + i; - prev_bfrate_normed[mgibfindex] = bfrate_raw[mgibfindex] * estimator_normfactor_over_H; + for (int i = 0; i < globals::bfestimcount; i++) { + const auto detailed_mgibfindex = nonemptymgi * globals::bfestimcount + i; + prev_bfrate_normed[detailed_mgibfindex] = bfrate_raw[detailed_mgibfindex] * estimator_normfactor_over_H; } } } @@ -1190,10 +1195,11 @@ auto get_bfrate_estimator(const int element, const int lowerion, const int lower if constexpr (!DETAILED_BF_ESTIMATORS_ON) { return -1; } else { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); const int allcontindex = get_bfcontindex(element, lowerion, lower, phixstargetindex); if (allcontindex >= 0) { - return prev_bfrate_normed[nonemptymgi * globals::nbfcontinua + allcontindex]; + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const int bfestimindex = globals::allcont[allcontindex].bfestimindex; + return (bfestimindex >= 0) ? prev_bfrate_normed[nonemptymgi * globals::bfestimcount + bfestimindex] : 0.; } printout("no bf rate for element Z=%d ion_stage %d lower %d phixstargetindex %d\n", get_atomicnumber(element), @@ -1257,8 +1263,8 @@ void reduce_estimators() MPI_Allreduce(MPI_IN_PLACE, nuJ.data(), nonempty_npts_model, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); if constexpr (DETAILED_BF_ESTIMATORS_ON) { - MPI_Allreduce(MPI_IN_PLACE, bfrate_raw, nonempty_npts_model * globals::nbfcontinua, MPI_DOUBLE, MPI_SUM, - MPI_COMM_WORLD); + MPI_Allreduce(MPI_IN_PLACE, bfrate_raw, grid::get_nonempty_npts_model() * globals::bfestimcount, MPI_DOUBLE, + MPI_SUM, MPI_COMM_WORLD); } if constexpr (MULTIBIN_RADFIELD_MODEL_ON) { @@ -1268,7 +1274,7 @@ void reduce_estimators() for (int nonemptymgi = 0; nonemptymgi < nonempty_npts_model; nonemptymgi++) { for (int binindex = 0; binindex < RADFIELDBINCOUNT; binindex++) { - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; MPI_Allreduce(MPI_IN_PLACE, &radfieldbins[mgibinindex].J_raw, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); MPI_Allreduce(MPI_IN_PLACE, &radfieldbins[mgibinindex].nuJ_raw, 1, MPI_DOUBLE, MPI_SUM, MPI_COMM_WORLD); MPI_Allreduce(MPI_IN_PLACE, &radfieldbins[mgibinindex].contribcount, 1, MPI_INT, MPI_SUM, MPI_COMM_WORLD); @@ -1306,12 +1312,12 @@ void do_MPI_Bcast(const int modelgridindex, const int root, int root_node_id) return; } - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); MPI_Bcast(&J_normfactor[nonemptymgi], 1, MPI_DOUBLE, root, MPI_COMM_WORLD); if constexpr (MULTIBIN_RADFIELD_MODEL_ON) { for (int binindex = 0; binindex < RADFIELDBINCOUNT; binindex++) { - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; if (globals::rank_in_node == 0) { MPI_Bcast(&radfieldbin_solutions[mgibinindex].W, 1, MPI_FLOAT, root_node_id, globals::mpi_comm_internode); MPI_Bcast(&radfieldbin_solutions[mgibinindex].T_R, 1, MPI_FLOAT, root_node_id, globals::mpi_comm_internode); @@ -1324,8 +1330,8 @@ 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::nbfcontinua], globals::nbfcontinua, MPI_FLOAT, root_node_id, - globals::mpi_comm_internode); + MPI_Bcast(&prev_bfrate_normed[nonemptymgi * globals::bfestimcount], globals::bfestimcount, MPI_FLOAT, + root_node_id, globals::mpi_comm_internode); } } @@ -1358,12 +1364,15 @@ void write_restart_data(FILE *gridsave_file) { const int nbfcontinua = globals::nbfcontinua; fprintf(gridsave_file, "%d\n", nbfcontinua); + const int bfestimcount = globals::bfestimcount; + fprintf(gridsave_file, "%d\n", bfestimcount); + for (int modelgridindex = 0; modelgridindex < grid::get_npts_model(); modelgridindex++) { if (grid::get_numassociatedcells(modelgridindex) > 0) { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); fprintf(gridsave_file, "%d\n", modelgridindex); - for (int i = 0; i < nbfcontinua; i++) { - fprintf(gridsave_file, "%a ", prev_bfrate_normed[nonemptymgi * nbfcontinua + i]); + for (int i = 0; i < bfestimcount; i++) { + fprintf(gridsave_file, "%a ", prev_bfrate_normed[nonemptymgi * bfestimcount + i]); } } } @@ -1379,13 +1388,13 @@ void write_restart_data(FILE *gridsave_file) { for (int modelgridindex = 0; modelgridindex < grid::get_npts_model(); modelgridindex++) { if (grid::get_numassociatedcells(modelgridindex) > 0) { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); assert_testmodeonly(nonemptymgi >= 0); fprintf(gridsave_file, "%d %la\n", modelgridindex, J_normfactor[nonemptymgi]); if constexpr (MULTIBIN_RADFIELD_MODEL_ON) { for (int binindex = 0; binindex < RADFIELDBINCOUNT; binindex++) { - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; fprintf(gridsave_file, "%la %la %a %a %d\n", radfieldbins[mgibinindex].J_raw, radfieldbins[mgibinindex].nuJ_raw, radfieldbin_solutions[mgibinindex].W, radfieldbin_solutions[mgibinindex].T_R, radfieldbins[mgibinindex].contribcount); @@ -1456,17 +1465,21 @@ void read_restart_data(FILE *gridsave_file) { assert_always(fscanf(gridsave_file, "%d\n", &gridsave_nbf_in) == 1); assert_always(gridsave_nbf_in == globals::nbfcontinua); + int gridsave_nbfestim_in = 0; + assert_always(fscanf(gridsave_file, "%d\n", &gridsave_nbfestim_in) == 1); + assert_always(gridsave_nbfestim_in == globals::bfestimcount); + for (int modelgridindex = 0; modelgridindex < grid::get_npts_model(); modelgridindex++) { if (grid::get_numassociatedcells(modelgridindex) > 0) { - const int nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); + const ptrdiff_t nonemptymgi = grid::get_modelcell_nonemptymgi(modelgridindex); int mgi_in = 0; assert_always(fscanf(gridsave_file, "%d\n", &mgi_in) == 1); assert_always(mgi_in == modelgridindex); - for (int i = 0; i < globals::nbfcontinua; i++) { + for (int i = 0; i < globals::bfestimcount; i++) { float bfrate_normed = 0; assert_always(fscanf(gridsave_file, "%a ", &bfrate_normed) == 1); - const int mgibfindex = nonemptymgi * globals::nbfcontinua + i; + const auto mgibfindex = nonemptymgi * globals::bfestimcount + i; #ifdef MPI_ON if (globals::rank_in_node == 0) #endif @@ -1505,7 +1518,7 @@ void read_restart_data(FILE *gridsave_file) { if constexpr (MULTIBIN_RADFIELD_MODEL_ON) { for (int binindex = 0; binindex < RADFIELDBINCOUNT; binindex++) { - const int mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; + const ptrdiff_t mgibinindex = nonemptymgi * RADFIELDBINCOUNT + binindex; float W = 0; float T_R = 0; assert_always(fscanf(gridsave_file, "%la %la %a %a %d\n", &radfieldbins[mgibinindex].J_raw, diff --git a/ratecoeff.cc b/ratecoeff.cc index bbdf776d8..c7420387c 100644 --- a/ratecoeff.cc +++ b/ratecoeff.cc @@ -709,8 +709,8 @@ auto get_spontrecombcoeff(int element, int ion, int level, int phixstargetindex, auto calculate_ionrecombcoeff(const int modelgridindex, const float T_e, const int element, const int upperion, const bool assume_lte, const bool collisional_not_radiative, const bool printdebug, - const bool lower_superlevel_only, const bool per_groundmultipletpop, const bool stimonly) - -> double + const bool lower_superlevel_only, const bool per_groundmultipletpop, + const bool stimonly) -> double // multiply by upper ion population (or ground population if per_groundmultipletpop is true) and nne to get a rate { const int lowerion = upperion - 1; diff --git a/rpkt.h b/rpkt.h index 55f5ea3e3..e32dc23a9 100644 --- a/rpkt.h +++ b/rpkt.h @@ -33,8 +33,8 @@ constexpr auto get_linedistance(const double prop_time, const double nu_cmf, con return CLIGHT * prop_time * (nu_cmf / nu_trans - 1); } -[[nodiscard]] inline auto get_ionestimindex_nonemptymgi(const int nonemptymgi, const int element, const int ion) - -> int { +[[nodiscard]] inline auto get_ionestimindex_nonemptymgi(const int nonemptymgi, const int element, + const int ion) -> int { assert_testmodeonly(ion >= 0); assert_testmodeonly(ion < get_nions(element) - 1); return nonemptymgi * get_includedions() + get_uniqueionindex(element, ion); diff --git a/stats.cc b/stats.cc index f584e2ec8..8ba380327 100644 --- a/stats.cc +++ b/stats.cc @@ -110,8 +110,8 @@ void increment_ion_stats_contabsorption(const struct packet *const pkt_ptr, cons } } -auto get_ion_stats(const int modelgridindex, const int element, const int ion, enum ionstattypes ionstattype) - -> double { +auto get_ion_stats(const int modelgridindex, const int element, const int ion, + enum ionstattypes ionstattype) -> double { assert_always(ion < get_nions(element)); assert_always(ionstattype < ION_STAT_COUNT); const int uniqueionindex = get_uniqueionindex(element, ion); diff --git a/tests/nebularonezone_1d_3dgrid_limitbfest_inputfiles/results_md5_final.txt b/tests/nebularonezone_1d_3dgrid_limitbfest_inputfiles/results_md5_final.txt new file mode 100644 index 000000000..9a425378d --- /dev/null +++ b/tests/nebularonezone_1d_3dgrid_limitbfest_inputfiles/results_md5_final.txt @@ -0,0 +1,20 @@ +77a8cefec72da8c0bcebdf4e04d2fcf5 absorption.out +aff44ecad68d986f04fcd6110331ba8b bflist.out +9cfd74fab2244cafd4bce62dc6fe2ad7 deposition.out +986056dc9a141a00e7b62e653b2b4f4b emission.out +8607d1185f06e29637171e28955097e1 emissiontrue.out +9aacf70462a00d805820041d07deaa02 gamma_light_curve.out +035228757b6e30f01ff6a04a6a94d364 gamma_spec.out +057b226c371f3819cba5e04bfea3d114 gammalinelist.out +47dd2e31cc98d5d1dda67520fb399429 grid.out +2c51d13a376994eb179c6e97d44f1dfd light_curve.out +9773becefdfe4afffb041b703210c67c linestat.out +98195ccb920b6831d137861c8c6cbfa7 modelgridrankassignments.out +ce06681bb20c3a895bc724caa303cb54 packets00_0000.out +195ea92baa7650d5ae6cab37393a1483 packets00_0001.out +7ff320b87eb952c5e24295c2e12a6258 spec.out +8e7163982f1aa938dc1505478b8c60d1 timesteps.out +463e7d1a82b94f8c0fa7083c1f9b59d7 job1/estimators_0000.out +ebc4df506dc9890bde3a762f581bd5ad job1/nlte_0000.out +1c4b1ccb2ea3bd56f3ef976e4dce870e job1/nonthermalspec_0000.out +ef9c7843dfc4b7dc5e6908270a6bd780 job1/radfield_0000.out \ No newline at end of file diff --git a/tests/nebularonezone_1d_3dgrid_limitbfest_inputfiles/results_md5_job0.txt b/tests/nebularonezone_1d_3dgrid_limitbfest_inputfiles/results_md5_job0.txt new file mode 100644 index 000000000..fe98f1ce9 --- /dev/null +++ b/tests/nebularonezone_1d_3dgrid_limitbfest_inputfiles/results_md5_job0.txt @@ -0,0 +1,19 @@ +f5e5d8fc1f4eb7d1d2e7624e6929954c absorption.out +aff44ecad68d986f04fcd6110331ba8b bflist.out +b1c057449d45d3eafc42a30c2530cb92 deposition.out +91195eb22d4b77ed93bfc5edc96fbd79 emission.out +6bdc6267ae27feab57a34c8a91f78c37 emissiontrue.out +59e10adcf5e9d1c05f935e5583ad9e1e gamma_light_curve.out +057b226c371f3819cba5e04bfea3d114 gammalinelist.out +47dd2e31cc98d5d1dda67520fb399429 grid.out +e0d00b6a32be200c5cf0e9962089c0f3 light_curve.out +9773becefdfe4afffb041b703210c67c linestat.out +98195ccb920b6831d137861c8c6cbfa7 modelgridrankassignments.out +353ee551590fbe95446d10f2554adb9c packets00_0000.out +58ae1da9bf04c21767d40fd722bf737e packets00_0001.out +65539a0900e70d998038852567575617 spec.out +8e7163982f1aa938dc1505478b8c60d1 timesteps.out +2a274951a4fff75b7962281cf7bd043c job0/estimators_0000.out +7d86d86b36f9208ba00af17f31846025 job0/nlte_0000.out +1c4b1ccb2ea3bd56f3ef976e4dce870e job0/nonthermalspec_0000.out +fc10e7d1fd915e9372566c536d9d4652 job0/radfield_0000.out \ No newline at end of file diff --git a/tests/setup_nebularonezone_1d_3dgrid_limitbfest.sh b/tests/setup_nebularonezone_1d_3dgrid_limitbfest.sh new file mode 100755 index 000000000..a2b0c6083 --- /dev/null +++ b/tests/setup_nebularonezone_1d_3dgrid_limitbfest.sh @@ -0,0 +1,45 @@ +#!/usr/bin/env zsh + +set -x + +runfolder=nebularonezone_1d_3dgrid_limitbfest_testrun + +rsync -av nebularonezone_1d_3dgrid_inputfiles/ nebularonezone_1d_3dgrid_limitbfest_testrun/ + +rsync -av nebularonezone_1d_3dgrid_limitbfest_inputfiles/ nebularonezone_1d_3dgrid_limitbfest_testrun/ + +if [ ! -f atomicdata_feconi.tar.xz ]; then curl -O https://theory.gsi.de/~lshingle/artis_http_public/artis/atomicdata_feconi.tar.xz; fi + +tar -xf atomicdata_feconi.tar.xz --directory nebularonezone_1d_3dgrid_limitbfest_testrun/ + +cp ../data/* nebularonezone_1d_3dgrid_limitbfest_testrun/ + +cp ../artisoptions_nltenebular.h nebularonezone_1d_3dgrid_limitbfest_testrun/artisoptions.h + +cd nebularonezone_1d_3dgrid_limitbfest_testrun + +sed -i'' -e 's/constexpr int MPKTS.*/constexpr int MPKTS = 1000000;/g' artisoptions.h + +sed -i'' -e 's/constexpr int GRID_TYPE.*/constexpr int GRID_TYPE = GRID_CARTESIAN3D;/g' artisoptions.h + +sed -i'' -e 's/constexpr int CUBOID_NCOORDGRID_X.*/constexpr int CUBOID_NCOORDGRID_X = 50;/g' artisoptions.h +sed -i'' -e 's/constexpr int CUBOID_NCOORDGRID_Y.*/constexpr int CUBOID_NCOORDGRID_Y = 50;/g' artisoptions.h +sed -i'' -e 's/constexpr int CUBOID_NCOORDGRID_Z.*/constexpr int CUBOID_NCOORDGRID_Z = 50;/g' artisoptions.h + +sed -i'' -e 's/constexpr int TABLESIZE.*/constexpr int TABLESIZE = 20;/g' artisoptions.h +sed -i'' -e 's/constexpr double MINTEMP.*/constexpr double MINTEMP = 2000.;/g' artisoptions.h +sed -i'' -e 's/constexpr double MAXTEMP.*/constexpr double MAXTEMP = 10000.;/g' artisoptions.h + +sed -i'' -e 's/constexpr int FIRST_NLTE_RADFIELD_TIMESTEP.*/constexpr int FIRST_NLTE_RADFIELD_TIMESTEP = 7;/g' artisoptions.h + +sed -i'' -e 's/constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP.*/constexpr int DETAILED_BF_ESTIMATORS_USEFROMTIMESTEP = 7;/g' artisoptions.h + +sed -i'' -e 's/constexpr bool SF_AUGER_CONTRIBUTION_ON.*/constexpr bool SF_AUGER_CONTRIBUTION_ON = false;/g' artisoptions.h + +sed -i'' -e 's/constexpr bool WRITE_PARTIAL_EMISSIONABSORPTIONSPEC.*/constexpr bool WRITE_PARTIAL_EMISSIONABSORPTIONSPEC = true;/g' artisoptions.h + +sed -i'' -e 's/ \/\/ return LEVEL_IS_NLTE.*/return LEVEL_IS_NLTE(element_z, ionstage, level);/g' artisoptions.h + +cd - + +set +x diff --git a/thermalbalance.cc b/thermalbalance.cc index 41518c48d..d733285fe 100644 --- a/thermalbalance.cc +++ b/thermalbalance.cc @@ -76,8 +76,8 @@ static auto integrand_bfheatingcoeff_custom_radfield(double nu, void *voidparas) return sigma_bf * (1 - nu_edge / nu) * radfield::radfield(nu, modelgridindex) * (1 - exp(-HOVERKB * nu / T_R)); } -static auto calculate_bfheatingcoeff(int element, int ion, int level, int phixstargetindex, int modelgridindex) - -> double { +static auto calculate_bfheatingcoeff(int element, int ion, int level, int phixstargetindex, + int modelgridindex) -> double { double error = 0.0; const double epsrel = 1e-3; const double epsrelwarning = 1e-1; diff --git a/vpkt.cc b/vpkt.cc index 0f5bb28fb..fedde4168 100644 --- a/vpkt.cc +++ b/vpkt.cc @@ -933,8 +933,8 @@ auto vpkt_call_estimators(struct packet *pkt_ptr, const enum packet_type realtyp } } -auto rot_angle(std::span n1, std::span n2, std::span ref1, std::span ref2) - -> double { +auto rot_angle(std::span n1, std::span n2, std::span ref1, + std::span ref2) -> double { // Rotation angle from the scattering plane // We need to rotate Stokes Parameters to (or from) the scattering plane from (or to) // the meridian frame such that Q=1 is in the scattering plane and along ref1