This PR adds support for the Grackle chemistry and cooling library to Enzo. Use of Grackle requires compiling with "make grackle-yes". All Grackle functionality should be inside #ifdef USE_GRACKLE, but please leave a comment if you find anything that is not.
As I mentioned in my email to enzo-dev a while ago, Grackle is a chemistry and cooling library that was originally derived from Enzo's MultiSpecies and Cloudy cooling modules. However, there have been a number of updates to the UV background and metal cooling functionality make this not only unique from the original implementation, but potentially a more attractive option for simulations without radiative transfer, such as:
UV backgrounds are treated via interpolating from data tables loaded from disk rather
than piece-wise polynomial functions, making it somewhat easier to add support for
new background models.
UV background and cooling data are contained within the same input file and are
more consistent than the currently available Cloudy cooling data and Enzo UV
background models. Currently, the Grackle distribution comes with UV background
and cooling data for two different models:
Unlike the original Cloudy cooling which required separate input files for cooling
before the UV background turns on and after, all data is contained in a single
table. This means one no longer has to run with one input file to the redshift
where the UV background starts, stop the simulation, and restart with another input
Also unlike the original Cloudy cooling module, Grackle supports the option to also
solve the primordial cooling via interpolation from a table. Thus, one is no longer
required to run with the MultiSpecies functionality in order to calculate the
primordial component. This simplified method is somewhat faster and requires fewer
baryon fields to be stored, lowering the ram and disk footprint.
This also enables the use of Enzo for the AGORA project with the required common physics package, a.k.a. Grackle.
What has been added:
support to solve chemistry cooling and calculate cooling time using Grackle functions
narrative documentation discussing Grackle in general, how to build Enzo with it, and a full description of parameters
sample cosmology simulation parameter file using Grackle in run/CosmologySimulation/AMRCosmology_Grackle
I ran a very small cosmological simulation using Cloudy cooling, Grackle with MultiSpecies 1, and Grackle with MultiSpecies 0 (fully-tabulated). Note, the input data are not the same for the Cloudy run and the Grackle runs. The Grackle data corresponds to a newer version of the Haardt & Madau background which starts at a higher redshift (z = 14.849 instead of z = 8.9). Below are projections of metallicity and phase plots of cell mass in bins of density and temperature for each of the three simulations. The strange centering in the projections is due to the fact that this is a pretty low-resolution nested simulation and we are viewing only the refinement region.
Grackle w/ MultiSpecies = 1
Grackle Fully-tabulated (note: phase plot y range is different)
The most notable different is the addition of a small amount of cold gas at roughly 1e3 K in the Grackle runs. I'm not entirely sure where exactly this comes from, but I don't think it is indicative of something incorrect.
I've tested this out, and once I actually follow the directions it works great! Thanks for all of this work, @Britton Smith . Would you be up for a G+ dev hangout in the near future to have folks give feedback on the structure of enzo <--> grackle?
First, let me say that I'm very excited that you're PR'ing this, as I'm eager to use it for my own sims. Thank you for your hard work!
I followed your installation instruction and was able to get grackle installed on TACC stampede without any problems. I then followed your instructions for installing and building grackle support in enzo. This was pretty straightforward too. Finally, I tried to run your AMRCosmology Test simulation in the run directory, but I'm getting failures. Essentially, it builds the ICs fine, runs ring ok, and then starts up the simulation, but it keeps repeating this line in STDERR:
"Mean molecular weight not converged! NaN NaN"
I found this in the Grackle Source for cool1d_multi_g.F, in that it cannot seem to converge on Mu when adding in the different metal species. But I'm not exactly sure why this is happening. I'm just running with all of your defaults. It looks like it gets past all of the initializations and is actually starting calculation to cause this error. Any ideas?
To follow up on this comment, I was able to get enzo-grackle installed and working on a Darwin box (mountain lion). Everything seemed to work fine after I made the environment variable modification noted below so that enzo could find the grackle libraries. But I still don't know what's going awry with grackle-enzo on Stampede as discussed above.
@chummels has correctly pointed out to me that the grackle documentation states that you can use your Enzo make file to build grackle, which is not true and seems to be the cause of these problems. I will correct this as soon as I can.
Using gdb with this, it appears the error occurs when c++ passes pointers to fortran programs, specifically in the call to solve_rate_cool_g.F, many of the values it gets passed are recognized, which yield a quick segfault. I'm not sure how to fix this--perhaps some additional compiler flag in the Makefile? I'm sort of out of ideas here.
My guess is that it's related to -fdefault-real-8 and -fdefault-double-8 in your makefile. Try seeing if it fails in 32 bit mode. Also, did you mean "many of the values it gets passed aren't recognized"? C++ and F don't do interface checking, they just toss a set of values from one to the other and hope for the best.
I get the "molecular weight not converged" errors (repeating into infinity) when actually doing the enzo grackle test run in the run directory. But when I just try to run the grackle example, that is when I get the errors above (e.g. NaN in edot...).
This is because the fortran routines that were taken out of Enzo to become the Grackle predate @Daniel Reynolds's precision fixes. I will update the source and documentation as soon as I have a free moment. I think I have this working with static libraries on NICS Kraken as well.
So as I noted in a previous comment, those flags don't work for me on stampede. However, I found a related solution. By using the flags:
I was able to get stampede to recognize the appropriate fortran precision and accept the C++ values correctly. Thanks for your insight, @MattT and your perseverance @Britton Smith ! I'll PR that Makefile for stampede to the grackle repo, and I now approve this PR, since I can now run the AMRCosmology_Grackle test problem on an OS X (mountain lion) box and Stampede.
N.B. There was never anything wrong with this PR, it was just a problem with grackle's compilation itself.
The only thing I would like to happen before this gets accepted is for me to update the Grackle documentation and make files to avoid the type of errors we all saw. I will leave a comment in this PR when that happens. Does anyone else have any comments they would like to make before this goes in?
@Britton Smith this is really great, and functionally I think this could be accepted and merged right away. The only remaining question I have comes from how the grackle compute_pressure, calculate_temperature, and compute_gamma compare to the methods in enzo itself. It seems like at first glance they might be equivalent, but could you comment on how you see enzo + grackle evolving if/when those methods diverge from each other? As I mentioned above, I think it may be useful to have a hangout so that we could talk about the API between enzo and grackle.
In the instances that are identical to the Grackle, I would eventually like to switch over to using the Grackle functions. However, in some of these, there are other cases that the Grackle does not explicitly cover. I think the best thing to do will be to keep the Enzo functions and have them call the Grackle functions inside or do the alternate calculation themselves. I would be interested in having a hangout to talk about the future of this. Is anyone other than @Sam Skillman interested?