Introduction

Density-based solvers are well suited to compute high-speed flows (transonic or supersonic). However, in cases where very low convective velocities dominate and the incompressibility limit is approached the large disparity arising between the smallest and largest eigenvalues of the system's characteristic matrix results in a progressive degradation of both the convergence properties and quality of results. These undesired effects can be already observed for Mach numbers below 0.2, independently of the implicitness of the solution algorithm. Pre-multiplication of the time derivative of the Euler or Navier-Stokes system of equations by a suitable matrix rescales the various characteristic velocities. Preconditioning does not only reduce the stiffness of the system but also significantly increase the accuracy at low flow speeds.

Face-flux description

The system of equations have been spatially discretised using a finite volume method. The fluxes at the cell faces are calculated by Roe's Flux-difference splitting methods [44] . Thus, the flux at face \(i+1/2\)

\begin{eqnarray} \mathbf{F}_{i+1/2} = \frac{1}{2}\left[\mathbf{F}\left(\mathbf{Q}^L_{i+1/2}\right) + \mathbf{F}\left(\mathbf{Q}^R_{i+1/2}\right) - |\mathbf{D}_Q|\left(\mathbf{Q}^L_{i+1/2} - \mathbf{Q}^R_{i+1/2}\right)\right], \end{eqnarray}

where \(\mathbf{Q}^L\) and \(\mathbf{Q}^R\) are the physical states at the left and right side of the considered faces. The matrix \(\mathbf{D}_Q\) contains the sum of the Roe-average flux Jacobians

\begin{eqnarray} \mathbf{D}_Q = \kappa_x \frac{\partial \mathbf{F}}{\partial \mathbf{Q}} + \kappa_y \frac{\partial \mathbf{G}}{\partial \mathbf{Q}} + \kappa_z \frac{\partial \mathbf{H}}{\partial \mathbf{Q}}, \end{eqnarray}

where \(\kappa_x,\kappa_y,\kappa_z\) are the directional components [44] . The modification of Eqn_RoeFlux after the preconditioning matrix has been applied yields:

"RoeFluxesLowMach"

\begin{eqnarray} \mathbf{F}_{i+1/2} = \frac{1}{2}\left[\mathbf{F}\left(\mathbf{Q}^L_{i+1/2}\right) + \mathbf{F}\left(\mathbf{Q}^R_{i+1/2}\right) - \mathbf{P}_Q^{-1}|\mathbf{P}_Q\mathbf{D}_Q| \left(\mathbf{Q}^L_{i+1/2} - \mathbf{Q}^R_{i+1/2}\right) \right] \end{eqnarray}

An essential aspect of the procedure is the determination of the term \(|\mathbf{P}_Q\mathbf{D}_Q|\). The complexity of this product is reduced by a variable transformation to primitive variables. Hence, using the decomposition

\begin{eqnarray} |\mathbf{P}_Q \mathbf{D}_Q| = \frac{\partial{\mathbf{Q}}}{\partial{\mathbf{U}}}|\mathbf{P}_U\mathbf{D}_U| \frac{\partial{\mathbf{U}}}{\partial{\mathbf{Q}}} = \frac{\partial{\mathbf{Q}}}{\partial{\mathbf{U}}}\mathbf{M}_U|\mathbf{\Lambda}|\mathbf{M}_U^{-1} \frac{\partial{\mathbf{U}}}{\partial{\mathbf{Q}}} \end{eqnarray}

the flux formulation RoeFluxesLowMach can be expressed as

"RoeFluxPreconditionedPrimitive"

\begin{eqnarray} \mathbf{F}_{j} = \frac{1}{2}\left[\mathbf{F}\left(\mathbf{Q}^L_{j}\right) + \mathbf{F}\left(\mathbf{Q}^R_{j}\right) - \mathbf{P}_Q^{-1} \frac{\partial{\mathbf{Q}}}{\partial{\mathbf{U}}}\mathbf{M}_U|\mathbf{\Lambda}|\mathbf{M}_U^{-1} \frac{\partial{\mathbf{U}}}{\partial{\mathbf{Q}}} \left(\mathbf{Q}^L_{j} - \mathbf{Q}^R_{j}\right) \right]_{j = i+1/2}. \end{eqnarray}

where \(\mathbf{\Lambda}\) is a diagonal matrix containing the eigenvalues of \(\mathbf{P}_U \mathbf{D}_U\):

