API tools faq
paste
Advertisement
mabruzzo

arg-reduced cool_multi_time_g-cpp.C

mabruzzo
Dec 4th, 2024
57
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 24.37 KB | None | 0 0
raw download clone embed print report
  1. #include <cstdio>
  2. #include <vector>
  3.  
  4. #include "grackle.h"
  5. #include "fortran_func_decls.h"
  6. #include "utils-cpp.hpp"
  7. #include "cool_multi_time_g-cpp.h"
  8.  
  9. //_// TODO: ADD ANY OTHER INCLUDE DIRECTIVES
  10.  
  11. extern "C" {
  12. void cool_multi_time_g(
  13.   gr_float* cooltime_data_, int* imetal, double* utem, double* uxyz,
  14.   double* urho, chemistry_data* my_chemistry, chemistry_data_storage* my_rates,
  15.   code_units* my_units, grackle_field_data* my_fields,
  16.   photo_rate_storage* my_uvb_rates
  17. )
  18. {
  19.  
  20.   //   SOLVE RADIATIVE COOLING/HEATING EQUATIONS
  21.  
  22.   //   written by: Yu Zhang, Peter Anninos and Tom Abel
  23.   //   date:
  24.   //   modified1: January, 1996 by Greg Bryan; adapted to KRONOS
  25.   //   modified2: October, 1996 by GB; moved to AMR
  26.   //   modified3: February, 2003 by Robert Harkness; iteration mask
  27.  
  28.   //   PURPOSE:
  29.   //     Solve the energy cooling equations.
  30.  
  31.   //   INPUTS:
  32.   //     is,ie   - start and end indicies of active region (zero-based!)
  33.  
  34.   //   PARAMETERS:
  35.  
  36.   // -----------------------------------------------------------------------
  37.  
  38.  
  39.   // Arguments
  40.  
  41.   // -- removed line (previously just declared arg types) --
  42.  
  43.   // -- removed line (previously just declared arg types) --
  44.   grackle::impl::View<gr_float***> cooltime(cooltime_data_, my_fields->grid_dimension[0], my_fields->grid_dimension[1], my_fields->grid_dimension[2]);
  45.   // -- removed line (previously just declared arg types) --
  46.   // -- removed line (previously just declared arg types) --
  47.   // -- removed line (previously just declared arg types) --
  48.   // -- removed line (previously just declared arg types) --
  49.   // -- removed line (previously just declared arg types) --
  50.   // -- removed line (previously just declared arg types) --
  51.   // -- removed line (previously just declared arg types) --
  52.   // -- removed line (previously just declared arg types) --
  53.   // -- removed line (previously just declared arg types) --
  54.   // -- removed line (previously just declared arg types) --
  55.   // -- removed line (previously just declared arg types) --
  56.   // -- removed line (previously just declared arg types) --
  57.   // -- removed line (previously just declared arg types) --
  58.   // -- removed line (previously just declared arg types) --
  59.   // -- removed line (previously just declared arg types) --
  60.   // -- removed line (previously just declared arg types) --
  61.   // -- removed line (previously just declared arg types) --
  62.   // -- removed line (previously just declared arg types) --
  63.   // -- removed line (previously just declared arg types) --
  64.   // -- removed line (previously just declared arg types) --
  65.   // -- removed line (previously just declared arg types) --
  66.   // -- removed line (previously just declared arg types) --
  67.   // -- removed line (previously just declared arg types) --
  68.   // opacity table
  69.   // -- removed line (previously just declared arg types) --
  70.   // -- removed line (previously just declared arg types) --
  71.   // grain temperature
  72.   std::vector<double> tSiM(my_fields->grid_dimension[0]);
  73.   std::vector<double> tFeM(my_fields->grid_dimension[0]);
  74.   std::vector<double> tMg2SiO4(my_fields->grid_dimension[0]);
  75.   std::vector<double> tMgSiO3(my_fields->grid_dimension[0]);
  76.   std::vector<double> tFe3O4(my_fields->grid_dimension[0]);
  77.   std::vector<double> tAC(my_fields->grid_dimension[0]);
  78.   std::vector<double> tSiO2D(my_fields->grid_dimension[0]);
  79.   std::vector<double> tMgO(my_fields->grid_dimension[0]);
  80.   std::vector<double> tFeS(my_fields->grid_dimension[0]);
  81.   std::vector<double> tAl2O3(my_fields->grid_dimension[0]);
  82.   std::vector<double> treforg(my_fields->grid_dimension[0]);
  83.   std::vector<double> tvolorg(my_fields->grid_dimension[0]);
  84.   std::vector<double> tH2Oice(my_fields->grid_dimension[0]);
  85.   // -- removed line (previously just declared arg types) --
  86.  
  87.   // Cloudy cooling data
  88.  
  89.   // -- removed line (previously just declared arg types) --
  90.   // -- removed line (previously just declared arg types) --
  91.   // -- removed line (previously just declared arg types) --
  92.  
  93.   // Parameters
  94.  
  95.   // Locals
  96.  
  97.   int i, j, k;
  98.   int t, dj, dk;
  99.   double comp1, comp2, energy;
  100.   gr_float factor;
  101.  
  102.   // Slice locals
  103.  
  104.   std::vector<long long> indixe(my_fields->grid_dimension[0]);
  105.   std::vector<double> t1(my_fields->grid_dimension[0]);
  106.   std::vector<double> t2(my_fields->grid_dimension[0]);
  107.   std::vector<double> logtem(my_fields->grid_dimension[0]);
  108.   std::vector<double> tdef(my_fields->grid_dimension[0]);
  109.   std::vector<double> p2d(my_fields->grid_dimension[0]);
  110.   std::vector<double> tgas(my_fields->grid_dimension[0]);
  111.   std::vector<double> tgasold(my_fields->grid_dimension[0]);
  112.   std::vector<double> mmw(my_fields->grid_dimension[0]);
  113.   std::vector<double> tdust(my_fields->grid_dimension[0]);
  114.   std::vector<double> metallicity(my_fields->grid_dimension[0]);
  115.   std::vector<double> dust2gas(my_fields->grid_dimension[0]);
  116.   std::vector<double> rhoH(my_fields->grid_dimension[0]);
  117.   std::vector<double> mynh(my_fields->grid_dimension[0]);
  118.   std::vector<double> myde(my_fields->grid_dimension[0]);
  119.   std::vector<double> gammaha_eff(my_fields->grid_dimension[0]);
  120.   std::vector<double> gasgr_tdust(my_fields->grid_dimension[0]);
  121.   std::vector<double> regr(my_fields->grid_dimension[0]);
  122.   std::vector<double> edot(my_fields->grid_dimension[0]);
  123.  
  124.   // Cooling/heating slice locals
  125.  
  126.   std::vector<double> ceHI(my_fields->grid_dimension[0]);
  127.   std::vector<double> ceHeI(my_fields->grid_dimension[0]);
  128.   std::vector<double> ceHeII(my_fields->grid_dimension[0]);
  129.   std::vector<double> ciHI(my_fields->grid_dimension[0]);
  130.   std::vector<double> ciHeI(my_fields->grid_dimension[0]);
  131.   std::vector<double> ciHeIS(my_fields->grid_dimension[0]);
  132.   std::vector<double> ciHeII(my_fields->grid_dimension[0]);
  133.   std::vector<double> reHII(my_fields->grid_dimension[0]);
  134.   std::vector<double> reHeII1(my_fields->grid_dimension[0]);
  135.   std::vector<double> reHeII2(my_fields->grid_dimension[0]);
  136.   std::vector<double> reHeIII(my_fields->grid_dimension[0]);
  137.   std::vector<double> brem(my_fields->grid_dimension[0]);
  138.   std::vector<double> cieco(my_fields->grid_dimension[0]);
  139.   std::vector<double> hyd01k(my_fields->grid_dimension[0]);
  140.   std::vector<double> h2k01(my_fields->grid_dimension[0]);
  141.   std::vector<double> vibh(my_fields->grid_dimension[0]);
  142.   std::vector<double> roth(my_fields->grid_dimension[0]);
  143.   std::vector<double> rotl(my_fields->grid_dimension[0]);
  144.   std::vector<double> gpldl(my_fields->grid_dimension[0]);
  145.   std::vector<double> gphdl(my_fields->grid_dimension[0]);
  146.   std::vector<double> hdlte(my_fields->grid_dimension[0]);
  147.   std::vector<double> hdlow(my_fields->grid_dimension[0]);
  148.   int dummy_iter_arg;
  149.          
  150.   // Iteration mask for multi_cool
  151.  
  152.   std::vector<gr_mask_type> itmask(my_fields->grid_dimension[0]);
  153.   std::vector<gr_mask_type> itmask_metal(my_fields->grid_dimension[0]);
  154.  
  155.   // \\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\\/////////////////////////////////
  156.   // =======================================================================
  157.   dk = my_fields->grid_end[2] - my_fields->grid_start[2] + 1;
  158.   dj = my_fields->grid_end[1] - my_fields->grid_start[1] + 1;
  159.  
  160.   // Convert densities from comoving to 'proper'
  161.  
  162.   if (my_units->comoving_coordinates == 1)  {
  163.  
  164.     factor = (gr_float)(std::pow(my_units->a_value,(-3)) );
  165.  
  166.      FORTRAN_NAME(scale_fields_g)(
  167.       &my_chemistry->primordial_chemistry, imetal, &my_chemistry->use_dust_density_field, &my_chemistry->metal_chemistry,
  168.       &my_chemistry->metal_abundances, &my_chemistry->dust_species, &my_chemistry->multi_metals, &my_chemistry->grain_growth, &my_chemistry->dust_sublimation, &factor,
  169.       &my_fields->grid_start[0], &my_fields->grid_end[0], &my_fields->grid_start[1], &my_fields->grid_end[1], &my_fields->grid_start[2], &my_fields->grid_end[2], &my_fields->grid_dimension[0], &my_fields->grid_dimension[1], &my_fields->grid_dimension[2],
  170.       my_fields->density, my_fields->e_density, my_fields->HI_density, my_fields->HII_density, my_fields->HeI_density, my_fields->HeII_density, my_fields->HeIII_density,
  171.       my_fields->HM_density, my_fields->H2I_density, my_fields->H2II_density, my_fields->DI_density, my_fields->DII_density, my_fields->HDI_density, my_fields->DM_density, my_fields->HDII_density, my_fields->HeHII_density,
  172.       my_fields->metal_density, my_fields->dust_density,
  173.       my_fields->CI_density, my_fields->CII_density, my_fields->CO_density, my_fields->CO2_density,
  174.       my_fields->OI_density, my_fields->OH_density, my_fields->H2O_density, my_fields->O2_density,
  175.       my_fields->SiI_density, my_fields->SiOI_density, my_fields->SiO2I_density,
  176.       my_fields->CH_density, my_fields->CH2_density, my_fields->COII_density, my_fields->OII_density,
  177.       my_fields->OHII_density, my_fields->H2OII_density, my_fields->H3OII_density, my_fields->O2II_density,
  178.       my_fields->Mg_density, my_fields->Al_density, my_fields->S_density, my_fields->Fe_density,
  179.       my_fields->SiM_dust_density, my_fields->FeM_dust_density, my_fields->Mg2SiO4_dust_density, my_fields->MgSiO3_dust_density, my_fields->Fe3O4_dust_density,
  180.       my_fields->AC_dust_density, my_fields->SiO2_dust_density, my_fields->MgO_dust_density, my_fields->FeS_dust_density, my_fields->Al2O3_dust_density,
  181.       my_fields->ref_org_dust_density, my_fields->vol_org_dust_density, my_fields->H2O_ice_dust_density,
  182.       my_fields->local_ISM_metal_density, my_fields->ccsn13_metal_density, my_fields->ccsn20_metal_density, my_fields->ccsn25_metal_density, my_fields->ccsn30_metal_density,
  183.       my_fields->fsn13_metal_density, my_fields->fsn15_metal_density, my_fields->fsn50_metal_density, my_fields->fsn80_metal_density,
  184.       my_fields->pisn170_metal_density, my_fields->pisn200_metal_density, my_fields->y19_metal_density);
  185.  
  186.   }
  187.  
  188.   // Loop over slices (in the k-direction)
  189.  
  190.   // parallelize the k and j loops with OpenMP
  191.   // flat j and k loops for better parallelism
  192.   //_// PORT: #ifdef _OPENMP
  193.   //_// PORT: !$omp parallel do schedule(runtime) private(
  194.   //_// PORT: !$omp&   i, j, k,
  195.   //_// PORT: !$omp&   comp1, comp2, energy,
  196.   //_// PORT: !$omp&   indixe,
  197.   //_// PORT: !$omp&   t1, t2, logtem, tdef, p2d,
  198.   //_// PORT: !$omp&   tgas, tgasold,
  199.   //_// PORT: !$omp&   tdust, metallicity, dust2gas, rhoH, mmw,
  200.   //_// PORT: !$omp&   mynh, myde, gammaha_eff, gasgr_tdust, regr, edot,
  201.   //_// PORT: !$omp&   ceHI, ceHeI, ceHeII,
  202.   //_// PORT: !$omp&   ciHI, ciHeI, ciHeIS, ciHeII,
  203.   //_// PORT: !$omp&   reHII, reHeII1, reHeII2, reHeIII,
  204.   //_// PORT: !$omp&   brem, cieco,
  205.   //_// PORT: !$omp&   hyd01k, h2k01, vibh, roth, rotl,
  206.   //_// PORT: !$omp&   gpldl, gphdl, hdlte, hdlow,
  207.   //_// PORT: !$omp&   itmask )
  208.   //_// PORT: #endif
  209.   for (t = 0; t<=(dk * dj - 1); t++) {
  210.     k = t/dj      + my_fields->grid_start[2]+1;
  211.     j = grackle::impl::mod(t,dj) + my_fields->grid_start[1]+1;
  212.  
  213.     for (i = my_fields->grid_start[0] + 1; i<=(my_fields->grid_end[0] + 1); i++) {
  214.       itmask[i-1] = MASK_TRUE;
  215.     }
  216.  
  217.     // Compute the cooling rate// Compute the cooling time on the slice
  218.     //   (the gamma used here is the same as used to calculate the pressure
  219.     //    in cool1d_multi_g)
  220.  
  221.     for (i = my_fields->grid_start[0] + 1; i<=(my_fields->grid_end[0] + 1); i++) {
  222.       energy = std::fmax(p2d[i-1]/(my_chemistry->Gamma-1.), tiny);
  223.       cooltime(i-1,j-1,k-1) = (gr_float)(energy/edot[i-1] );
  224.     }
  225.  
  226.   }
  227.   //_// PORT: #ifdef _OPENMP
  228.   //_// PORT: !$omp end parallel do
  229.   //_// PORT: #endif
  230.  
  231.   // Convert densities back to comoving from 'proper'
  232.  
  233.   if (my_units->comoving_coordinates == 1)  {
  234.  
  235.     factor = (gr_float)(std::pow(my_units->a_value,3) );
  236.  
  237.      FORTRAN_NAME(scale_fields_g)(
  238.       &my_chemistry->primordial_chemistry, imetal, &my_chemistry->use_dust_density_field, &my_chemistry->metal_chemistry,
  239.       &my_chemistry->metal_abundances, &my_chemistry->dust_species, &my_chemistry->multi_metals, &my_chemistry->grain_growth, &my_chemistry->dust_sublimation, &factor,
  240.       &my_fields->grid_start[0], &my_fields->grid_end[0], &my_fields->grid_start[1], &my_fields->grid_end[1], &my_fields->grid_start[2], &my_fields->grid_end[2], &my_fields->grid_dimension[0], &my_fields->grid_dimension[1], &my_fields->grid_dimension[2],
  241.       my_fields->density, my_fields->e_density, my_fields->HI_density, my_fields->HII_density, my_fields->HeI_density, my_fields->HeII_density, my_fields->HeIII_density,
  242.       my_fields->HM_density, my_fields->H2I_density, my_fields->H2II_density, my_fields->DI_density, my_fields->DII_density, my_fields->HDI_density, my_fields->DM_density, my_fields->HDII_density, my_fields->HeHII_density,
  243.       my_fields->metal_density, my_fields->dust_density,
  244.       my_fields->CI_density, my_fields->CII_density, my_fields->CO_density, my_fields->CO2_density,
  245.       my_fields->OI_density, my_fields->OH_density, my_fields->H2O_density, my_fields->O2_density,
  246.       my_fields->SiI_density, my_fields->SiOI_density, my_fields->SiO2I_density,
  247.       my_fields->CH_density, my_fields->CH2_density, my_fields->COII_density, my_fields->OII_density,
  248.       my_fields->OHII_density, my_fields->H2OII_density, my_fields->H3OII_density, my_fields->O2II_density,
  249.       my_fields->Mg_density, my_fields->Al_density, my_fields->S_density, my_fields->Fe_density,
  250.       my_fields->SiM_dust_density, my_fields->FeM_dust_density, my_fields->Mg2SiO4_dust_density, my_fields->MgSiO3_dust_density, my_fields->Fe3O4_dust_density,
  251.       my_fields->AC_dust_density, my_fields->SiO2_dust_density, my_fields->MgO_dust_density, my_fields->FeS_dust_density, my_fields->Al2O3_dust_density,
  252.       my_fields->ref_org_dust_density, my_fields->vol_org_dust_density, my_fields->H2O_ice_dust_density,
  253.       my_fields->local_ISM_metal_density, my_fields->ccsn13_metal_density, my_fields->ccsn20_metal_density, my_fields->ccsn25_metal_density, my_fields->ccsn30_metal_density,
  254.       my_fields->fsn13_metal_density, my_fields->fsn15_metal_density, my_fields->fsn50_metal_density, my_fields->fsn80_metal_density,
  255.       my_fields->pisn170_metal_density, my_fields->pisn200_metal_density, my_fields->y19_metal_density);
  256.  
  257.   }
  258.  
  259.   return;
  260. }
  261.  
  262. }  // extern "C"
  263.  
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!