Riemann Problem for Linear Hyperbolic Systems

In this post, I will be solving the linear system

\displaystyle \frac{\partial \vec{q}}{\partial t} + A \frac{\partial \vec{q}}{\partial x} =\vec{0}

with the Riemann initial conditions

\displaystyle \vec{q}(x,0) = \left\{ \begin{array}{lr} \displaystyle \vec{q}_l &, x <0 \\ \displaystyle \vec{q}_r &, x>0\end{array}

where x \in \mathbb{R} and A is an n by n matrix which has distinct real eigenvalues.

Since A has real distinct eigenvalues, we know that the matrix can be decomposed as

\displaystyle A = R \Lambda R^{-1}

where columns in R correspond to A‘s eigenvectors and \Lambda is a diagonal matrix containing the eigenvalues

\displaystyle R = \left [ r_1 | r_2 | r_3|\cdots|r_n \right]
\displaystyle \Lambda = \left( \begin{array}{ccccc} \lambda_1 &0 &\ldots & \ldots & 0 \\ 0 & \lambda_2 &0 &\ldots & 0\\ & & \ddots & & \\ & & & \ddots & \\ 0 & \ldots & \ldots & 0 & \lambda_n \end{array} \right)

Additionally, because we have distinct real eigenvalues, we know that R is a unitary matrix (R^*R=I).

Because R is unitary, we can uniquely decompose a vector into

\displaystyle \vec{u} = \sum_{p=1}^n ( \vec{r}_p^* \vec{u} ) \vec{r}_p

Using this fact, it is beneficial to write \vec{q}_l and \vec{q}_r as

\displaystyle \vec{q}_l = \sum_{p=1}^n \alpha_p \vec{r}_p
\displaystyle \vec{q}_r = \sum_{p=1}^n \beta_p \vec{r}_p

where \alpha_p = \vec{r}_p^* \vec{q}_l and \beta_p = \vec{r}_p^* \vec{q}_r .

Now, we will transform the original pde using \vec{v} = R^{-1} \vec{q} to

\displaystyle R^{-1} \left ( \frac{\partial \vec{q}}{\partial t} + R \Lambda R^{-1} \frac{\partial \vec{q}}{\partial x} \right ) =\vec{0}
\displaystyle \frac{\partial \vec{v}}{\partial t} + \Lambda \frac{\partial \vec{v}}{\partial x} \right ) =\vec{0}

Likewise, we can rewrite the initial conditions as

\displaystyle \vec{v}(x,0) = R^{-1} \vec{q}(x,0)
\displaystyle = \left\{ \begin{array}{lr} \displaystyle R^{-1} \sum_{p=1}^n \alpha_p \vec{r}_p &, x <0 \\ \displaystyle R^{-1} \sum_{p=1}^n \beta_p \vec{r}_p &, x>0\end{array}
\displaystyle = \left\{ \begin{array}{lr} \displaystyle \sum_{p=1}^n \alpha_p \vec{e}_p &, x <0 \\ \displaystyle \sum_{p=1}^n \beta_p \vec{e}_p &, x>0\end{array}

where e_p is a unit vector.

We can conclude that the original system decouples to n linear advection equations

\displaystyle \frac{\partial \vec{v}^{(p)}}{\partial t} + \lambda_p \frac{\partial \vec{v}^{(p)}}{\partial x} =\vec{0}
\displaystyle \vec{v}^{(p)}(x,0) = \left\{ \begin{array}{lr} \displaystyle \alpha_p &, x <0 \\ \displaystyle \beta_p &, x>0\end{array}

which has the solution

\displaystyle \vec{v}^{(p)}(x,t) = \vec{v}^{(p)}(x-\lambda_p t,0) = \left\{ \begin{array}{lr} \displaystyle \alpha_p &, x/t < \lambda_p \\ \displaystyle \beta_p &, x/t >\lambda_p \end{array}

Therefore, we can solve the original equation by using the transformation \vec{q} = R \vec{v} which gives us

\displaystyle \vec{q}(x,t) = R \vec{v}(x,t) = \sum_{p=1}^n \vec{v}^{(p)}(x,t) \vec{r}_p
\displaystyle = \sum_{\lambda_p > x/t} \alpha_p \vec{r}_p + \sum_{\lambda_p < x/t} \beta_p \vec{r}_p