Move ions[ion].Alpha_sp to block allocation ion_alpha_sp

Brackets
artis-mcrt · Aug 30, 2024 · db80d03 · db80d03
1 parent 0353028
commit db80d03
Show file tree

Hide file tree

Showing 4 changed files with 11 additions and 16 deletions.
diff --git a/globals.h b/globals.h
@@ -108,7 +108,6 @@ struct Ion {
   int ncoolingterms;
   int uniquelevelindexstart;
   int groundcontindex;
-  float *Alpha_sp;
   double ionpot;  // Ionisation threshold to the next ionstage
 };
 

diff --git a/input.cc b/input.cc
@@ -823,9 +823,6 @@ void read_atomicdata_files() {
       globals::elements[element].ions[ion].groundcontindex = -1;
       globals::elements[element].ions[ion].first_nlte = -1;
 
-      globals::elements[element].ions[ion].Alpha_sp = static_cast<float *>(calloc(TABLESIZE, sizeof(float)));
-      assert_always(globals::elements[element].ions[ion].Alpha_sp != nullptr);
-
       globals::elements[element].ions[ion].levels = static_cast<EnergyLevel *>(calloc(nlevelsmax, sizeof(EnergyLevel)));
       assert_always(globals::elements[element].ions[ion].levels != nullptr);
 

diff --git a/ltepop.cc b/ltepop.cc
@@ -28,23 +28,20 @@ struct nneSolutionParas {
   bool force_lte;
 };
 
-auto interpolate_ions_spontrecombcoeff(const int element, const int ion, const double T) -> double {
-  assert_testmodeonly(element < get_nelements());
-  assert_testmodeonly(ion < get_nions(element));
-  assert_testmodeonly(T >= MINTEMP);
-
+auto interpolate_ions_spontrecombcoeff(const int uniqueionindex, const double T) -> double {
   const int lowerindex = floor(log(T / MINTEMP) / T_step_log);
+  assert_testmodeonly(lowerindex >= 0);
   if (lowerindex < TABLESIZE - 1) {
     const int upperindex = lowerindex + 1;
     const double T_lower = MINTEMP * exp(lowerindex * T_step_log);
     const double T_upper = MINTEMP * exp(upperindex * T_step_log);
 
-    const double f_upper = globals::elements[element].ions[ion].Alpha_sp[upperindex];
-    const double f_lower = globals::elements[element].ions[ion].Alpha_sp[lowerindex];
+    const double f_upper = globals::ion_alpha_sp[(uniqueionindex * TABLESIZE) + upperindex];
+    const double f_lower = globals::ion_alpha_sp[(uniqueionindex * TABLESIZE) + lowerindex];
 
     return f_lower + ((f_upper - f_lower) / (T_upper - T_lower) * (T - T_lower));
   }
-  return globals::elements[element].ions[ion].Alpha_sp[TABLESIZE - 1];
+  return globals::ion_alpha_sp[(uniqueionindex * TABLESIZE) + TABLESIZE - 1];
 }
 
 // use Saha equation for LTE ionization balance
@@ -93,7 +90,7 @@ auto phi_ion_equilib(const int element, const int ion, const int modelgridindex,
     std::abort();
   }
 
-  const double Alpha_sp = interpolate_ions_spontrecombcoeff(element, ion, T_e);
+  const double Alpha_sp = interpolate_ions_spontrecombcoeff(uniqueionindex, T_e);
 
   // const double Col_rec = calculate_ionrecombcoeff(modelgridindex, T_e, element, ion + 1, false, true, false, false,
   // false);

diff --git a/ratecoeff.cc b/ratecoeff.cc
@@ -645,21 +645,23 @@ void read_recombrate_file() {
 }
 
 void precalculate_ion_alpha_sp() {
+  globals::ion_alpha_sp.resize(get_includedions() * TABLESIZE);
   for (int iter = 0; iter < TABLESIZE; iter++) {
     const float T_e = MINTEMP * exp(iter * T_step_log);
     for (int element = 0; element < get_nelements(); element++) {
       const int nions = get_nions(element) - 1;
       for (int ion = 0; ion < nions; ion++) {
-        const int nlevels = get_ionisinglevels(element, ion);
+        const auto uniqueionindex = get_uniqueionindex(element, ion);
+        const int nionisinglevels = get_ionisinglevels(element, ion);
         double zeta = 0.;
-        for (int level = 0; level < nlevels; level++) {
+        for (int level = 0; level < nionisinglevels; level++) {
           const auto nphixstargets = get_nphixstargets(element, ion, level);
           for (int phixstargetindex = 0; phixstargetindex < nphixstargets; phixstargetindex++) {
             const double zeta_level = get_spontrecombcoeff(element, ion, level, phixstargetindex, T_e);
             zeta += zeta_level;
           }
         }
-        globals::elements[element].ions[ion].Alpha_sp[iter] = zeta;
+        globals::ion_alpha_sp[(uniqueionindex * TABLESIZE) + iter] = zeta;
       }
     }
   }