What I'm looking for:

Let $\vec{W}$ be the vector of conserved variables for a 1-dimensional, adiabatic, (special) relativistic, electrically neutral fluid. (Yes, something that simple!) I'm looking for a paper that derives the form of the matrix $A$ that linearizes the evolution equation. That is, $$ \partial_t \vec{W} \approx A \partial_x \vec{W}. $$ Alternatively, the matrix $B$ that does the same for the conserved variables will work: $$ \partial_t \vec{U} \approx B \partial_x \vec{U}. $$

These matrices are useful in fluid computations for several reasons. I need more than just the eigenvalues (used in certain solvers) - I also need the eigenfunctions. I've found plenty of papers that treat the nonrelativistic case (for one of many, many examples, see the appendices of Stone et al. 2008, ApJS 178 137) both with and without magnetism. I've found some papers that just quote a few eigenfunctions for relativistic MHD, but these are often the ones that are interesting only with magnetism in play.

I'm looking for whatever paper derives these matrices in the rather simple case I'm dealing with. An answer that gives the derivation would be nice, but I'm also trying to locate the relevant literature. In particular, I would like a paper addresses the physical reliability/usefulness of the linear approximation.

Background:

There are three primitive variables defining my fluid: rest-mass density $\rho$, velocity $v$, and pressure $p$. By convention, these variables are combined into a vector $\vec{W} = (\rho, v, p)^\mathrm{T}$.

Many approaches to evolving fluids deal with the equations in flux-conservative form. Here, the conserved variables are \begin{align} D & = \gamma\rho && \text{(lab-frame density),} \\ M & = Dh\gamma v && \text{(relativistic momentum),} \\ E & = Dh\gamma - p && \text{(relativistic energy).} \end{align} Here I define \begin{align} \gamma & = \frac{1}{\sqrt{1-v^2}} && \text{(standard Lorentz factor),} \\ h & = 1 + \frac{\Gamma}{\Gamma-1} \left(\frac{p}{\rho}\right) && \text{(enthalpy),} \end{align} where the ratio of specific heats $\Gamma$ is assumed to be constant. These are often combined as $\vec{U} = (D, M, E)^\mathrm{T}$.

Along with the vector of conserved quantities, we can define the vector of fluxes $$ \vec{F} = \begin{pmatrix} Dv \\ Mv + p \\ M \end{pmatrix}. $$ Then we have the relation $$ \partial_t \vec{U} + \partial_x \vec{F} = 0, $$ which is used as the basis for most Riemann solvers and many fluid codes in general.