Numerical Relativistic Hydrodynamics

1. i.e. The Cauchy problem from the viewpoint of an (Astro)Physicist… Numerical Relativistic Hydrodynamics Luciano Rezzolla SISSA, International School for Advanced Studies, Trieste INFN, Department of Physics, University of Trieste Cargese, 06/08/02

2. Plan of the Talk • A look at the eqs: when matter is not “on the way”… • Valencia Formulation of GR Hydrodynamics • The discretization problem • Numerical Methods • Artificial Viscosity Methods • The Riemann problem and HRSC Methods • Godunov and his idea (basic physics  num. methods) • A quick primer: reconstruct, solve & build… • Conclusions and lookout

3. A look at the equations: 3+1, Tmn and all of that… Foliate the spacetime in t=const spatial hypersurfaces St Let n be the unit timelike4-vectororthogonaltoSt such that We then fill the spacetime with a generic fluid such that: u: fluid’s 4-velocity, p:isotropic pressure, r : rest-mass density, e : specific internal energy density, e=r( 1+e ): energy density

4. Once the fluid variables have been specified, a stress energy tensor can be constructed non-ideal part In practice, it’s more convenient to use 3-velocities vectors v which can be readily defined once the ”normal observers”n have been introduced in the spacetime foliations so that

5. The set of equations to solve in 3D and without symmetries is In addition, we need to specify suitable gauge conditions for the lapse a and the components of the shift vector b i To best hard-wire the conservative nature of the equations, they can be written in a first-order, “flux-conservative form”.

6. Replace the ”primitive variables”,in terms of “conserved variables”: As a result: where and Valencia formulation where is the vector of conserved variables “fluxes” “sources”

7. Discretizating the problem… Let’s restrict to a simpler but instructive problem: a homogeneous, flux- conservative differential equation for the scalaru=u(x,t) in one dimension Its generic, finite-difference form is (1st-order in time) where and “some approximation to the average flux at j” Any finite-difference form of(1)must representin the most accurate way.Different forms of calculating lead to different evolution schemes (Forward-Time-Centred-Space, Lax, Runge-Kutta, etc…)

8. A generic problem arises when a Cauchy problem described by a set of continuous PDEs is solved in a discretized form: the numerical solution is, at best,piecewise constant. This is particularly problematic when discretizing hydrodynamical eqs in compressible fluids, whose nonlinear properties generically produce, in a finite time, nonlinear waves with discontinuities (ie shocks, rarefaction waves, etc) even from smooth initial data!

9. Because the occurrence of discontinuities is a fundamental property of the hydrodynamical equations, any finite-difference scheme must be able to handle them in a satisfactory way. • Possible solutions to the discontinuities problem: • 1st order accurate schemes • generally fine, but very inaccurate across discontinuities (eccessive diffusion, e.g. Lax method) • 2nd order accurate schemes • generally introduce oscillations across discontinuities • 2nd order accurate schemes with artificial viscosity • mimic Nature (more later) • Godunov Methods • discontinuities are not eliminated, rather they are exploited! (more later)

10. Artificial Viscosity Methods (von Neumann & Richtmyer 1950) They are largely used methods that introduce an “artificial” viscosityQsuch that the spurious oscillations across the shocks are damped. They resemble Nature in that the kinetic energy build-up on the largest wavenumber modes are converted into internal energy. These additional terms are introduced only in the vicinity of the shocks where and k is an adjustable coefficient • PROS • easy to implement • computationally efficients • CONS • problem dependent • inaccurate for ultrarelativ. flows

11. High Resolution Shock Capturing Methods Based on a simple, yet brilliant idea by Godunov (’59). An example of how basic physics can boost research in computational physics. Core idea: a piecewise constant description of hydrodynamical quantities will produce a collection of local Riemann problems whose solution can be found exactly. where is the exact solution of the Riemann problem with initial data “left” state “right” state

12. What is exactly a Riemann problem?… It’s the evolution of a fluid initially composed of two states with different and constant values of velocity, pressure and density. If the problem is linear, it can be handled analytically after rewriting the flux conservative equation as whereA(u)is the Jacobian matrix of const. coefficients. In this way, (2) is written as a set of i linear equations for the characteristic variables withLthe diagonal matrix of the eigenvaluesli . The solution is andare the right eigenvectors ofA