\begin{eqnarray} \mathbf{\Lambda} = diag\left(\lambda_1,\lambda_2,\lambda_2,\lambda_4,\lambda_5\right). \end{eqnarray}

given by

\begin{eqnarray} \lambda_{1,2,3} &=& \mathbf{u}\cdot\kappa\\ \lambda_{4,5} &=& \frac {1}{2}\left((1 + \beta^2)\lambda_1 \pm \sqrt{\left(1 - \beta^2\right)^2 \lambda^2_1 + \beta^2 a^2 \|\kappa\|}\right) \end{eqnarray}

The matrix \(\mathbf{P}_U \mathbf{D}_U\) is given by

\begin{eqnarray} \mathbf{P}_U \mathbf{D}_U = \begin{pmatrix} \mathbf{u}\cdot\kappa & \kappa_x \rho \beta^2 & \kappa_y \rho \beta^2 & \kappa_z \rho \beta^2 & \mathbf{u}\cdot\kappa\frac{\beta^2 - 1}{a^2} \\[2mm] 0 & \mathbf{u}\cdot\kappa & 0 & 0 & \kappa_x \rho^{-1} \\[2mm] 0 & 0 & \mathbf{u}\cdot\kappa & 0 & \kappa_y \rho^{-1} \\[2mm] 0 & 0 & 0 & \mathbf{u}\cdot\kappa & \kappa_z \rho^{-1} \\[2mm] 0 & \beta^2 a^2 \rho \kappa_x & \beta^2 a^2 \rho \kappa_y & \beta^2 a^2 \rho \kappa_z & \beta^2 \mathbf{u}\cdot\kappa \end {pmatrix} \end{eqnarray}

where the preconditioning matrix reads

\begin{eqnarray} \mathbf{P}_U = \begin{pmatrix} 1 & 0 & 0 & 0 & \frac{\beta^2 - 1}{a^2} \\[2mm] 0 & 1 & 0 & 0 & 0 \\[2mm] 0 & 0 & 1 & 0 & 0 \\[2mm] 0 & 0 & 0 & 1 & 0 \\[2mm] 0 & 0 & 0 & 0 & \beta^2 \end {pmatrix}\qquad \mathbf{P}^{-1}_U = \begin{pmatrix} 1 & 0 & 0 & 0 & -\frac{\beta^2 -1}{a^2 \beta^2} \\[2mm] 0 & 1 & 0 & 0 & 0 \\[2mm] 0 & 0 & 1 & 0 & 0 \\[2mm] 0 & 0 & 0 & 1 & 0 \\[2mm] 0 & 0 & 0 & 0 & \frac{1}{\beta^2} \end {pmatrix} \end{eqnarray}

The left and right eigenvector matrices obtained using a scaling proposed by Hirsch read

\begin{eqnarray} \mathbf{M}^{-1}_U = \begin{pmatrix} \hat{\kappa}_x & 0 & \hat{\kappa}_z & -\hat{\kappa}_y & -\frac{\hat{\kappa}_x}{a^2}\\[2mm] \hat{\kappa}_y & -\hat{\kappa}_z & 0 & \hat{\kappa}_x & -\frac{\hat{\kappa}_y}{a^2}\\[2mm] \hat{\kappa}_z & \hat{\kappa}_y & -\hat{\kappa}_x & 0 & -\frac{\hat{\kappa}_z}{a^2}\\[2mm] 0 & \hat{\kappa}_x & \hat{\kappa}_y & \hat{\kappa}_z & -\frac{\hat{\lambda}_x - \hat{\lambda}_4}{\beta^2 a^2 \rho}\\[2mm] 0 & -\hat{\kappa}_x & -\hat{\kappa}_y & -\hat{\kappa}_z & \frac{\hat{\lambda}_1 -\hat{\lambda}_5}{\beta^2 a^2 \rho} \end{pmatrix}\qquad \mathbf{M}_{U} = \begin{pmatrix} \hat{\kappa}_x & \hat{\kappa}_y & \hat{\kappa}_z & \frac{\rho \beta^2}{\hat{\lambda}_4 - \hat{\lambda}_5} & \frac{\rho \beta^2}{\hat{\lambda}_4 - \hat{\lambda}_5}\\[2mm] 0 & -\hat{\kappa}_z & \hat{\kappa}_y & \hat{\kappa}_x \frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\hat{\lambda}_4 - \hat{\lambda}_5} & \hat{\kappa}_x \frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\hat{\lambda}_4 - \hat{\lambda}_5}\\[2mm] \hat{\kappa}_z & 0 & -\hat{\kappa}_x & \hat{\kappa}_y \frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\hat{\lambda}_4 - \hat{\lambda}_5} & \hat{\kappa}_y \frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\hat{\lambda}_4 - \hat{\lambda}_5} \\[2mm] -\hat{\kappa}_y & \hat{\kappa}_x & 0 & \hat{\kappa}_z \frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\hat{\lambda}_4 - \hat{\lambda}_5} & \hat{\kappa}_z \frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\hat{\lambda}_4 - \hat{\lambda}_5} \\[2mm] 0 & 0 & 0 & \frac{\beta^2 a^2 \rho}{\hat{\lambda}_4 - \hat{\lambda}_5} & \frac{\beta^2 a^2 \rho}{\hat{\lambda}_4 - \hat{\lambda}_5} \end{pmatrix} \end{eqnarray}

