[santacruzgalaxy-list] Grackle proposal
Junhwan Choi (최준환)
choi.junhwan at gmail.com
Thu Mar 20 08:07:12 PDT 2014
Hi Romain,
For me, it took about two weeks to make it works.
But, I spend considerable time to find a way to pass the system unit
to grackle properly.
There are a few outstanding issues on this adjustment.
1. Grackle use C++ and wrap the FORTRAN cooling routine. It is not
trivial to call this C++ wrapped the FORTRAN but there is example in
the grackle example directory.
2. In the gadget case, I pass one particle by one particle to the
grackle in order to calculate the cooling. The reason is that the
gadget cooling and star formation routine is based on particle by
particle. Indeed, it is least efficient way but so far I have not
experienced significant slow down.
[If I use tabulated H and He cooling rates from grackle, it takes
almost the same time with the original gadget cooling.]
3. For unit adjustment, we can set as follow:
my_units.comoving_coordinates = 0;
my_units.density_units = udensity;
my_units.length_units = ulength;
my_units.time_units = utime;
my_units.a_units = 1.0;
[udensity, ulength, and utime is the conversion factor from system
physical units to cgs.]
In this case, you can pass the density, internal energy, and chemistry
value in system physical unit.
Besides these issues, the incorporating Grackle should not be very difficult.
Best,
Junhwan
On Thu, Mar 20, 2014 at 2:35 AM, Romain Teyssier
<romain.teyssier at gmail.com> wrote:
> That's great.
> A very important question: how long did it take him to succeed ?
> Since in my opinion, many "adjustments" need to be done, I would be curious to
> have an estimate on the time required that have this working.
>
> Cheers,
> Romain
>
> On 20 Mar 2014, at 00:33, Ken Nagamine <kn at physics.unlv.edu> wrote:
>
>> Dear all,
>>
>> Let me just mention that, on the SPH side, Junhwan Choi has already succeeded in
>> implementing the Grackle package into Gadget3-UNLV version.
>> We still need to do more checks, but at least the phase diagram
>> looks reasonable compared to the previous TREECOOL table implementation of original Gadget.
>> Junhwan can explain the details of his implementation method if it can be of help for this discussion.
>>
>> cheers,
>> Ken
>>
>>
>> On Mar 20, 2014, at 6:44 AM, Romain Teyssier wrote:
>>
>>>
>>> On 19 Mar 2014, at 22:35, Matthew Turk <matthewturk at gmail.com> wrote:
>>>
>>>> On Wed, Mar 19, 2014 at 5:26 PM, Romain Teyssier
>>>> <romain.teyssier at gmail.com> wrote:
>>>>> Hi Matt,
>>>>>
>>>>>
>>>>> For example:
>>>>>
>>>>> int solve_chemistry(chemistry_data &my_chemistry, ???
>>>>> code_units &my_units, OK set to 1
>>>>> gr_float a_value, gr_float dt_value, why a ?
>>>>> gr_int grid_rank, gr_int *grid_dimension, ????
>>>>> gr_int *grid_start, gr_int *grid_end, ????
>>>>> gr_float *density, gr_float *internal_energy, OK
>>>>> gr_float *x_velocity, gr_float *y_velocity, gr_float *z_velocity, ???
>>>>> gr_float *HI_density, gr_float *HII_density, gr_float *HM_density, too much
>>>>> gr_float *HeI_density, gr_float *HeII_density, gr_float *HeIII_density, too much
>>>>> gr_float *H2I_density, gr_float *H2II_density, too much
>>>>> gr_float *DI_density, gr_float *DII_density, gr_float *HDI_density, too much
>>>>> gr_float *e_density, gr_float *metal_density); OK
>>>>>
>>>>> What are the *grid related variables ?
>>>>
>>>> Good question! These, like the code units variables, are designed to
>>>> minimize the overhead of any simulation that wants to put its code
>>>> into it. For instance, if one had a patch/block-based code (like
>>>> Enzo, FLASH, Nyx, etc) or an Octree code where the 2x2x2 zones were
>>>> included in (2+2*NGZ, 2+2*NGZ, 2+2*NGZ) blocks and one didn't want to
>>>> spend time solving the ghost zones, they can be masked out. These are
>>>> the variables:
>>>>
>>>> grid_start => array of size (grid_rank) indicating the indices to
>>>> start at in each rank of dimensionality
>>>> grid_end => array of size (grid_rank) indicating the indices to *end*
>>>> at in each rank of dimensionality
>>>> grid_dimension => size of the block of data in each dimension
>>>> grid_rank => dimensionality; if you're supplying a pencil beam, this would be 1.
>>>>
>>>> While these typically are better for block or patch based solvers, I
>>>> think that the primary goal -- *minimizing* the overhead to using a
>>>> new code -- is met with them.
>>>
>>> Well I clearly disagree.
>>> It is brain damaging to say the least.
>>>
>>>>
>>>> As I mentioned in our previous email, the velocities are supplied for
>>>> the Sobolev approximation. You can see the full discussion on the
>>>> public Grackle mailing list here:
>>>>
>>>> https://groups.google.com/forum/#!searchin/grackle-cooling-users/velocity/grackle-cooling-users/Dr77TM2te9g/225RNzoZADEJ
>>>>
>>>>>
>>>>> If chemistry_data is an external variable, this is a pain because this variable need to be declared in the main code.
>>>>> Same thing for code_units.
>>>>
>>>> The alternative is that it be a global defined in a different
>>>> namespace, loading in by the dlloader at runtime. Not sure that's
>>>> substantially different.
>>>>
>>>
>>> Well I think this should be completely hidden to the user.
>>> This is a very substantial difference.
>>>
>>>> Code units, if you want to convert to CGS, will just be 1.0 for all of
>>>> them, which can be constant.
>>>>
>>>
>>> Yes but you still have to create a code_unit type, and set the variable to 1.
>>> Unless you change the calling sequence.
>>>
>>>>> This is all too ENZO specifics.
>>>>
>>>> I'm not sure I believe that argument. Perhaps the names of the
>>>> variables share too much with Enzo naming schemes, but I think it's
>>>> rather intuitive to think of the 3D dataset as a base starting point
>>>> and then reducing overhead by supplying 1-dimensional arrays.
>>>>
>>>
>>> For a purely local process like cooling and chemistry, the natural data type is 1D.
>>> Cells, particles and what not, could then be sent as a 1D array to the chemistry solver.
>>> 3D structures are relevant only for patch based codes.
>>>
>>>
>>>> -Matt
>>>>
>>>>>
>>>>> Romain
>>>>>
>>>>>
>>>>>
>>>>> On 19 Mar 2014, at 22:08, Matthew Turk <matthewturk at gmail.com> wrote:
>>>>>
>>>>>> Hi all,
>>>>>>
>>>>>> On Wed, Mar 19, 2014 at 4:58 PM, Romain Teyssier
>>>>>> <romain.teyssier at gmail.com> wrote:
>>>>>>>
>>>>>>> I totally agree with the 3 routines that need to be used.
>>>>>>> Nothing else should be required on the code's side, except passing the required variables.
>>>>>>>
>>>>>>>
>>>>>>>> grackle_init_run()
>>>>>>>> called once per simulation
>>>>>>>>
>>>>>>>
>>>>>>> The parameters that should be passed here are
>>>>>>> - some cosmological model parameters
>>>>>>> - UV background model (spectrum, reionization redshift...)
>>>>>>> - starting redshift to compute the initial temperature and ionisation state
>>>>>>>
>>>>>>>> grackle_init_step(aexp,...)
>>>>>>>> called once per each time-step, uses the value of the scale
>>>>>>>> factor and other parameters as needed
>>>>>>>
>>>>>>> In principle at this stage, only aexp is required for cosmo runs.
>>>>>>> For non cosmo runs, the cooling tables do not need to be modified.
>>>>>>>
>>>>>>>>
>>>>>>>> grackle_solve_chemistry(dt,den,tem,...)
>>>>>>>> called for each resolution element one or more times per one
>>>>>>>> global time-step, uses gas properties to update internal energy.
>>>>>>>> Input gas properties may be required to have specific units and
>>>>>>>> be of specific data type. It should internally sense if it
>>>>>>>> runs within an OpenMP construct and support OpenMP
>>>>>>>> parallelization.
>>>>>>>>
>>>>>>>
>>>>>>> I would rather give the possibility of passing an array of cell.
>>>>>>> Users could also set the array size to 1 to deal with cells one by one.
>>>>>>> The required information could be
>>>>>>> nH, T or Tovermu or specific energy, metallicity in some units, dt_hydro, ncell
>>>>>>> On output, Delta T or Delta Tovermu or Delta specific energy
>>>>>>> I suggest we use fixed units (cgs or mks) for all input and output variables.
>>>>>>>
>>>>>>> These 3 routines are the only one that the user should care about.
>>>>>>>
>>>>>>
>>>>>> This is awfully similar to the existing API:
>>>>>>
>>>>>> https://bitbucket.org/brittonsmith/grackle/src/642cd133535a2dd3c57626b768c7c8c107096ac7/src/clib/grackle.h?at=default
>>>>>>
>>>>>> The only necessary functions, according to
>>>>>> http://grackle.readthedocs.org/en/latest/Integration.html , are:
>>>>>>
>>>>>> initialize_chemistry_data
>>>>>> solve_chemistry or solve_chemistry for tabular data
>>>>>>
>>>>>> Looking at the example (non-equilibrium) executable:
>>>>>>
>>>>>> https://bitbucket.org/brittonsmith/grackle/src/642cd133535a2dd3c57626b768c7c8c107096ac7/src/example/example.C?at=default
>>>>>>
>>>>>> those are the only *functional* routines. Everything else is to
>>>>>> provide additional information, *not* to actually do any computation.
>>>>>> I think that we have actually met these standards. All of the units
>>>>>> are provided so as the *avoid* any boilerplate -- but if you want to
>>>>>> supply in CGS, you can set the input unit conversions to 1.0.
>>>>>>
>>>>>> To compute the equilibrium tables:
>>>>>>
>>>>>> https://bitbucket.org/brittonsmith/grackle/src/642cd133535a2dd3c57626b768c7c8c107096ac7/src/example/table_example.C?at=default
>>>>>>
>>>>>> it's even simpler. The only *active* routine is called on line 120,
>>>>>> "solve_chemistry". I am genuinely being earnest when I ask, from the
>>>>>> .h file linked above, which arguments to the functions would you like
>>>>>> to see removed? I've included the full function signatures below. If
>>>>>> you want to get *out* various things, you can utilize the other
>>>>>> functions -- like the cooling time function and so on.
>>>>>>
>>>>>> -Matt
>>>>>>
>>>>>> int initialize_chemistry_data(chemistry_data &my_chemistry,
>>>>>> code_units &my_units, gr_float a_value);
>>>>>>
>>>>>> int solve_chemistry(chemistry_data &my_chemistry,
>>>>>> code_units &my_units,
>>>>>> gr_float a_value, gr_float dt_value,
>>>>>> gr_int grid_rank, gr_int *grid_dimension,
>>>>>> gr_int *grid_start, gr_int *grid_end,
>>>>>> gr_float *density, gr_float *internal_energy,
>>>>>> gr_float *x_velocity, gr_float *y_velocity, gr_float *z_velocity,
>>>>>> gr_float *HI_density, gr_float *HII_density, gr_float *HM_density,
>>>>>> gr_float *HeI_density, gr_float *HeII_density, gr_float *HeIII_density,
>>>>>> gr_float *H2I_density, gr_float *H2II_density,
>>>>>> gr_float *DI_density, gr_float *DII_density, gr_float *HDI_density,
>>>>>> gr_float *e_density, gr_float *metal_density);
>>>>>>
>>>>>> int solve_chemistry(chemistry_data &my_chemistry,
>>>>>> code_units &my_units,
>>>>>> gr_float a_value, gr_float dt_value,
>>>>>> gr_int grid_rank, gr_int *grid_dimension,
>>>>>> gr_int *grid_start, gr_int *grid_end,
>>>>>> gr_float *density, gr_float *internal_energy,
>>>>>> gr_float *x_velocity, gr_float *y_velocity,
>>>>>> gr_float *z_velocity,
>>>>>> gr_float *metal_density);
>>>>>>
>>>>>>> Cheers,
>>>>>>> Romain
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>>> The API can also define its own data types with natutal conversion from
>>>>>>>> standard types, for example the following code should be valid:
>>>>>>>>
>>>>>>>> float f;
>>>>>>>> double d;
>>>>>>>> gr_float gr_f;
>>>>>>>>
>>>>>>>> gr_f = f;
>>>>>>>> gr_f = d;
>>>>>>>>
>>>>>>>> An analogous API should be provided for F77.
>>>>>>>>
>>>>>>>> Then in a code the API will be implemented as follows:
>>>>>>>>
>>>>>>>> Begin_code
>>>>>>>>
>>>>>>>> grackle_init_run()
>>>>>>>>
>>>>>>>> Loop_over_timesteps(aexp)
>>>>>>>>
>>>>>>>> grackle_init_step(aexp,...)
>>>>>>>>
>>>>>>>> #pragma omp parallel loop
>>>>>>>> Loop_over_resolution_elements(elem)
>>>>>>>> {
>>>>>>>> gr_float dt = code_time_step*time_unit;
>>>>>>>> gr_float den = code_density(elem)*den_unit;
>>>>>>>> gr_float tem = code_temperature(elem)*tem_unit;
>>>>>>>> grackle_solve_chemistry(dt,den,tem,...)
>>>>>>>> }
>>>>>>>>
>>>>>>>> End_loop
>>>>>>>>
>>>>>>>> End_code
>>>>>>>>
>>>>>>>> Let's converge on the API, and then the Grackle team will be able to
>>>>>>>> write a couple of wrappers that will suit everyone.
>>>>>>>>
>>>>>>>> Given the wrappers, if Grackle is installed as an external library and
>>>>>>>> provides a proper Linus-style installer (that will handle all HDF5 and
>>>>>>>> other dependencies), then it should be trivial to integrate it in any code.
>>>>>>>>
>>>>>>>> Nick
>>>>>>>>
>>>>>>>> To unsubscribe from this group and stop receiving emails from it, send an email to santacruzgalaxy-list+unsubscribe at ucsc.edu.
>>>>>>>> _______________________________________________
>>>>>>>> santacruzgalaxy-list mailing list
>>>>>>>> santacruzgalaxy-list at ucsc.edu
>>>>>>>> https://lists.ucsc.edu/mailman/listinfo/santacruzgalaxy-list
>>>>>>>
>>>>>>> --
>>>>>>> You received this message because you are subscribed to the Google Groups "Grackle Cooling Library Users" group.
>>>>>>> To unsubscribe from this group and stop receiving emails from it, send an email to grackle-cooling-users+unsubscribe at googlegroups.com.
>>>>>>> For more options, visit https://groups.google.com/d/optout.
>>>
>>>
>>> To unsubscribe from this group and stop receiving emails from it, send an email to santacruzgalaxy-list+unsubscribe at ucsc.edu.
>>> _______________________________________________
>>> santacruzgalaxy-list mailing list
>>> santacruzgalaxy-list at ucsc.edu
>>> https://lists.ucsc.edu/mailman/listinfo/santacruzgalaxy-list
>>
>
> --
> You received this message because you are subscribed to the Google Groups "Grackle Cooling Library Users" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to grackle-cooling-users+unsubscribe at googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.
To unsubscribe from this group and stop receiving emails from it, send an email to santacruzgalaxy-list+unsubscribe at ucsc.edu.
More information about the santacruzgalaxy-list
mailing list