Frontiers in Computational Astrophysics: The Origin of Stars, Planets and Galaxies
An ESF exploratory workshop
Here we describe several tests that we would like participants to try in readiness for the meeting. Please attempt what you can. At the meeting we will spend some time discussing results and to refine the tests that we think highlight the positive and negative characteristics of different codes. We may also come up with more appropriate tests at the meeting. For now you should look at standard diagnostics, such as temperature, pressure or density maps through central slices etc. We can decide on detailed diagnostic comparisons later. We hope to eventually publish a joint paper summarising our efforts.
We decided to set up the tests using particles and provide a code for transforming particle distributions to an accurate realisation on a grid. This tool currently outputs binary files that can be read using the HDF tools. We have given particle IC's at several different resolutions which you can run directly. If you are transforming to a grid, please use the highest resolution particle IC's as an input to the conversion tool and you can specify the size of the output grid. If you run the simulations at several resolutions we can try to determine at which resolution different codes agree - please run these at the highest resolution you can in the time we have left.
The particle initial conditions are provided in both standard tipsy RPC-binary (no endianness problems) format and in ascii format. TIPSY is a nice particle visualization tool written by Neal Katz and Thomas Quinn. There are a few common conversions to other formats to be found following the above link. If you need to write your own conversion macro, the structure of the tipsy ascii file is described here. One odditiy is the hsmooth field for the gas particles which, by its name, should contain the smoothing length for SPH actually contains the gravitational softening for the particle (the same as the eps field for dark and star particles). The reason for this is that the gravitational softening is a more fundamental parameter, while the smoothing length is typically determined from the 32 nearest neighbors in our case. TIPSY gas particles should really have both fields, but there are always historical constraints with data formats.
tipgrid & smooth
# tipgrid –
Will give you a brief description of the cmdline arguments, and I will try to put together a man page for this tool in the next week. Also, there will likely by several updates to this code to strengthen the types of visualizations it can produce (right now just projected density) adding slices, max-in-column and other variables (temperature, entropy, pressure).
The code TIPGRID creates HDF5-ready files in binary, but can also produce ascii files if requested. This conversion to grid is broken up into 3 steps:
Assume the provided input particle file (in standard TIPSY binary format) is called
test3.bin in what follows.
# smooth -s 32g -px 2000 -py 2000 -pz 8000 hsmooth -o test3 < test3.bin
# hsmtoeps test3.bin test3.hsm > test3.eps
# tipgrid -v -n 64 64 256 -l 2000 -xpd 1 -ypd 2 -p -dmin 0 -tmin 0 -msol 2.3262e5 -m 4 -o test3 < test3.eps
-asciiflag to the tipgrid command line. This produces an additional 5 files as well as a gzipped tar file containing these 5 ascii files. We provide this .tgz (gzipped tar file) for the grid intial conditions below. The x-dimension in these ascii files is again the most rapidly changing. To extract the individual ascii files from the tar files given below just type the following.
# tar xfvz test1-128.gasc.tgz
# h5import test3.grid.d -c test3.cfg.d test3.grid.t -c test3.cfg.t test3.grid.vx -c test3.cfg.vx test3.grid.vy -c test3.cfg.vy test3.grid.vz -c test3.cfg.vz -o test3.h5
# h5import test3.plane.d -c test3.pcfg.d -o test3.proj.h5
# h5topng -c jet test3.proj.h5
# eog/xv test3.proj.h5.png
For each of the four tests we provide this sequence of commands to generate the grid initial condition.
Units used in the simulations are based on a G=1 system, which provides a convenient system when the length and velocity units are in kpc and km/s respectively; namely that the time unit is very nearly Gyrs. This system of units is thus very practical in astrophysical applications. For converting to other unit systems the following should be helpful: the time unit is given by the gravitational timescale t_unit = 1/sqrt(G_physical_value * density_unit); the energy per unit mass unit is given by energy_per_unit_mass_unit = G_physical * mass_unit / length_unit. Thus the gas constant in system units is given by gas_const = (length_unit/(G_physical_value*mass_unit))*k_boltzman/mass_hydrogen(or amu).
To convert temperature to thermal energy U = gas_const/(mean_molar_mass)/(gamma - 1)*Temperature. The mass_unit = 2.3262e5 solar masses, the length_unit = 1 kpc, the velocity_unit = 1 km/s, the time_unit = 0.9777 Gyrs. For SPH codes one can convert masses, positions, velocities, softening radii, and temperatures to other system units if this is required using the above relations. For GRID codes we provide the 5 Euler equation quantities in these units (density, velocity in x, y and z, and thermal energy) as 5 separate one-dimension arrays which can be mapped to 3-dimensional arrays using A3D(ix,iy,iz) = A1D((iz-1)*(nyDim*nxDim) + (iy-1)*(nxDim) + ix), i.e., x-dimension changes most rapidly. We will also specify the grid dimensions in these cases.
As mentioned, for very high resolution grid realizations one should use the provided tool to produce an accurate grid realization from the SPH particle initial conditions. Plots produced in the above system of units are preferred, as this will make comparisons more straight forward.
We are meant to solve only the Euler equations so we will neglect physical viscosity in all the tests.In addition all tests should be performed as adiabatic. What we mean by adiabatic is what is the standard practice in the field of structure/galaxy formation. Therefore for adiabatic we mean the gas can undergo heating and cooling by adiabatic compression/decompression + irreversible heating by shocks (achieved through an artificial viscosity term in the thermal energy equation for the case SPH codes like ours). In other words we allow conversion of kinetic energy into thermal energy. Entropy will be generated in shocks so adiabatic here is not equivalent to isentropic.
Our code employs a spline kernel for the gravitational softening (the same type of function used in the calculation of the SPH smoothing length). The gravitational force between two particles becomes keplerian at twice the softening length (see for example the Appendix in Hernquist & Katz 1989 (ApJ, 70, 419) for a detailed description of the function). Therefore softening lengths quoted in this page are always given assuming a spline kernel. Slightly different lengths will have to be specified with different softening functions in order for the gravitational force resolution to be comparable (for example for a Plummer softening function the length should be about 30% lower than the spline).