where the eigenvalues and metric components have been normalised by

\begin{eqnarray} \hat{\kappa}_m = \frac{\kappa_m}{\|\kappa\|},\quad m\in \{x,y,z\}\\[2mm] \hat{\lambda}_n = \frac{\lambda_n}{\|\kappa\|},\quad n\in \{1,2,3,4,5\} \end{eqnarray}

In terms of the conservative variables the right and left eigenvector matrices read

\begin{eqnarray} \scaletolinewidth{ \mathbf{M}_Q = \frac{\partial \mathbf{Q}}{\partial \mathbf{Q}_{U}} \mathbf{M}_U = \begin{pmatrix} \hat{\kappa}_x & \hat{\kappa}_y & \hat{\kappa}_z & \frac{\rho \beta^2}{\hat{\lambda}_4 - \hat{\lambda}_5} & \frac{\rho \beta^2}{\hat{\lambda}_4 - \hat{\lambda}_5}\\[3mm] u\hat{\kappa}_x & u\hat{\kappa}_y - \rho\hat{\kappa}_z & u\hat{\kappa}_z + \rho\hat{\kappa}_y & \frac{u\rho\beta^2 + \rho\hat{\kappa}_x\left(\hat{\lambda}_1 - \hat{\lambda}_5\right)}{\hat{\lambda}_4 - \hat{\lambda}_5} & \frac{u\rho\beta^2 + \rho\hat{\kappa}_x\left(\hat{\lambda}_1 - \hat{\lambda}_4\right)}{\hat{\lambda}_4 - \hat{\lambda}_5}\\[3mm] v\hat{\kappa}_x + \rho\hat{\kappa}_z & v\hat{\kappa}_y & v\hat{\kappa}_z - \rho\hat{\kappa}_x & \frac{v\rho\beta^2 + \rho\hat{\kappa}_y\left(\hat{\lambda}_1 - \hat{\lambda}_5\right)}{\hat{\lambda}_4 - \hat{\lambda}_5} & \frac{v\rho\beta^2 + \rho\hat{\kappa}_y\left(\hat{\lambda}_1 - \hat{\lambda}_4\right)}{\hat{\lambda}_4 - \hat{\lambda}_5}\\[3mm] w\hat{\kappa}_x - \rho\hat{\kappa}_y & w\hat{\kappa}_y + \rho\hat{\kappa}_x & w\hat{\kappa}_z & \frac{w\rho\beta^2 + \rho\hat{\kappa}_z\left(\hat{\lambda}_1 - \hat{\lambda}_5\right)}{\hat{\lambda}_4 - \hat{\lambda}_5} & \frac{w\rho\beta^2 + \rho\hat{\kappa}_z\left(\hat{\lambda}_1 - \hat{\lambda}_4\right)}{\hat{\lambda}_4 - \hat{\lambda}_5}\\[3mm] K_2\hat{\kappa}_x + \rho v \hat{\kappa}_z - \rho w \hat{\kappa}_y & K_2\hat{\kappa}_y + \rho u \hat{\kappa}_z - \rho w \hat{\kappa}_x & K_2\hat{\kappa}_z + \rho u \hat{\kappa}_y - \rho v \hat{\kappa}_x & \frac{K_2\rho\beta^2 + \hat{\lambda}_1\left(\hat{\lambda}_1 - \hat{\lambda}_5\right)}{\hat{\lambda}_4 - \hat{\lambda}_5} + \frac{\beta^2 a^2 \rho}{\left(\hat{\lambda}_4 - \hat{\lambda}_5\right) \left(\gamma - 1\right)} & \frac{K_2\rho\beta^2 + \hat{\lambda}_1\left(\hat{\lambda}_1 - \hat{\lambda}_4\right)}{\hat{\lambda}_4 - \hat{\lambda}_5} + \frac{\beta^2 a^2 \rho}{\left(\hat{\lambda}_4 - \hat{\lambda}_5\right)\left(\gamma - 1\right)}\\ \end{pmatrix} } \end{eqnarray}

