The Ohio State University viscous hydrodynamics code for relativistic heavy-ion collisions
I am not the author of this software: it was originally written by Dr. Huichao Song and further modified and enhanced by Zhi Qiu, Chun Shen, and Jia Liu. Alternate versions may be found at https://github.com/chunshen1987/VISHNew and https://github.com/JiaLiu11/iEBE/tree/hydro/EBE-Node/VISHNew.
Some references:
This is a stripped down and cleaned up version of the Ohio State University viscous 2+1D hydrodynamics code (also known elsewhere as "VISHNew").
I started from Jia Liu's version including shear and bulk viscosity at https://github.com/JiaLiu11/iEBE/tree/hydro/EBE-Node/VISHNew, however this repository does not share git history with Jia's because I used git filter-branch to remove a number of large files from old commits.
I made this version with one task in mind: running hydro as part of an event-by-event hybrid model, i.e. to be preceded by a Monte Carlo initial condition model and followed by particlization and hadronic transport models.
To compile and install, follow the standard CMake sequence:
mkdir build && cd build
cmake .. [-DCMAKE_INSTALL_PREFIX=<prefix>]
make install
This will place the compiled binary osu-hydro in <prefix>/bin/ and the configuration file osu-hydro.conf and equation of state data eos.dat in <prefix>/share/osu-hydro/.
At runtime, osu-hydro searches for osu-hydro.conf and eos.dat in the following locations:
- The current directory.
- User data directory
$XDG_DATA_HOME/osu-hydro/, with$XDG_DATA_HOMEdefaulting to~/.local/share/if unset. - System directory
/usr/local/share/osu-hydro/. - System directory
/usr/share/osu-hydro/.
If either of the files cannot be found, the program will exit with an error message.
Thus, osu-hydro may be executed from anywhere provided it can locate the data files.
However, it always expects the initial condition data files (described below) to be in the current working directory, and it will also write output files in the current directory.
The EOS is generated at build time by the script eos/eos.py, which connects a hadron resonance gas EOS at low temperature to a lattice EOS (HotQCD, http://inspirehep.net/record/1307761) at high temperature.
See the help output (eos.py --help) for more options.
This depends on frzout to compute the hadron resonance gas part.
The config file osu-hydro.conf has brief comments for each parameter.
Note that the code expects the parameters in a specific order, so do not change the order of parameters in the config file.
Parameters may also be passed on the command line with key=value syntax (case-insensitive), e.g.
osu-hydro edec=0.30
The shear viscosity is parametrized as the modified linear ansatz
(eta/s)(T) = min + slope*(T - T0)*(T/T0)^crv
for T > T0, where T0 is the temperature at which eta/s reaches its minimum value.
Previously, T0 was fixed to the HotQCD EoS transition temperature Tc = 0.154 GeV, and that is now its default value.
The parameter min is the minimum value at temperature T0, slope is the slope above T0, and crv controls the curvature above T0.
The min and slope must be non-negative; crv may be negative but probably not below -1, since this would cause decreasing (eta/s)(T).
For T < T0 (HRG phase), eta/s is constant, set by the parameter etas_hrg.
The eta/s parameters may be set in the config file or on the command line with keys etas_t0, etas_hrg, etas_min, etas_slope, etas_crv.
The bulk viscosity is parametrized as a Cauchy distribution with tunable peak location, max value, and width:
(zeta/s)(T) = max / (1 + ((T - T0)/width)^2)
The location (temperature T0), max, and width may be set in the config file or on the command line with keys zetas_t0, zetas_max, zetas_width.
The computational grid is always square and centered at the origin.
The size may be set in the config file or on the command line by the nls parameter.
LS ("lattice size") means the number of cells from the origin to the grid edge, not counting the cell precisely at the origin; the full grid is (2*LS + 1) × (2*LS + 1).
The grid cell size is set in the config file or on the command line by the dxy parameter.
Ensure the cell size is small enough to resolve the system geometry.
For Monte Carlo initial conditions with Gaussian nucleons of width w, dxy = 0.15*w is a good compromise between computation time and precision.
Example: With the default DXY = 0.1 fm and LS = 200, the grid is 401 × 401 with cells centered at -20.0, -19.9, ..., 0.0, 0.1, ..., 20.0 fm.
WARNING:
The grid size is defined by the center of the outermost grid cell, which is inconsistent with programs that use the cell edge, such as trento and freestream.
The two definitions differ by half a cell width.
For example, the default 401 × 401 grid has its outermost cells centered at ±20.0 fm, and the cells are 0.1 fm wide, so the physical extent of the grid is a half-cell larger: ±20.05 fm.
To generate a matching grid with trento, use options --grid-max 20.05 --grid-step 0.1.
The timestep is set in the config file or on the command line by the dt parameter.
The SHASTA algorithm used in this code requires the timestep to be less than one-half the spatial step.
However, a smaller timestep around one-quarter to one-fifth the spatial step is a good compromise between computation time and precision.
By default (with option InitialURead = 1), the initial energy density, flow, and shear viscous tensor must be provided.
Create the following binary data files:
ed.dat- energy densityu{1,2}.dat- flow velocity in x and y directionspi{11,12,22}.dat- shear tensor components xx, xy, yy (the remaining components are computed internally)- The bulk pressure is computed internally as the deviation from ideal pressure, Π = e/3 - p(e).
If InitialURead = 0, only the energy density is required.
Alternatively, set IEin = 1 and provide the entropy density sd.dat.
All files must contain an array of doubles in binary format.
They can be easily created in Python using numpy.ndarray.tofile, e.g., if ed is a numpy array containing the energy density, then ed.tofile('ed.dat') writes it as a binary file.
WARNING: The input data files must have precisely the expected grid size, as described above in grid settings. This is not checked by the code.
During time evolution, the code prints out a line for each timestep
timestep_number tau max_energy_density max_temperature
Evolution stops after the max energy density drops below the configured decoupling energy density.
The only output file is the hypersurface data surface.dat, a binary file containing an N × 16 array of doubles, where N is the number of freeze-out volume elements and the 16 columns are
| Index | Quantity | Index | Viscous pressures [GeV/fm3] |
|---|---|---|---|
| 0 | τ [fm] | 8 | πtt |
| 1 | x [fm] | 9 | πtx |
| 2 | y [fm] | 10 | πty |
| 3 | dσt [fm2] | 11 | πxx |
| 4 | dσx [fm2] | 12 | πxy |
| 5 | dσy [fm2] | 13 | πyy |
| 6 | vx | 14 | πzz |
| 7 | vy | 15 | Π |
The file can be read in Python with
surface = np.fromfile('surface.dat', dtype='float64').reshape(-1, 16)Note: The surface normal vectors dσμ are output in covariant form, which is standard. From 84fc226 (July 2016) until 9995ee1 (April 2017), they were contravariant (dσμ).