Many numerical methods for hyperbolic PDEs are based on the use of Riemann solvers. Such solvers are essential for accurately capturing shock waves. There are a range of such solvers available for the most well-studied systems (e.g., exact solvers, Roe solvers, HLL solvers). How should I decide which to use?
For the numerical solution of hyperbolic PDEs the use of Riemann solvers are essential components of conservative shock capturing methods for accurate simulation of wave problems which may have shocks (discontinuous jumps in conserved variables). To be able to obtain accurate solutions to such problems, we need to use proper upwinding techniques -- the Riemann solver is responsible for this. The Riemann solver seek an accurate solution to the interface problem between cells (fx. in Finite Volumes) or elements (fx. in Discontinuous Galerkin Finite Element Methods). The solution of this interface problem is based on solution of either side of the interface and seeks to use this as basis for accurate reconstruction of the (numerical) flux (in terms of conserved variables) across the interface. Goal is to mimic the underlying physics as accurately as possible/necessary to obtain accurate solutions.
There are two standard approaches to the solution of such (local to the interface) Riemann problems, namely, exact and approximate Riemann solvers. For many PDEs there are no exact closed-form solution available in which case we have to resort to approximate Riemann solvers. In practice, it may also be (too) expensive to exactly solve the Riemann problems in which case it may be more practical to resort to approximate Riemann solvers. For the same reason Lax-Freidrichs type fluxes are widely used as a simple means.
Essentially, the choice between Riemann solvers has to do with how accurately one seeks to take represent the wave speeds of the solution and the resulting efficiency.
It is problem dependent. The Riemann problem is based on data from either side of cell interfaces. To reconstruct the flux at the interface based on this data we must know information about the full wave structure of the hyperbolic PDE in question. This makes the Riemann problem problem-dependent and therefore also the choice of Riemann solver. In short, exact solvers seek to take into account the full wave structure, the Roe solver is based on local approximation (by linearization and special averaging) of the local wave structure, the HLL solver is based on estimating two dominant wave speeds in the local wave structure and then impose conservation by satisfying Rankine-Hugoniot condition to hold across shocks or contact discontinuities.
Thus, the choice between specific solvers, exact solvers or approximate Roe/HLL/etc solvers depends on striking a balance between accuracy (in mimicking the underlying physics of the model equations) and efficiency needs. In the end -- as I see it -- in practical application it is often efficiency requirements that dictates the use of approximate Riemann solvers (fx. of the Lax-Friedrichs type).
A good exposition to the subject is given by E. F. Toro in his textbook "Riemann solvers and numerical methods for fluid dynamics", Springer.
I've been lead to believe that for low order numerics one needs high quality Riemann solvers and for high order numerics one can use low quality Riemann solvers. Intuitively there's some number of FLOPs necessary to capture the physics below which one is hosed.
And, yes, there's also zero content in this answer from an evaluation metric perspective...