13. rarefaction wave shock front Solution at the time n+1 of the two Riemann problems at the cell boundaries xj+1/2 and xj-1/2 Spacetime evolution of the two Riemann problems at the cell boundaries xj+1/2 and xj-1/2.Each problem leads to a shock wave and a rarefaction wave moving in opposite directions Initial data at the time n for the two Riemann problems at the cell boundaries xj+1/2 and xj-1/2 cell boundaries where fluxes are required

14. A quick primer.The numerical solution of a Riemann problems is based on three basic steps:first reconstructthen solvefinally convert and build…

15. A quick primer: first reconstruct Higher accuracy is reached with a better representation of the solution. ”Reconstructing” the initial data for the Riemann problem at the cell boundaries can be made with a number of algorithms; the most interesting are the TVDs(minmod, MC, Superbee) for which the solutionu=u(x)is linear reconst. Wherefis a coefficient based on the slope of uand varying from0(near a discontinuity) to1(in smooth regions of the solution). slope limiter

16. then solve… While in principle one could solve “exactly” the Riemann problem at each cell interface, this is almost never done in practice. Rather “approximate”Riemann solvers are used in place of the exact ones and a number of them are available (Roe, Marquina, HLLE, etc…). As an example, Roe’s approximate Riemann solver can be calculated as where wR and wL are the values of the primitive variables at the right/left sides of the i-th interface and are the eigenvalues and right eigenvectors of the Jacobian matrix. The coefficients measure the jumps of the characteristic variables across the characteristic field

17. Finally, average, convert and build… Once the solution in terms of the conserved variables D, Sj and E has been obtained, it is necessary to return to the primitive variables after inverting numerically the set of equations Once the primitive variables have been calculated, the stress-energy tensor can be reconstructed and used on the right hand side of the Einstein equations. This is repeated at each grid point and for each time level…

18. Conclusions and lookout • We now have a good (but improvable!) mathematical and numerical framework for highly accurate, hydrodynamical simulations in GR • The use of HRSC methods provides us with the opportunity/privilege of tapping numerical experties from a much wider community (CFD) • New “ideas” are making this framework even more robust and efficient • A wish-list: i.e. what would be nice to have and how mathematicians can be of great help: • Formulations of EFEs better suited to numerical calculations • More efficient mathematical framework for the solution of the Riemann problem • More efficient mathematical framework for relativistic non-ideal hydrodynamics

19. New in relativistic hydrodynamics • Predicting the “waves” (LR & O. Zanotti, JFM, 02) Uses the invariant expression for the relative normal velocity vx12between the two initial states of a RP and predicts the wave-pattern produced. limiting values vx12is a monotonic function of the pressure across the contact discontinuityp3and assumes a different form according to the wave pattern produced (2S, SR, 2R). For each initial data, the comparison of vx12with the relevant limits predicts the wave pattern that will be produced. This allows for a simple logic, hence a more compact and efficient algorithm. 2 shocks branch 1 shock, 1 rw 2 rw’s

20. Riemann Problem in Multidimensional flows (Pons et al., JFM 01) The special relativistic exact Riemann problem is investigated in multidimensional flows. It appears evident that the tangential velocities, introduce major differences and couple quantities on either side of the nonlinear waves In relativistic hydrodynamics the tangential velocity vector is discontinuos across a shock while preserving the direction This is of great importance for the implementation of exact Riemann solvers in multidimensional codes

21. New relativisitic hydrodynamical effects (LR & Zanotti, 02) This is an example in which: num. methods  basic physics The coupling of the states on either side of a nonlinear waves via the tangentialvelocities gives rise to new relativistic effects. While maintaining the initial conditions (p, r, vx)not altered, a change in the tangential velocities produces a shift from one wave-pattern to another.

22. Consider Sod’s problem: r1=1.0, p1=1.0, r2=0.1, p2=0.125 zero tangential velocities nonzero tangential velocities vt1=0.0 vt2=0.9 vx1=0.5, vx2=0.0 vt1=0.9 vt2=0.0 vx1=0.0, vx2=0.5