\begin{eqnarray} \scaletolinewidth{ \mathbf{M}^{-1}_Q = \mathbf{M}^{-1}_U \frac{\partial \mathbf{Q}_{U}}{\partial \mathbf{Q}} = \begin{pmatrix} b\hat{\kappa}_x - \hat{\kappa}_z\frac{v}{\rho} + \hat{\kappa}_y\frac{w}{\rho} & \frac{u}{a^2}\left(\gamma - 1\right)\hat{\kappa}_x & \frac{\hat{\kappa}_z}{\rho} + \frac{v}{a^2}\left(\gamma - 1\right)\hat{\kappa}_x & -\frac{\hat{\kappa}_y}{\rho} + \frac{w}{a^2}\left(\gamma - 1\right)\hat{\kappa}_x & -\frac{\left(\gamma - 1\right)}{a^2}\hat{\kappa}_x\\[2mm] b\hat{\kappa}_y + \hat{\kappa}_z\frac{u}{\rho} - \hat{\kappa}_x\frac{w}{\rho} & -\frac{\hat{\kappa}_z}{\rho} + \frac{u}{a^2}\left(\gamma - 1\right)\hat{\kappa}_y & \frac{v}{a^2}\left(\gamma - 1\right)\hat{\kappa}_y & \frac{\hat{\kappa}_x}{\rho} + \frac{w}{a^2}\left(\gamma - 1\right)\hat{\kappa}_y & -\frac{\left(\gamma - 1\right)}{a^2}\hat{\kappa}_y \\[2mm] b\hat{\kappa}_z - \frac{u}{\rho}\hat{\kappa}_y + \frac{v}{\rho}\hat{\kappa}_x & \frac{\hat{\kappa}_y}{\rho} + \frac{u}{a^2}\left(\gamma - 1\right)\hat{\kappa}_z & -\frac{\hat{\kappa}_x}{\rho} + \frac{v}{a^2}\left(\gamma - 1\right)\hat{\kappa}_z & \frac{w}{a^2}\left(\gamma - 1\right)\hat{\kappa}_z & -\frac{\left(\gamma - 1\right)}{a^2}\hat{\kappa}_z \\[2mm] -\frac{\hat{\lambda}_1}{\rho} - K_3\frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\beta^2a^2\rho} & \frac{\hat{\kappa}_x}{\rho} + u\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\beta^2a^2\rho} & \frac{\hat{\kappa}_y}{\rho} + v\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\beta^2a^2\rho} & \frac{\hat{\kappa}_z}{\rho} + w\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\beta^2a^2\rho} & -\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_4}{\beta^2a^2\rho}\\[2mm] \frac{\hat{\lambda}_1}{\rho} + K_3\frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\beta^2a^2\rho} & -\frac{\hat{\kappa}_x}{\rho} - u\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\beta^2a^2\rho} & -\frac{\hat{\kappa}_y}{\rho} - v\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\beta^2a^2\rho} & -\frac{\hat{\kappa}_z}{\rho} - w\left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\beta^2a^2\rho} & \left(\gamma - 1\right)\frac{\hat{\lambda}_1 - \hat{\lambda}_5}{\beta^2a^2\rho}\\[2mm] \end{pmatrix} } \end{eqnarray}

The variables \(K_2\) and \(K_3\) are defined by

\begin{eqnarray} K_2 &:=& \left(\frac{q^2}{2} - \frac{\omega^2_R}{2}\right)\\ K_3 &:=& \left(\gamma - 1\right)\left(\frac{q^2}{2} + \frac{\omega^2_R}{2}\right)\\ b &:=& \left(1 - \frac{K_3}{a^2}\right)\\ \end{eqnarray}