math]†‡§¶∥††‡‡

# Artificial diffusion for convective and acoustic low Mach number flows II:
Application to Liou-Steffen, Zha-Bilgen and Toro-Vasquez convection-pressure
flux splittings

Joshua Hope-Collins***Oxford Thermofluids Institute, University of Oxford,
Osney Mead, Oxford, OX2 0ES, England, United Kingdom †††Department of
Mathematics, Imperial College London, London, SW7 2AZ, England, United Kingdom
‡‡‡Corresponding author. Email address<EMAIL_ADDRESS>
Luca di Mare†

(August, 2024)

###### Abstract

Liou-Steffen splitting (AUSM) schemes are popular for low Mach number
simulations, however, like many numerical schemes for compressible flow they
require careful modification to accurately resolve convective features in this
regime. Previous analyses of these schemes usually focus only on a single
discrete scheme at the convective limit, only considering flow with acoustic
effects empirically, if at all. In our recent paper _Hope-Collins & di Mare,
2023_ we derived constraints on the artificial diffusion scaling of low Mach
number schemes for flows both with and without acoustic effects, and applied
this analysis to Roe-type finite-volume schemes. In this paper we form
approximate diffusion matrices for the Liou-Steffen splitting, as well as the
closely related Zha-Bilgen and Toro-Vasquez splittings. We use the constraints
found in _Hope-Collins & di Mare, 2023_ to derive and analyse the required
scaling of each splitting at low Mach number. By transforming the diffusion
matrices to the entropy variables we can identify erroneous diffusion terms
compared to the ideal form used in _Hope-Collins & di Mare, 2023_. These terms
vanish asymptotically for the Liou-Steffen splitting, but result in spurious
entropy generation for the Zha-Bilgen and Toro-Vasquez splittings unless a
particular form of the interface pressure is used. Numerical examples for
acoustic and convective flow verify the results of the analysis, and show the
importance of considering the resolution of the entropy field when assessing
schemes of this type.

Keywords: Low Mach number flows; Computational fluid dynamics; Asymptotic
analysis; Numerical dissipation; Flux vector splittings; AUSM;

Highlights:

  * •

Formulation of approximate diffusion matrices for convection-pressure flux
splittings

  * •

Identification of required asymptotic scaling of diffusion coefficients

  * •

Considers both convective and acoustic flow features

  * •

Resolution of the entropy field is particularly sensitive to the choice of
splitting

  * •

Numerical results show excellent agreement with analysis

###### Contents

  1. 1 Low Mach number limits of the Euler equations
  2. 2 Artificial diffusion for low Mach number schemes
    1. 2.1 Entropy variables form
    2. 2.2 Conserved variables form
  3. 3 Convection-pressure flux-vector splittings at low Mach number
    1. 3.1 Conserved variables form of convection-pressure flux-vector splittings
      1. 3.1.1 Liou-Steffen
      2. 3.1.2 Zha-Bilgen
      3. 3.1.3 Toro-Vasquez
    2. 3.2 Entropy variables form
    3. 3.3 Limit equations of flux-vector splittings
      1. 3.3.1 Roe form
      2. 3.3.2 Liou-Steffen splitting
      3. 3.3.3 Zha-Bilgen splitting
      4. 3.3.4 Toro-Vasquez splitting
  4. 4 Comparison of existing schemes
  5. 5 Numerical examples
    1. 5.1 One dimensional examples
      1. 5.1.1 Isolated soundwave
    2. 5.2 Two dimensional examples
      1. 5.2.1 Circular cylinder
      2. 5.2.2 Isentropic vortex
  6. 6 Conclusions

## Introduction

Low Mach number flows occur in many applications of engineering and scientific
interest, and the design of numerical schemes that can accurately simulate
these flows has been an active research area for several decades [61, 18]. It
is well-known that collocated density-based schemes for compressible flow may
fail to accurately resolve convective flow features at low Mach number due to
the incorrect scaling of the artificial diffusion terms in the asymptotic
limit $M\to 0$ [20, 60]. A great many schemes which overcome this failing have
been developed, including central-difference schemes [62, 61, 65] and flux-
difference splitting schemes [20, 19, 57, 46], as well as convection-pressure
flux-vector splitting schemes [36, 32, 54] which are the focus of the current
work.

Convection-pressure flux-vector splittings separate the physical flux of the
Euler equations into two terms: one describing the transport of some set of
convected quantities, and the other containing the contribution of the
pressure. This is in contrast with other flux-vector splitting methods where
the flux vector is split into a forward travelling wave and a backward
travelling wave [55, 28]. A numerical flux function can be constructed from a
convection-pressure flux-vector splitting by finding suitable discretisations
for each term, and fluxes of this type are attractive for a number of reasons.
Separating the flux allows different modelling approaches to be used for each
term, which is especially useful at low Mach number where the velocity and
pressure play increasingly distinct roles [42]. Splitting in this manner also
gives a balance of accuracy and computational efficiency: intermediate contact
waves can be accurately resolved unlike standard flux-vector splitting
schemes, but the computational expense scales with the size $n$ of the system
being solved instead of $n^{2}$ for flux-difference splittings, which is
especially useful for multi-component flows [37, 33].

The most commonly used convection-pressure flux-vector splitting is the Liou-
Steffen splitting, named advection upstream splitting method (AUSM) in the
original papers from Liou & Steffen 1993 and Liou 1996 [37, 33]. This
splitting has been widely applied, including to low Mach number [32],
hypersonic [25], multicomponent and multiphase [43, 26] flows. The earliest
low Mach number schemes were inspired by the ideas of low Mach number
preconditioning [14, 36, 13], which can be combined with flux-difference
splitting schemes to overcome the accuracy problem at low Mach number [60,
20]. Liou 2006 [32] integrated these ideas into the AUSM+-up for all speeds
scheme, which is still in widespread use today. Of the more recent Liou-
Steffen splitting schemes suitable for low Mach number, SLAU and its
descendants are particularly popular [54, 27, 25]. Two other convection-
pressure flux-vector splittings are the Zha-Bilgen [69] and Toro-Vasquez [59]
splittings, which differ from the Liou-Steffen splitting only in the treatment
of the energy equation. The Zha-Bilgen splitting was derived by splitting the
exact flux Jacobian into convective and acoustic components in [49, 55]. Zha &
Bilgen 1993 [69] developed a numerical scheme for this splitting, and as such
it is usually referred to by their names. The Toro-Vasquez splitting was
proposed by Toro & Vasquez 2012 [59] in the context of developing a general
framework for the discretisation of convection-pressure flux-vector
splittings. Over the last decade, the Zha-Bilgen and Toro-Vasquez splittings
have gained some attention for low Mach number flow [45, 4, 22, 56, 31, 7],
and it has been shown that accurate schemes for this regime can be produced
from these splittings. We cover more detail on the literature for low Mach
number schemes using convection-pressure flux-vector splittings in section 4.

In our previous paper _“Artificial diffusion for convective and acoustic low
Mach number flows I: Analysis of the modified equations, and application to
Roe-type schemes”_ [21] we analysed the form and behaviour of artificial
diffusion in numerical schemes at low Mach number. Three diffusion scalings
naturally arise - one for purely convective flow, one for purely acoustic
flow, and one for flow with either or both convective and acoustic features.
Our approach used the modified equations, which make the artificial diffusion
terms of the discrete numerical scheme explicit in the continuous PDE and
enable the behaviour of these terms to be studied independently of the
specific discretisation they arose from. The modified equations have
previously been used to study low Mach number numerical schemes in [60, 11].
By applying this analysis to Roe-type schemes, we showed that many existing
low Mach number schemes from the literature adhere to one of the three
diffusion scalings, and display the expected behaviour. While the behaviours
of these schemes have been shown and studied previously, this approach
demonstrated that many properties of low Mach number schemes can be predicted
and understood using the continuous equations, and that considering both
convective and acoustic flow is necessary to provide a complete picture of
this behaviour.

In this paper, we apply this approach to investigate the behaviour of
convection-pressure flux-vector splitting schemes at low Mach number,
specifically the Liou-Steffen, Zha-Bilgen, and Toro-Vasquez splittings. By
forming explicit diffusion matrices for each splitting, we can clearly see the
differences in the forms of the diffusion, and by transforming to the entropy
variables we can make direct comparison to the results of [21]. Many previous
studies of low Mach number numerical schemes have focused on Roe-type schemes
- as we did in [21] \- and the findings are extended to convection-pressure
flux-vector splitting schemes by identifying the diffusion terms that these
schemes have in common with the Roe-type schemes [11, 48, 30]. While we also
do this here, working in the entropy variables means that we can additionally
show the effect of the diffusion terms which are different between the
convection-pressure flux-vector splitting schemes and the Roe-type schemes.
The particular form of these terms turns out to be crucial to the success of
these schemes at low Mach number.

The structure of the remainder of the paper is as follows. In section 1 we
very briefly cover the necessary background for the low Mach number limits of
the Euler equations, and in section 2 we recap some of the key findings of
[21] which will be used later. Section 3 contains the main content of this
study. The approximate diffusion matrix for each splitting is shown in the
conserved variables, then transformed to the entropy variables. The form in
the entropy variables is discussed, and the limiting form of the modified
equations for each splitting at both the convective and acoustic limits are
found. In section 4 we look at a number of existing schemes from the
literature, compare them to the generic forms used in section 3, and identify
which of the three diffusion scalings they use. The results of several
numerical examples are shown in section 5 and display excellent agreement with
the analysis of section 3 for all three splittings for both convective and
acoustic flow.

##  1 Low Mach number limits of the Euler equations

Here we give a brief recap of the necessary background for the low Mach number
asymptotic limits of the Euler equations. More details are given in [42, 11,
21]. We use the non-dimensionalisation:

$\begin{gathered}\rho=\frac{\tilde{\rho}}{\rho_{\infty}},\quad\underline{u}=\frac{\underline{\tilde{u}}}{u_{\infty}},\quad\underline{x}=\frac{\underline{\tilde{x}}}{L_{\infty}},\quad
E=\frac{\tilde{E}}{p_{\infty}/\rho_{\infty}},\quad
H=\frac{\tilde{H}}{p_{\infty}/\rho_{\infty}},\quad
p=\frac{\tilde{p}}{p_{\infty}},\quad
t=\frac{\tilde{t}}{L_{\infty}/u_{\infty}}\end{gathered}$ (1)

where tildes indicate local dimensional quantities, $\infty$ indicates
reference dimensional quantities, and underlining indicates vector-valued
quantities. $E$ and $H$ are the specific total energy and enthalpy
respectively. Using this non-dimensionalisation, the Euler equations in the
entropy variables $d\underline{w}=(dp,du,dv,ds)^{T}$ where $ds=dp-a^{2}d\rho$
are:

$\displaystyle\partial_{t}p+\underline{u}\boldsymbol{\cdot}\nabla p+\gamma
p\nabla\boldsymbol{\cdot}\underline{u}=0$ (2)
$\displaystyle\rho\partial_{t}\underline{u}+M^{\scalebox{0.4}[1.0]{$-$}2}\nabla
p+\rho\underline{u}\boldsymbol{\cdot}\nabla\underline{u}=0$
$\displaystyle\partial_{t}s+\underline{u}\boldsymbol{\cdot}\nabla s=0$

where $M=\frac{\sqrt{\gamma}u_{\infty}}{a_{\infty}}$ is the reference Mach
number. Two single scale limits will be shown here for purely convective and
purely acoustic flow. A multiple time (space) scale, single space (time) scale
limit which contains both convective and acoustic flow features also exists,
see [21] and [42, 3] for more detail. To find the convective limit, all
variables in (2) are expanded as power series of the Mach number:

$\psi(\underline{x},t,M)=\psi^{(0)}(\underline{x},t)+M\psi^{(1)}(\underline{x},t)+M^{2}\psi^{(2)}(\underline{x},t)+\mathcal{O}(M^{3})$
(3)

Grouping terms of equal order, we find the following relations for the
pressure and velocity:

$\displaystyle\nabla\,p^{(0,1)}=0$ (4)
$\displaystyle\rho^{(0)}\partial_{t}\underline{u}^{(0)}+\nabla
p^{(2)}+\underline{\rho u}^{(0)}\boldsymbol{\cdot}\nabla\underline{u}^{(0)}=0$
$\displaystyle
d_{t}p^{(0,1)}+\gamma(p\nabla\boldsymbol{\cdot}\underline{u})^{(0,1)}=0$

The pressure has only $\mathcal{O}(M^{2})$ spatial variations, and at steady-
state the divergence of the zeroth and first order velocities vanish, with the
system (4) approaching the incompressible regime as $M\to 0$. The same
relation is found for every order of the entropy equation:

$\partial_{t}s^{(n)}+(\underline{u}\boldsymbol{\cdot}\nabla s)^{(n)}=0$ (5)

The acoustic limit can be found by replacing the convective reference
timescale $L_{\infty}/u_{\infty}$ with the acoustic reference timescale
$L_{\infty}/a_{\infty}$ in the non-dimensionalisation (1). The acoustic non-
dimensional time is then $\tau=a_{\infty}\tilde{t}/L_{\infty}=t/M$. Repeating
the expansion with this timescale we find:

$\displaystyle\nabla\,p^{(0)}=\partial_{\tau}p^{(0)}=0$ (6)
$\displaystyle\rho^{(0)}\partial_{\tau}\underline{u}^{(0)}+\nabla p^{(1)}=0$
$\displaystyle d_{\tau}p^{(1)}+\gamma
p^{(0)}\nabla\boldsymbol{\cdot}\underline{u}^{(0)}=0$

The pressure now has $\mathcal{O}(M^{1})$ spatial variations, the zeroth order
velocity divergence is non-zero, and the system (6) approaches the equations
for low Mach number acoustics. The leading order relation from the entropy
equation is

$\partial_{\tau}s^{(0)}=0$ (7)

which shows that the (leading order) entropy does not vary on the acoustic
timescale. The asymptotic scaling of quantities at each of the two regimes are
shown in table 1, and will be used later to find the limit equations for
different numerical schemes in each regime.

| Convective | Acoustic
---|---|---
$p$ | $p^{(0)}\sim\mathcal{O}(M^{0})$
$u$ | $u^{(0)}\sim\mathcal{O}(M^{0})$
$\rho$ | $\rho^{(0)}\sim\mathcal{O}(M^{0})$
$\partial_{x}p$ | $\partial_{x}p^{(2)}\sim\mathcal{O}(M^{2})$ | $\partial_{x}p^{(1)}\sim\mathcal{O}(M)$
$\partial_{x}u$ | $\partial_{x}u^{(0)}\sim\mathcal{O}(M^{0})$
$\partial_{x}v$ | $\partial_{x}v^{(0)}\sim\mathcal{O}(M^{0})$
$\partial_{x}s$ | $\partial_{x}s^{(0)}\sim\mathcal{O}(M^{0})$
$\partial_{t}p$ | $\partial_{t}p^{(0)}\sim\mathcal{O}(M^{0})$ | $\partial_{\tau}p^{(1)}\sim\mathcal{O}(M^{0})$
$\partial_{t}u$ | $\partial_{t}u^{(0)}\sim\mathcal{O}(M^{0})$ | $\partial_{\tau}u^{(0)}\sim\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}1})$
$\partial_{t}v$ | $\partial_{t}v^{(0)}\sim\mathcal{O}(M^{0})$
$\partial_{t}s$ | $\partial_{t}s^{(0)}\sim\mathcal{O}(M^{0})$
Table 1: Low Mach number scaling of the various terms in equations (4), (6)
for convective or acoustic variations

##  2 Artificial diffusion for low Mach number schemes

In [21], the $x$-split form of the modified equations in the entropy variables
was used to determine the appropriate scaling of the artificial diffusion at
low Mach number. Although convective and acoustic flow features can only
properly be distinguished with a truly multi-dimensional analysis [11, 1],
many finite-volume/difference schemes use a dimension splitting procedure, and
many low Mach number behaviours of these schemes can be demonstrated with the
1D-split equations - with the understanding that the final scheme will be the
sum over all dimensions. The modified equations make the artificial diffusion
terms which are inherent in many numerical schemes explicit in the continuous
equations [67], which allows analysis of the behaviour of these terms without
restriction to a specific discrete form. The modified equations in the
symmetry variables were used in [60] to demonstrate the behaviour of the
artificial diffusion of the standard and preconditioned Roe scheme, and the
modified equations of the barotropic Euler equations were used in [11] to
investigate the requirements for accurate schemes at the convective limit. The
entropy variables are useful at low Mach number because in the 1D-split form
the vorticity (transverse-velocity) and entropy equations are equivalent to
scalar advection equations, which require only the usual upwind diffusion
using the convective velocity scale. This reduces the degrees of freedom we
need to consider to only those in the pressure and (normal-)velocity
equations, which are the distinguishing feature between the different low Mach
number limits and schemes.

In the first part of this section, the form of the artificial diffusion in the
entropy variables used in [21] will be introduced, and the required scaling at
each of the three low Mach number limits will be recalled and briefly
discussed. In the second part of the section, the form of the artificial
diffusion transformed to the conserved variables will be presented, which will
be useful later for comparison with the diffusion of the convection-pressure
flux-vector splitting schemes.

###  2.1 Entropy variables form

The form of the modified equations in the non-dimensional entropy variables
used in [21] to investigate the behaviour of the artificial diffusion
is:§§§The form written here is equivalent to equation (14) with the diffusion
coefficients from equation (41) in [21].

$\displaystyle\partial_{t}p+u\partial_{x}p+\gamma p\partial_{x}u$
$\displaystyle=\phantom{\rho}\mu_{u}|u|\partial_{xx}p+M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p}{\rho|v|}\mu_{11}\partial_{xx}p\pm\gamma p\mu_{12}\partial_{xx}u$ (8)
$\displaystyle\rho\partial_{t}u+M^{\scalebox{0.4}[1.0]{$-$}2}\partial_{x}p+\rho
u\partial_{x}u$
$\displaystyle=\rho\mu_{u}|u|\partial_{xx}u\phantom{\frac{\gamma
p}{\rho|v|}}\pm
M^{\scalebox{0.4}[1.0]{$-$}2}\mu_{21}\partial_{xx}p+\rho|v|\mu_{22}\partial_{xx}u$
$\displaystyle\partial_{t}v+u\partial_{x}v$
$\displaystyle=\phantom{\rho}\mu_{u}|u|\partial_{xx}v$
$\displaystyle\partial_{t}s+u\partial_{x}s$
$\displaystyle=\phantom{\rho}\mu_{u}|u|\partial_{xx}s$

The terms with the $\mu_{u}$ diffusion coefficient are the convective
upwinding, so $\mu_{u}\sim\mathcal{O}(M^{0})$, and $|v|$ is some velocity
scale. The diffusion terms with the $\mu_{ij}$, $i,j=1,2$ coefficients are
elements of the artificial diffusion Jacobian $\underline{\underline{A}}$. The
form of the Jacobian elements used here was chosen because it simplifies the
diffusion coefficients once we transform to the conserved variables below. The
precise form of $\mu_{ij}$, $i,j=1,2$ is particular to each discrete scheme.
The sign of the off-diagonal terms may be negative so long as the diffusion
Jacobian remains positive (semi-)definite. See [21] for a discussion of the
requirements on the form of the off-diagonal diffusion components.

In [21] the method of Turkel 1999 [60] was used to find the limiting forms of
(8) for each of the convective and acoustic limits by enforcing the variations
in table 1. The Mach number scaling of the coefficients $\mu_{ij}$ was then
chosen for each of the two limits such that the artificial diffusion terms had
the same Mach number scaling as the dominant physical terms on the left hand
side. This ensures that the artificial diffusion terms do not vanish
asymptotically - which would result in loss of stability - nor do they
dominate asymptotically - which would result in loss of accuracy. The
$M^{\scalebox{0.4}[1.0]{$-$}2}$ coefficients in $A_{11}$ and $A_{21}$ are
included so that at the convective limit all $\mu_{ij}$ are independent of the
Mach number, as is appropriate for the incompressible regime, while keeping
the correct scaling for the full Jacobian elements $A_{ij}$. The convective
and acoustic diffusion scalings were then combined into a mixed diffusion
scaling which is suitable for both convective and acoustic flow features. The
scaling of the coefficients for the convective
($\underline{\underline{A}}^{c}$), acoustic ($\underline{\underline{A}}^{a}$)
and mixed ($\underline{\underline{A}}^{m})$ diffusion Jacobians are shown in
table 2. Each of these three diffusion scalings has its own strengths and
weaknesses. A brief description is given here, with more detail given in [21]
and references therein.

| $\underline{\underline{A}}^{c}$ | $\underline{\underline{A}}^{a}$ | $\underline{\underline{A}}^{m}$
---|---|---|---
$\begin{matrix}\begin{pmatrix}A_{11}&A_{12}\\\ A_{21}&A_{22}\end{pmatrix}\end{matrix}$ | $\mathcal{O}\begin{matrix}\begin{pmatrix}M^{\scalebox{0.4}[1.0]{$-$}2}&M^{0}\\\ M^{\scalebox{0.4}[1.0]{$-$}2}&M^{0}\end{pmatrix}\end{matrix}$ | $\mathcal{O}\begin{matrix}\begin{pmatrix}M^{\scalebox{0.4}[1.0]{$-$}1}&M^{0}\\\ M^{\scalebox{0.4}[1.0]{$-$}2}&M^{\scalebox{0.4}[1.0]{$-$}1}\end{pmatrix}\end{matrix}$ | $\mathcal{O}\begin{matrix}\begin{pmatrix}M^{\scalebox{0.4}[1.0]{$-$}1}&M^{0}\\\ M^{\scalebox{0.4}[1.0]{$-$}2}&M^{0}\end{pmatrix}\end{matrix}$
$\begin{matrix}\begin{pmatrix}\mu_{11}&\mu_{12}\\\ \mu_{21}&\mu_{22}\end{pmatrix}\end{matrix}$ | $\mathcal{O}\begin{matrix}\begin{pmatrix}M^{0}&M^{0}\\\ M^{0}&M^{0}\end{pmatrix}\end{matrix}$ | $\mathcal{O}\begin{matrix}\begin{pmatrix}M&M^{0}\\\ M^{0}&M^{\scalebox{0.4}[1.0]{$-$}1}\end{pmatrix}\end{matrix}$ | $\mathcal{O}\begin{matrix}\begin{pmatrix}M&M^{0}\\\ M^{0}&M^{0}\end{pmatrix}\end{matrix}$
Suitable for convective flow | Yes | No | Yes1
Suitable for acoustic flow | No | Yes | Yes2
Timestep limit | $\mathcal{O}(M^{2})$ | $\mathcal{O}(M^{1})$ | $\mathcal{O}(M^{1})$
Table 2: Details of the convective $\underline{\underline{A}}^{c}$, acoustic
$\underline{\underline{A}}^{a}$, and mixed $\underline{\underline{A}}^{m}$
diffusion scalings. The mixed scheme is susceptible to 1pressure chequerboard
instabilities at the convective limit and 2velocity instabilities at the
acoustic limit.

The convective diffusion scaling has a large
$\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}2})$ pressure diffusion coefficient
$A_{11}$ to compensate for the $\mathcal{O}(M^{2})$ scaling of $\partial_{x}p$
at the convective limit. This term prevents pressure chequer-board
instabilities, which result from the same pressure-velocity decoupling that
occurs with collocated schemes for incompressible flow. This increased
pressure diffusion also damps out any acoustic features in the flow, making
the convective scheme unsuitable for the acoustic or mixed limits. The
pressure diffusion also increases the spectral radius of the diffusion matrix
by a factor of $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}1})$ compared to the
physical flux Jacobian and the other two diffusion matrix scalings. This means
that an implicit time integration strategy is required for time-accurate
simulations, as the timestep restriction scale as $\Delta
t\sim\mathcal{O}(M^{2})$ instead of the usual $\Delta t\sim\mathcal{O}(M)$ for
stability of the explicit scheme.

The acoustic diffusion scaling has a smaller pressure diffusion coefficient,
so can resolve acoustic flow features. However, the velocity diffusion
$A_{22}$ is a factor of $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}1})$ larger
than for the convective scaling, which destroys the accuracy for convective
flow features and makes this scaling unsuitable for the convective or mixed
limits. Most schemes for compressible flow designed for higher Mach number
regimes approach this scaling at low Mach number, and hence are inaccurate at
low Mach number unless at least this term is modified.

The mixed scaling has neither the $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}2})$
pressure diffusion of the convective diffusion nor the
$\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}1})$ velocity diffusion of the
acoustic scaling. As such it does not overdamp acoustic nor convective flow
features and is accurate for the convective, acoustic, and mixed limits,
although this flexibility does not come without cost. The pressure diffusion
vanishes asymptotically at the purely convective limit, meaning that this
scheme is susceptible to pressure chequerboard instabilities. The velocity
diffusion vanishes asymptotically at the purely acoustic limit, meaning that
there is an undamped grid mode which can be excited by non-linearities,
initial or boundary conditions, poor mesh quality, and acoustic
discontinuities such as shockwaves.

###  2.2 Conserved variables form

Transforming the modified equations in the non-dimensional entropy variables
(8) to the dimensional conserved variables, we can construct a numerical
interface flux $\underline{\hat{f}}$ over a face with normal $\underline{n}$:

$\underline{\hat{f}}=\\{\underline{\tilde{f}}\\}-\frac{1}{2}\Bigg{[}\mu_{u}|\tilde{U}|\begin{matrix}\begin{pmatrix}\Delta\tilde{\rho}\\\
\Delta\tilde{\rho}\tilde{\underline{u}}\\\
\Delta\tilde{\rho}\tilde{E}\end{pmatrix}\end{matrix}+\delta\tilde{U}\begin{matrix}\begin{pmatrix}\tilde{\rho}\\\
\tilde{\rho}\tilde{\underline{u}}\\\
\tilde{\rho}\tilde{H}\end{pmatrix}\end{matrix}+\delta\tilde{p}\begin{matrix}\begin{pmatrix}0\\\
\underline{n}\\\
\tilde{U}\end{pmatrix}\end{matrix}\Bigg{]}=\\{\underline{\tilde{f}}\\}-\underline{\tilde{f}}^{d}$
(9)
$\begin{split}\delta\tilde{U}&=\frac{\mu_{11}}{\tilde{\rho}|\tilde{v}|}\Delta\tilde{p}+\frac{\tilde{U}}{|\tilde{U}|}\mu_{12}\Delta\tilde{U}\\\
\delta\tilde{p}&=\frac{\tilde{U}}{|\tilde{U}|}\mu_{21}\Delta\tilde{p}+\tilde{\rho}|\tilde{v}|\mu_{22}\Delta\tilde{U}\end{split}$
(10)

Where
$\\{(\boldsymbol{\cdot})\\}=\frac{1}{2}((\boldsymbol{\cdot})_{L}+(\boldsymbol{\cdot})_{R})$
and
$\Delta(\boldsymbol{\cdot})=(\boldsymbol{\cdot})_{R}-(\boldsymbol{\cdot})_{L}$
are the central average and the jump between the left/right states of the
interface and $U=\underline{u}\boldsymbol{\cdot}\underline{n}$ is the
interface-normal velocity. The first term $\\{\underline{\tilde{f}}\\}$ is the
central approximation to the physical flux, and the other terms are the
artificial diffusion $\underline{\tilde{f}}^{d}$. The first of these is the
convective upwinding on the conserved variables. The interface velocity
perturbation $\delta\tilde{U}$ comes from the diffusion on the pressure in the
entropy variables, and the interface pressure perturbation $\delta\tilde{p}$
comes from the diffusion on the velocity in the entropy variables. The switch
$\tilde{U}/|\tilde{U}|$ ensures that the off-diagonal terms maintain a
diffusive character. Matrix Rusanov fluxes such as the Roe scheme [47] can be
arranged in the form of (9) [38, 68] to compare existing Roe-type schemes at
low Mach number to the scalings in table 2, as was done in [30] and [21]. This
form will be referred to as the Roe form for the rest of the paper.

##  3 Convection-pressure flux-vector splittings at low Mach number

The form (9,10) is useful for analysing interface fluxes as it explicitly
separates the non-diffusive central component from the stabilising artificial
diffusion components, however the artificial diffusion in convection-pressure
flux-vector splitting schemes is usually implicit in the formulation. In this
section, we construct approximate forms for three families of convection-
pressure flux-vector splittings in the conserved variables that explicitly
separate the artificial diffusion terms, so that we can make comparison to
(9,10). We then transform the splittings to the entropy variables to make
direct comparison to (8) and the results of [21].

###  3.1 Conserved variables form of convection-pressure flux-vector
splittings

Convection-pressure flux-vector splittings split the physical flux vector
$\underline{f}$ into two components:

$\underline{f}=U\underline{\varphi}+\underline{P}$ (11)

where $U$ is the local convecting velocity, $\underline{\varphi}$ is some
vector of convected quantities, and $\underline{P}$ is some vector
representing the pressure flux. The choice of $\underline{\varphi}$ and
$\underline{P}$ specify the flux-splitting. For the three flux-splittings
considered in this paper, the pressure flux can be written as
$\underline{P}=p\underline{\Pi}$, where $p$ is the local pressure. A numerical
interface flux can be constructed from the flux-splitting (11) with the form:

$\underline{\hat{f}}=U_{1/2}\underline{\varphi}+\underline{P}_{1/2}$ (12)

where $U_{1/2}$ is an interface velocity and $\underline{P}_{1/2}$ is an
interface pressure flux. We assume that the interface velocity and pressure
can be rearranged (possibly approximately) into the form:

$U_{1/2}=\\{U\\}-\delta
U,\quad\underline{P}_{1/2}=\\{\underline{P}\\}-\delta\underline{P}$ (13)

where the interface velocity and pressure perturbations $\delta U$ and
$\delta\underline{P}$ are diffusion terms similar to those in (10). The
interface flux (12) can then be rearranged into a form that separates the
central and diffusive components by assuming $U_{1/2}>0$ and
$\underline{\varphi}$ is taken as the upwind value based on $U_{1/2}$, as is
usually the case. Adding
$(U_{R}\underline{\varphi}_{R}-U_{R}\underline{\varphi}_{R})$ to the right
hand side, we obtain:

$\underline{\hat{f}}=\\{\underline{f}\\}-\Big{(}\frac{U_{R}}{2}\Delta\underline{\varphi}+\delta
U\underline{\varphi}_{L}+\delta\underline{P}\Big{)}$ (14)

The terms in brackets are all artificial diffusion terms. Assuming smooth
flow, the differences between left and right values are small, and we can
approximate the flux as:

$\underline{\hat{f}}=\\{\underline{f}\\}-\Big{(}\frac{U}{2}\Delta\underline{\varphi}+\delta
U\underline{\varphi}+\delta\underline{P}\Big{)}$ (15)

where $U$ and $\underline{\varphi}$ are some consistent approximation of the
local values. Comparing to equations (9,10) in section 2, it can be seen that
equation (15) closely approximates the Roe form of the artificial diffusion.
The first diffusion term is the convective upwinding, the second term is
diffusion due to the interface velocity perturbation, and the third term is
due to the interface pressure flux perturbation. The interface velocity
perturbation $\delta U$ is assumed to have the same form as in the Roe
diffusion (10), however the form of interface pressure flux perturbation
$\delta\underline{P}$ can be different for each splitting. The recasting of
the convection term into central and upwinding components was shown for AUSM
in the original paper from Liou & Steffen [37], and the splitting of the
pressure term was also used by Liou [33] and Edwards & Liou [14, 35, 12]. The
CUSP scheme [23] uses explicit diffusion terms with a similar form. Now we
have an explicit expression for the artificial diffusion we can substitute the
specific forms of $\underline{\varphi}$ and $\underline{P}$ for each flux
splitting.

####  3.1.1 Liou-Steffen

For the Liou-Steffen splitting we have:

$\underline{\varphi}=\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\
\rho
H\end{pmatrix}\end{matrix},\quad\underline{P}=p\underline{\Pi}=p\begin{matrix}\begin{pmatrix}0\\\
n_{x}\\\ n_{y}\\\ 0\end{pmatrix}\end{matrix},\quad\delta\underline{P}=\delta
p\underline{\Pi}$ (16)

where $\delta p$ is the same interface pressure perturbation as the Roe
diffusion (10). Substituting (16) into (15), the approximate form of the
diffusion is:

$\underline{\hat{f}}=\\{\underline{f}\\}-\frac{U}{2}\begin{matrix}\begin{pmatrix}\Delta\rho\\\
\Delta\rho u\\\ \Delta\rho v\\\ \Delta\rho H\end{pmatrix}\end{matrix}-\delta
U\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\ \rho
H\end{pmatrix}\end{matrix}-\delta p\begin{matrix}\begin{pmatrix}0\\\ n_{x}\\\
n_{y}\\\ 0\end{pmatrix}\end{matrix}$ (17)

The first term, the convective upwinding, is almost exactly the same as the
Roe diffusion, except that the energy equation is upwinded on jumps in the
total enthalpy, not total energy. For low Mach number flows, the zeroth order
pressure is a factor of $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}2})$ larger
than the kinetic energy, so to leading order the total energy and total
enthalpy differ by a factor of $\gamma$. The second term, the velocity
perturbation diffusion, is exactly that of the Roe diffusion, which is
particularly important for preventing pressure instabilities in purely
convective flows. The third term, the pressure perturbation diffusion,
naturally has the form $\delta p\underline{\Pi}$ and matches the Roe form for
the momentum equations but provides no diffusion on the energy equation,
unlike the Roe form.

####  3.1.2 Zha-Bilgen

For the Zha-Bilgen splitting we have:

$\underline{\varphi}=\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\
\rho
E\end{pmatrix}\end{matrix},\quad\underline{P}=p\underline{\Pi}=p\begin{matrix}\begin{pmatrix}0\\\
n_{x}\\\ n_{y}\\\
U\end{pmatrix}\end{matrix},\quad\delta\underline{P}=\begin{matrix}\begin{pmatrix}0\\\
\delta p\,n_{x}\\\ \delta p\,n_{y}\\\ \delta(pU)\end{pmatrix}\end{matrix},$
(18)

leading to the approximate form:

$\underline{\hat{f}}=\\{\underline{f}\\}-\frac{U}{2}\begin{matrix}\begin{pmatrix}\Delta\rho\\\
\Delta\rho u\\\ \Delta\rho v\\\ \Delta\rho E\end{pmatrix}\end{matrix}-\delta
U\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\ \rho
E\end{pmatrix}\end{matrix}-\begin{matrix}\begin{pmatrix}0\\\ \delta
p\,n_{x}\\\ \delta p\,n_{y}\\\ \delta(pU)\end{pmatrix}\end{matrix}$ (19)

The convective upwinding term is exactly that of the Roe form, so can be
expected to give good resolution for the convective features. The velocity
perturbation term matches the Roe term in the mass and momentum equations, but
uses the total energy scale in the energy equation instead of the total
enthalpy. The third term matches the Roe form except possibly in the energy
equation, which depends on the specific expression for $\delta(pU)$. If we use
the expression:

$\delta(pU)=U\delta p$ (20)

then we can find a first form for the diffusion in terms of $\delta U$ and
$\delta p$:

$\underline{\hat{f}}=\\{\underline{f}\\}-\frac{U}{2}\begin{matrix}\begin{pmatrix}\Delta\rho\\\
\Delta\rho u\\\ \Delta\rho v\\\ \Delta\rho E\end{pmatrix}\end{matrix}-\delta
U\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\ \rho
E\end{pmatrix}\end{matrix}-\delta p\begin{matrix}\begin{pmatrix}0\\\ n_{x}\\\
n_{y}\\\ U\end{pmatrix}\end{matrix}$ (21)

where the pressure perturbation term now has the form $\delta
p\underline{\Pi}$, and is exactly that of the Roe scheme. A second form is
found if the $\delta(pU)$ term is instead written as:

$\delta(pU)=p\delta U+U\delta p$ (22)

Then the total diffusion on the energy equation is:

$\begin{split}\underline{\hat{f}}^{d}_{\rho E}&=U\Delta\rho E+(\rho E+p)\delta
U+U\delta p\\\ &=U\Delta\rho E+\rho H\delta U+U\delta p\end{split}$ (23)

which matches that in the Roe-type scheme (9).

####  3.1.3 Toro-Vasquez

For the Toro-Vasquez splitting we have:

$\underline{\varphi}=\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\
\rho
k\end{pmatrix}\end{matrix},\quad\underline{P}=p\underline{\Pi}=p\begin{matrix}\begin{pmatrix}0\\\
n_{x}\\\ n_{y}\\\ \dfrac{\gamma
U}{\gamma-1}\end{pmatrix}\end{matrix},\quad\delta\underline{P}=\begin{matrix}\begin{pmatrix}0\\\
\delta p\,n_{x}\\\ \delta p\,n_{y}\\\
\dfrac{\gamma\delta(pU)}{\gamma-1}\end{pmatrix}\end{matrix},$ (24)

leading to the approximate form:

$\underline{\hat{f}}=\\{\underline{f}\\}-\frac{U}{2}\begin{matrix}\begin{pmatrix}\Delta\rho\\\
\Delta\rho u\\\ \Delta\rho v\\\ \Delta\rho k\end{pmatrix}\end{matrix}-\delta
U\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\ \rho
k\end{pmatrix}\end{matrix}-\begin{matrix}\begin{pmatrix}0\\\ \delta
p\,n_{x}\\\ \delta p\,n_{y}\\\
\dfrac{\gamma\delta(pU)}{\gamma-1}\end{pmatrix}\end{matrix}$ (25)

The Toro-Vasquez diffusion has the least similarity to the Roe form. The
diffusion on the continuity and momentum equations match, but the energy
equation diffusion differs in all three terms. As stated above, the ratio of
the kinetic energy to the internal energy is $\mathcal{O}(M^{2})$, so the
convective upwind diffusion and velocity perturbation diffusion on the energy
equation in this form of the Toro-Vasquez splitting will be significantly
smaller than in the Liou-Steffen and Zha-Bilgen splittings. If we again use
the first expression for $\delta(pU)$, equation (20), then the first form of
the diffusion in terms of $\delta U$ and $\delta p$ is:

$\underline{\hat{f}}=\\{\underline{f}\\}-\frac{U}{2}\begin{matrix}\begin{pmatrix}\Delta\rho\\\
\Delta\rho u\\\ \Delta\rho v\\\ \Delta\rho k\end{pmatrix}\end{matrix}-\delta
U\begin{matrix}\begin{pmatrix}\rho\\\ \rho u\\\ \rho v\\\ \rho
k\end{pmatrix}\end{matrix}-\delta p\begin{matrix}\begin{pmatrix}0\\\ n_{x}\\\
n_{y}\\\ \dfrac{\gamma U}{\gamma-1}\end{pmatrix}\end{matrix}$ (26)

The third diffusion term now has the form $\delta p\underline{\Pi}$, and
matches the Roe form, up to a factor $\frac{\gamma}{\gamma-1}$ in the energy
equation. A second form of the Toro-Vasquez diffusion can be found similarly
to the second Zha-Bilgen form by using (22). This results in the total
diffusion on the energy equation having the form:

$\begin{split}\underline{\hat{f}}^{d}_{\rho E}&=U\Delta\rho k+\Big{(}\rho
k+\frac{\gamma p}{\gamma-1}\Big{)}\delta U+\frac{\gamma U\delta
p}{\gamma-1}\\\ &=U\Delta\rho k+\rho H\delta U+\frac{\gamma U\delta
p}{\gamma-1}\end{split}$ (27)

where the velocity perturbation now matches the Roe form, but the convective
upwinding and the pressure perturbation terms remain unchanged.

In this section it has been shown how convection-pressure flux-vector
splittings can be approximately reformulated to make the artificial diffusion
terms explicit, allowing comparison with the findings of the first paper. This
reformulation requires a number of assumptions which limit the scope of its
applicability. The first two assumptions are that the convection term in (12)
can be cast as $U_{1/2}\underline{\varphi}$, and that the interface velocity
and pressure flux can be cast as the sum of a central approximation and a
diffusive term (13), which hold true for many schemes in the literature. In
this study we are interested in the asymptotic scaling of the artificial
diffusion coefficients $\mu_{ij}$, and not their specific form. This means
that, for the purpose of our analysis, the reformulation into (12,13) can be
approximate and these first two assumptions are not particularly limiting even
when they do not hold exactly. The next assumptions are on the form of the
pressure diffusion term $\delta\underline{P}$. Firstly, that this term has the
same sparsity as the pressure flux $\underline{P}$. While this may seem like a
logical assumption, we will see later that it does not apply to all schemes in
the literature. The second assumption is the simplification
$\delta(pU)=U\delta p$ (20), which allows the pressure diffusion term of the
Zha-Bilgen and Toro-Vasquez splittings to be written as
$\delta\underline{P}=\delta p\underline{\Pi}$, as it naturally is for the
Liou-Steffen splitting and the Roe flux. However, (20) will not necessarily
hold for all schemes even if the sparsity assumption holds, in which case the
second forms of the Zha-Bilgen and Toro-Vasquez splittings may be more
applicable. Lastly, we note that the only assumption which is specific to low
Mach number flow is that the flow is smooth, which we used in the step from
(14) to (15).

The assumptions on the diffusion terms enforce only light requirements, and
result in diffusion forms which closely match that of the Roe scheme so allow
better comparison to the findings of [21]. The approximate diffusion terms of
all of these forms match the Roe form on the continuity and momentum
equations, although they all differ to some extent on the energy equation. The
Liou-Steffen and both of the Zha-Bilgen forms match the Roe form of the
artificial diffusion exactly in at least one of the terms (all three in the
case of the second Zha-Bilgen form (23)), and approximately match in the
others. On the other hand, the Toro-Vasquez splitting shows little resemblance
to the Roe form in the energy equation, other than in the second form (27)
where the velocity perturbation term matches.

###  3.2 Entropy variables form

The artificial diffusion fluxes (17), (21), and (26) are now transformed to
the non-dimensional entropy variables, and the resulting modified equations
are compared to (8) and the findings of the previous paper. The Roe type
diffusion (9,10) in the entropy variables is:

$\mu_{u}|u|\begin{matrix}\begin{pmatrix}\phantom{\rho}\partial_{xx}p\\\
\rho\partial_{xx}u\\\ \rho\partial_{xx}v\\\
\phantom{\rho}\partial_{xx}s\end{pmatrix}\end{matrix}+\delta
U\begin{matrix}\begin{pmatrix}\gamma p\\\ 0\\\ 0\\\
0\end{pmatrix}\end{matrix}+\delta p\begin{matrix}\begin{pmatrix}0\\\ 1\\\ 0\\\
0\end{pmatrix}\end{matrix}$ (28)

which is exactly the right hand side of (8) written in terms of $\delta U$ and
$\delta p$. The diffusion terms found from the Liou-Steffen splitting (17)
are:

$\mu_{u}|u|\begin{matrix}\begin{pmatrix}\gamma\partial_{xx}p\\\
\rho\partial_{xx}u\\\ \rho\partial_{xx}v\\\
\partial_{xx}s+(\gamma-1)\partial_{xx}p\end{pmatrix}\end{matrix}+\delta
U\begin{matrix}\begin{pmatrix}\gamma p\\\ 0\\\ 0\\\
0\end{pmatrix}\end{matrix}+\delta
p\begin{matrix}\begin{pmatrix}-M^{2}u(\gamma-1)\\\ 1\\\ 0\\\
-M^{2}u(\gamma-1)\end{pmatrix}\end{matrix}$ (29)

As expected, the second term in (29) exactly matches that of the Roe scheme.
The convective upwinding is slightly different to the Roe form - increased by
a factor of $\gamma$ on the pressure, and the upwinding on the entropy field
also has a contribution from the pressure variations. Interestingly, the
upwinding terms on the entropy field are equivalent to temperature diffusion:
$ds+(\gamma-1)dp=\rho R\gamma dT$. The pressure perturbation diffusion in the
third term, which should act only on the velocity field, now also has anti-
diffusion terms on the pressure and entropy fields, although these terms are
reduced by a factor of $M^{2}$.

For the first form of the Zha-Bilgen splitting (21) we have:

$\mu_{u}|u|\begin{matrix}\begin{pmatrix}\phantom{\rho}\partial_{xx}p\\\
\rho\partial_{xx}u\\\ \rho\partial_{xx}v\\\
\phantom{\rho}\partial_{xx}s\end{pmatrix}\end{matrix}+\delta
U\begin{matrix}\begin{pmatrix}p\\\ 0\\\ 0\\\
-(\gamma-1)p\end{pmatrix}\end{matrix}+\delta
p\begin{matrix}\begin{pmatrix}0\\\ 1\\\ 0\\\ 0\end{pmatrix}\end{matrix}$ (30)

As expected, the first and third terms in (30) exactly match those in the Roe
scheme. However, the $\delta U$ diffusion is a factor of $\gamma$ smaller on
the pressure field, and has introduced an anti-diffusion term on the entropy
field. By inspection, the second Zha-Bilgen form (23) will appear the same as
the Roe form (28), with the correct pressure diffusion $\delta U$ term.
Finally, for the first form of the Toro-Vasquez splitting (26) we have:

$\mu_{u}|u|\begin{matrix}\begin{pmatrix}0\\\ \rho\partial_{xx}u\\\
\rho\partial_{xx}v\\\
\partial_{xx}s-\partial_{xx}p\end{pmatrix}\end{matrix}+\delta
U\begin{matrix}\begin{pmatrix}0\\\ 0\\\ 0\\\ -\gamma
p\end{pmatrix}\end{matrix}+\delta p\begin{matrix}\begin{pmatrix}M^{2}u\\\ 1\\\
0\\\ M^{2}u\end{pmatrix}\end{matrix}$ (31)

None of the three terms in (31) match those in (28). The upwinding and $\delta
U$ term provide no diffusion on the pressure field, which only has (erroneous)
diffusion from the $\delta p$ term scaled by $M^{2}$. The upwinding on the
entropy field has a contribution from the pressure variations, similar to the
Liou-Steffen form. This time, the upwinding on the entropy field is equivalent
to density diffusion: $ds-dp=-a^{2}d\rho$.¶¶¶The negative sign does not
indicate anti-diffusion because $\frac{\partial s}{\partial\rho}<0$. The only
effect of the $\delta U$ term is an anti-diffusion term on the entropy
equation, similar to that in first form of the Zha-Bilgen splitting. By
inspection, the second Toro-Vasquez form will appear as (31) but with the
velocity perturbation form of the Roe scheme, which removes the pressure anti-
diffusion term from the entropy equation and restores it to the pressure
equation.

Because all forms are consistent with the Roe form on the mass and momentum
equations in the conserved variables, they all have the correct diffusion on
the velocity and vorticity equations in the entropy variables. However, the
diffusion on the pressure and entropy fields does not fully match the Roe form
for any of the flux-splitting forms, apart from the second Zha-Bilgen form.

###  3.3 Limit equations of flux-vector splittings

In this section the limiting forms of the modified equations of each diffusion
scheme are found for the convective and acoustic asymptotic regimes. The limit
equations are found with the following steps:

  1. 1.

Replace the right hand side of the modified equations (8) with one of the
diffusion schemes (28-31).

  2. 2.

Set the scaling of $\mu_{11}$, $\mu_{22}$ according to either the convective
or mixed diffusion scaling in table 2.

  3. 3.

Enforce the scaling of the physical quantities and their gradients according
to either the convective or acoustic regime scaling in table 1.

  4. 4.

Retain only the leading order terms in each equation.

The limit equations provide less detail than full asymptotic expansions of the
modified equations would, but are significantly more compact and still give
useful insight into the asymptotic behaviour of schemes at low Mach number, as
shown in [21]. The velocity equations for all three flux splittings match the
Roe form which means that with the acoustic diffusion scaling they will all be
inaccurate for convective features. For this reason only the limit equations
of convective and mixed diffusion scalings will considered. To reduce the size
of the relations, all schemes are taken to have asymptotically diagonal
diffusion Jacobians i.e. $\mu_{12}\sim\mu_{21}\sim o(M^{0})$, which is almost
always the case for convection-pressure flux-vector splittings at low Mach
number.

####  3.3.1 Roe form

The limit equations for the ‘ideal’ modified equations (8) are included here
for the flux-vector-splitting schemes to be compared against. These limit
equations are equivalent to those presented in section 3 of [21], except with
the specific form of the diffusion coefficients from the present equation (8).
Because (8) are the modified equations for Roe-type schemes in the entropy
variables, we shall refer to them as the Roe form of the modified equations.
Their main features are described here; see [21] for more in-depth discussion.
The limit equations for the Roe form (28) using the convective diffusion
scaling and enforcing convective variations from table 1 are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(2)}$ (32)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}v^{(0)}+u^{(0)}\partial_{x}v^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}v^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}$

The limit equations for the Roe form (28) using the convective diffusion
scaling and enforcing acoustic variations from table 1 are:

$\displaystyle 0$ $\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$ (33)
$\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}v^{(0)}$ $\displaystyle=0$
$\displaystyle\partial_{\tau}s^{(0)}$ $\displaystyle=0$

The left-hand-side of the convective limit equations match those of the single
scale convective asymptotic expansion (4) by construction. The diffusion on
the velocity equation is well balanced, and the vorticity and entropy
equations appear as scalar advection-diffusion equations with upwind
diffusion. Note that at steady state the divergence does not disappear, but is
equal to the pressure diffusion term. As mentioned previously, this term damps
pressure oscillations from the solution, preventing chequerboard modes. In the
acoustic limit, this term overwhelms the physical terms in the pressure
equation, causing it to become a Poisson-like relation for the acoustic
pressure $p^{(1)}$ which quickly damps out any acoustic waves from the
solution.

The limit equations for the Roe form (28) using the mixed diffusion scaling
and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (34)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}v^{(0)}+u^{(0)}\partial_{x}v^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}v^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}$

The limit equations for the Roe form (28) using the mixed diffusion scaling
and enforcing acoustic variations are:

$\displaystyle\partial_{\tau}p^{(1)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$ (35)
$\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}v^{(0)}$ $\displaystyle=0$
$\displaystyle\partial_{\tau}s^{(0)}$ $\displaystyle=0$

Because the only difference between the convective and mixed diffusion
scalings is in the pressure diffusion term, the limit equations for velocity,
vorticity, and entropy with the mixed diffusion scaling are identical to those
with the convective diffusion scaling, with only the pressure equation
differing. At the convective limit the pressure diffusion disappears from the
continuity equation. This allows for divergence-free solutions, consistent
with the incompressible limit but, as stated previously, makes the scheme
susceptible to pressure chequerboard instabilities. At the acoustic limit, the
pressure diffusion no longer overwhelms the physical terms. The left-hand-side
of the pressure and velocity equations are the equations for low Mach number
acoustics. The diffusion on the velocity vanishes asymptotically, making the
mixed scaling susceptible to velocity instabilities at the acoustic limit. The
leading order entropy does not vary on the acoustic timescale, consistent with
the approximation of isentropic acoustics.

Because the diffusion on the velocity and vorticity equations match the Roe
scheme for all three flux-vector splittings, these limit equations will match
those of the Roe form. In the following sections, the limit equations for the
vorticity will be omitted, but the limit equations for the velocity will be
retained so that the pressure-velocity subsystem can be understood as a whole.

####  3.3.2 Liou-Steffen splitting

The limit equations for the Liou-Steffen form (29) using the convective
diffusion scaling and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(2)}$ (36)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}$

The limit equations for the Liou-Steffen form (29) using the convective
diffusion scaling and enforcing acoustic variations are:

$\displaystyle 0$ $\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$ (37)
$\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}s^{(0)}$ $\displaystyle=0$

The limit equations for the Liou-Steffen form (29) using the mixed diffusion
scaling and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (38)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}$

The limit equations for the Liou-Steffen form (29) using the mixed diffusion
scaling and enforcing acoustic variations are:

$\displaystyle\partial_{\tau}p^{(1)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$ (39)
$\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}s^{(0)}$ $\displaystyle=0$

The erroneous diffusion from the $\delta p$ terms vanishes asymptotically from
the pressure and entropy equations for both the convective and mixed diffusion
scalings at both the convective and acoustic limits. The pressure diffusion in
the entropy equation from the convective upwinding term also vanishes at both
limits. As a result, the limit equations for the Liou-Steffen form exactly
match those of the Roe form. This should not be surprising given the success
of low Mach number schemes based on this splitting with both convective and
mixed diffusion scalings.

####  3.3.3 Zha-Bilgen splitting

The limit equations for the first Zha-Bilgen form (30) using the convective
diffusion scaling and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(2)}$
(40) $\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}-M^{\scalebox{0.4}[1.0]{$-$}2}(\gamma-1)\frac{p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(2)}$

The limit equations for the Zha-Bilgen form (30) using the convective
diffusion scaling and enforcing acoustic variations are:

$\displaystyle 0$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$
(41) $\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}s^{(0)}$
$\displaystyle=-M^{\scalebox{0.4}[1.0]{$-$}2}(\gamma-1)\frac{p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$

The pressure diffusion is a factor of $\gamma$ smaller than for the Roe and
Liou-Steffen forms, although if $\gamma\sim\mathcal{O}(1)$ then this should
not have a significant effect, and can easily be compensated for in the
definition of $\mu_{11}$. Other than this difference, the pressure, velocity,
and vorticity equations match the Roe and Liou-Steffen forms, so can be
expected to perform similarly. On the other hand, the entropy equation retains
the pressure anti-diffusion term in both limits, which could cause erroneous
entropy generation from leading order pressure variations.
The limit equations for the Zha-Bilgen form (30) using the mixed diffusion
scaling and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (42)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}$

The limit equations for the Zha-Bilgen form (30) using the mixed diffusion
scaling and enforcing acoustic variations are:

$\displaystyle\partial_{\tau}p^{(1)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=M^{\scalebox{0.4}[1.0]{$-$}2}\frac{p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$
(43) $\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}s^{(0)}$ $\displaystyle=0$

With the mixed diffusion scaling the pressure anti-diffusion term in the
entropy equation vanishes asymptotically at both limits, so the Zha-Bilgen
form matches the Roe form at the convective limit, and the only difference at
the acoustic limit is the factor of $\gamma$ in the pressure equation
diffusion.

The limit equations for the second Zha-Bilgen form (23) will be identical to
the limit equations for the Roe form, which will remove the erroneous pressure
anti-diffusion term from the entropy equation at the convective limit.

####  3.3.4 Toro-Vasquez splitting

The limit equations for the first Toro-Vasquez form (31) using the convective
diffusion scaling and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (44)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}-M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(2)}$

The limit equations for the Toro-Vasquez form (31) using the convective
diffusion scaling and enforcing acoustic variations are:

$\displaystyle\partial_{\tau}p^{(1)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (45)
$\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}s^{(0)}$
$\displaystyle=-M^{\scalebox{0.4}[1.0]{$-$}2}\frac{\gamma
p^{(0)}}{\rho^{(0)}|v|}\mu_{11}\partial_{xx}p^{(1)}$

The lack of pressure diffusion on the pressure equation means that, even with
the convective diffusion scaling, the divergence relation is undamped and
therefore susceptible to pressure chequerboard instabilities. At the acoustic
limit, the absence of this pressure diffusion means that the pressure equation
does not become a Poisson-like relation for the acoustic pressure $p^{(1)}$
which would damp out any acoustic waves. The erroneous $\delta p$ terms in the
pressure and entropy equations vanish asymptotically at both limits, as does
the pressure diffusion in the entropy equation from the convective upwinding
term, leaving the correct entropy upwinding. However, the pressure anti-
diffusion term from $\delta U$ remains on the entropy equation, as for the
Zha-Bilgen form.

The limit equations for the Toro-Vasquez form (31) using the mixed diffusion
scaling and enforcing convective variations are:

$\displaystyle\partial_{t}p^{(0)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (46)
$\displaystyle\rho^{(0)}\partial_{t}u^{(0)}+\partial_{x}p^{(2)}+\rho
u^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=\rho^{(0)}\mu_{u}|u|\partial_{xx}u^{(0)}+\rho^{(0)}|v|\mu_{22}\partial_{xx}u^{(0)}$
$\displaystyle\partial_{t}s^{(0)}+u^{(0)}\partial_{x}s^{(0)}$
$\displaystyle=\mu_{u}|u|\partial_{xx}s^{(0)}$

The limit equations for the Toro-Vasquez form (31) using the mixed diffusion
scaling and enforcing acoustic variations are:

$\displaystyle\partial_{\tau}p^{(1)}+\gamma p^{(0)}\partial_{x}u^{(0)}$
$\displaystyle=0$ (47)
$\displaystyle\rho^{(0)}\partial_{\tau}u^{(0)}+\partial_{x}p^{(1)}$
$\displaystyle=0$ $\displaystyle\partial_{\tau}s^{(0)}$ $\displaystyle=0$

With the mixed diffusion scaling, the limit equations at the convective limit
match the Roe form, however at the acoustic limit all diffusion terms vanish
asymptotically, which leaves the scheme completely undamped. It should be
noted that the pressure diffusion in the pressure equation of the Roe form
(28) with the mixed diffusion scaling still stabilises the zeroth and first
order pressures $p^{(0,1)}$ even if it disappears from the limit equations at
the convective limit (34). This can be seen from an asymptotic expansion of
the modified or discrete equations, even if it is not evident from the limit
equations [21]. The pressure diffusion in the pressure equation is completely
absent from the first Toro-Vasquez form, which means that even though the
limit equations (46) match those of the Roe scheme, it is potentially
susceptible to instabilities which do not occur for the Roe scheme.

The limit equations for the second Toro-Vasquez form (27) are, like those for
the second Zha-Bilgen form, identical to the limit equations for the Roe form.
This not only removes the erroneous pressure anti-diffusion term from the
entropy limit equations, it also restores the pressure diffusion terms to the
pressure limit equations, which will stabilise the convective diffusion
scaling at the convective limit, the mixed diffusion scaling at the acoustic
limit, and will provide the proper damping on the acoustic variations for the
convective diffusion scaling at the acoustic limit.

One last point to mention is that that erroneous $\delta p$ terms in pressure
and entropy equations of the Liou-Steffen and Toro-Vasquez forms (29) and (31)
vanish by a factor of $M^{2}$ from the limit equations, so these terms are
unlikely to have any visible effect on the results found with these schemes.
On the other hand, the erroneous $\delta U$ terms in the entropy equations of
the first Zha-Bilgen and Toro-Vasquez forms (30) and (31) only vanish by a
factor of $M$ from the limit equations for the mixed diffusion scaling (42,43)
and (46,47). For moderately low Mach number, the pressure anti-diffusion terms
in the entropy equations of the first Zha-Bilgen and Toro-Vasquez forms with
mixed diffusion scaling may leave small but noticeable residues.

In this section we have shown that transforming to the entropy variables and
examining the limiting forms of the modified equations provides additional
insights into the structure of flux-vector-splittings at low Mach number. The
limit equations for the Liou-Steffen splitting exactly match those of the Roe
form, despite only matching one of three diffusion terms exactly in the
conserved variables. The Zha-Bilgen splitting matches the Roe form on the
pressure, velocity, and vorticity equations, up to a factor of $\gamma$ on the
pressure diffusion. However, with the convective diffusion scaling, a pressure
anti-diffusion term remains in the entropy equations, which may lead to
spurious entropy generation. The Toro-Vasquez splitting has no diffusion on
the pressure equation, which may lead to instabilities for both the convective
and acoustic limits for any diffusion scaling. Like for the Zha-Bilgen
splitting, the Toro-Vasquez splitting with the convective diffusion scaling
has a pressure anti-diffusion term on the entropy equation at both the
convective and acoustic limits. The second forms of the Zha-Bilgen and Toro-
Vasquez splittings match the limit equations of the Roe scheme.

##  4 Comparison of existing schemes

Scheme | Splitting | Scaling |  Comments
---|---|---|---
$\sigma_{u}\sim\mathcal{O}(1)$ | $\sigma_{a}\sim\mathcal{O}(1)$
Zha-Bilgen 1993 [69] | Zha-Bilgen | Acoustic |  Original Zha-Bilgen transonic scheme.
AUSM/AUSM+ 1993/1996 [37, 33] | Liou-Steffen | Acoustic |  Original Liou-Steffen transonic scheme.
LDFSS 1998/1999/2001 [14, 36, 13] | Liou-Steffen | Convective |  AUSM with elements of flux-preconditioning.
AUSM+-up 2006 [32] | Liou-Steffen | Convective |  Discrete asymptotic expansion to find scaling of diffusion coefficients.
SLAU 2009/11 [54] | Liou-Steffen | Mixed |  Numerically shows acoustic capability.
Toro-Vasquez 2012 [59] | Toro-Vasquez | Acoustic |  Original Toro-Vasquez transonic scheme.
Li & Gu 2013 [30] | Liou-Steffen | Mixed |  Modified AUSM+-up to mixed scheme according to their guidelines.
Sachdev et al. 2012 [48] | Liou-Steffen | Convective | Mixed |  Adaptive schemes for AUSM+-up and SLAU.
SLAU-WS 2013 [50] | Liou-Steffen | Mixed/convective |  Modifies the pressure diffusion in SLAU to increase by a factor of $M$ in the presence of pressure ‘wiggles’ independently of the timestep.
TV-MAS 2017 [56] | Toro-Vasquez | Mixed |  $U_{1/2}$ only adds upwind diffusion. $\underline{P}_{1/2}$ adds Zha-Bilgen $\delta U$ diffusion to the density and energy equations only, and adds Toro-Vasquez $\delta p$ diffusion to the momentum equations with mixed scaling and to the energy equation with the acoustic scaling.
TV-MAS2 2018 [31] | Toro-Vasquez | Mixed |  Reduces $\delta p$ term in energy equation of TVMAS to the mixed scaling.
E-AUSMPWAS 2018 [45] | Zha-Bilgen | Mixed |  Uses $\delta(pU)=p\delta U$, resulting in second Zha-Bilgen form of $\delta U$ and Liou-Steffen form of $\delta p$.
AUPM 2018 [4] | Zha-Bilgen | Mixed |  Second Zha-Bilgen form.
ZB-FVS-Corr 2020 [22] | Zha-Bilgen | Mixed |  Second Zha-Bilgen form.
AUSM-M 2020 [5] | Liou-Steffen | Mixed |  Reduces $\mu_{22}$ diffusion through shear layer.
TVAP 2021 [7] | Toro-Vasquez | Mixed |  $U_{1/2}$ only adds upwind diffusion. $\underline{P}_{1/2}$ adds Liou-Steffen $\delta U$ diffusion to the density and energy equations only, and adds Liou-Steffen $\delta p$ diffusion to momentum equations.
Table 3: Diffusion scaling of a number of existing low Mach number convection-
pressure flux-vector splitting schemes. Some schemes have different scalings
when the timestep is calculated using the convective CFL ($\sigma_{u}$) and
$M_{u}\sim\mathcal{O}(M)$, or using the acoustic CFL ($\sigma_{a}$) and
$M_{u}\sim\mathcal{O}(1)$.

In this section we compare the current analysis to previous studies and
schemes, and identify the diffusion scaling of a number of existing schemes
from the literature, which are shown in table 3 along with the diffusion
scaling they use. Almost all of the studies on convection-pressure flux-vector
splittings at low Mach number have focused on the Liou-Steffen splitting,
which reflects its popularity over the last 30 years, and all of these studies
were carried out in the conserved variables.∥∥∥Liou 2017 [34] swapped the
energy equation for the entropy equation in a number of schemes including AUSM
in the context of the overheating problem in high speed receding flow.

Liou & Edwards [14, 35, 36, 13] discuss the use of ideas from low Mach number
flux preconditioning for improving the performance of Liou-Steffen splitting
schemes at low Mach number, showing that the diffusion coefficients should
become independent of the Mach number as the incompressible regime is
approached. This is equivalent to the convective diffusion scaling in table 2,
which flux-difference splitting schemes using flux preconditioning also
approach [60, 20]. An overview of these developments, and of the AUSM family
of schemes more generally, is given by Edwards in [13, 12]. In the development
of AUSM+-up for all speeds [32] Liou goes into detail - including discrete
asymptotic analysis - on the design of an AUSM scheme for low Mach number
independently from the preconditioned system (although inspired by the earlier
work mentioned), and also arrives at the convective diffusion scaling.

Dellacherie 2010 and Dellacherie et al 2016 [11, 10] analysed the requirements
on Godunov schemes to accurately resolve convective flow features at low Mach
number, showing that the reduction from
$\mu_{22}\sim\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}1})$ for the acoustic
diffusion scaling to $\mu_{22}\sim\mathcal{O}(M^{0})$ for the mixed or
convective diffusion scalings is essential. It is also shown that AUSM/AUSM+
[37, 33] are unsuitable for low Mach number flow, having the acoustic
diffusion scaling on the pressure perturbation term, whereas AUSM+-up has the
correct scaling on $\mu_{22}$ so is suitable for low Mach number flow as
designed. By studying existing low Mach number Roe-type schemes, Li & Gu 2013
[30] draw out a set of guidelines for the design of such schemes. By applying
these guidelines to the AUSM+-up scheme, they modified the velocity
perturbation to change the diffusion scaling from convective to mixed. More
detail on how the approach used in the current paper and [21] relates to the
works of Dellacherie and Li & Gu can be found in [21].

The works of Venkateswaran and co-workers [64, 65, 44, 48] were, to the
authors’ knowledge, the earliest examples of the mixed diffusion scaling in
the literature, and some of the only studies to explicitly identify the mixed
diffusion scaling as a combination of the convective and acoustic diffusion
scalings which is suitable for both limits. In one of their more recent
studies, Sachdev et al. [48] show the two instabilities that the mixed
diffusion scaling is susceptible to: pressure chequerboards at the convective
limit and a velocity grid mode at the acoustic limit. They show that AUSM+-up
and SLAU [54] have the convective and mixed diffusion scalings respectively,
with the expected behaviour, and go on to show how they can both be modified
into adaptive schemes i.e. approaching the convective diffusion scaling as
$\Delta t\to\infty$ and approaching the mixed diffusion scaling as $\Delta
t\to 0$. This adaptive behaviour is achieved using an “unsteady Mach number”
$M_{u}=L_{\infty}/(\Delta ta_{\infty})$. Shima 2013 [50] takes a different
approach to resolving this problem, damping out pressure instabilities by
switching the diffusion scaling near pressure “wiggles” i.e. locations where
the ratio of fourth to second derivatives of the pressure field is large.

From table 3 it is clear that the Liou-Steffen splitting is the most commonly
used, although the Zha-Bilgen and Toro-Vasquez splittings have recently gained
more attention for low Mach number flow. Just as for Roe-type schemes, the
early Liou-Steffen schemes were rooted in the low Mach number preconditioning
to obtain the correct scaling, which results in the convective diffusion
scaling. The almost total switch to the mixed diffusion scaling also happened
around the same time as for Roe-type schemes c.2008 - Postdam et al. 2007
[44], Thornber et al. 2008 [57, 58] and Li & Gu 2008 [29] for Roe schemes, and
Shima & Kitamura 2009/11 (SLAU) [53, 54] for Liou-Steffen schemes. This means
that, to the authors’ knowledge, no Zha-Bilgen or Toro-Vasquez schemes using
the convective diffusion scaling exist in the literature. Of these earlier
Roe-type fluxes with the mixed diffusion scaling, only Potsdam et al. had the
explicit aim of creating a scheme capable of resolving both convective and
acoustic low Mach number features, whereas aeroacoustic simulations were one
of the original target applications for SLAU.

The AUSM+-up scheme is one of the most popular convection-pressure flux-vector
splitting schemes for low Mach number. As stated earlier the diffusion
coefficients become independent of the Mach number as $M\to 0$, which means
that it has the convective diffusion scaling. As such, it can accurately
resolve convective flow features in low Mach number flow [32, 11] and is
robust against the pressure chequerboard instabilities that can occur for
collocated schemes for low Mach number or incompressible flow. On the other
hand, it damps acoustic variations very quickly [40, 41], and suffers from a
spectral radius which scales as $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}2})$
instead of the $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}1})$ scaling of the
mixed and acoustic diffusion scalings and the Jacobian of the physical flux
vector [21]. The $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}2})$ spectral radius
scaling of the convective diffusion scaling has long been known in the context
of Roe-type schemes since it was first proven by Birken & Meister 2005 [2].
This spectral radius means that the stability bound on the timestep for an
explicit scheme reduces prohibitively as $\Delta t\sim\mathcal{O}(M^{2})$, a
factor of $M$ faster than usual.

Birken & Meister noted that other schemes with the convective diffusion
scaling would also suffer from these problems, and gave the specific example
of the AUSMDV scheme [66]. That these findings also apply to AUSM+-up appears
to have gone unrecognised (or at least unmentioned) in the more recent
literature, although evidence of it can be seen in a number of studies. In the
original paper [32] it can be seen from figure 22 that the convergence of the
implicit scheme stalls after at least four orders of magnitude decrease in the
residual, whereas convergence to machine zero is achieved when low Mach number
preconditioning is used. This is consistent with the expectation that the
spectral radius and, more importantly, the condition number of the system is
reduced from $\mathcal{O}(M^{\scalebox{0.4}[1.0]{$-$}2})$ to
$\mathcal{O}(M^{0})$ by the preconditioning. Several papers have reported
stability issues with AUSM+-up at very low Mach number when using an explicit
timestepping strategy [39, 25, 6]. The response to these issues has been to
prevent the Mach number used in the calculation of the artificial diffusion
falling below some cut-off Mach number $M_{co}>M_{\infty}$. Keeping $M_{co}$
above the actual Mach number of the flow has two effects. Firstly, it reduces
the $\mu_{11}$ pressure diffusion, which reduces the severe scaling of the
spectral radius and alleviates the stability restrictions as intended.
Secondly, it increases the $\mu_{22}$ velocity diffusion, which results in a
more diffusive solution approaching that found with the acoustic diffusion
scaling. Li & Gu [30] mention that the AUSM+-up scheme suffers from the cut-
off Mach number problem just as the preconditioned Roe scheme does, but do not
explicitly connect this back to the spectral radius scaling and the results of
Birken & Meister. More recently, Chen et al. 2018 [6] showed empirically that
the scaling of the maximum stable timestep of AUSM+-up is what would be
expected of the convective diffusion scaling,******This can be seen from table
1 and figure 5 in [6] which show linear scaling of the maximum timestep with
the Mach number. This scaling is found by holding the speed of sound constant
and varying the convective velocity to control the Mach number, so $\Delta
t\sim\mathcal{O}(M)$ is equivalent to $\Delta t\sim\mathcal{O}(M^{2})$ when
holding the convective velocity constant and varying the speed of sound.
Incidentally, the authors also empirically found this timestep scaling for
AUSM+-up, independently of [6]. Investigating if and how it was connected to
the results of Birken & Meister was one of the original motivations for [21]
and the current study. but again do not remark on the connection to the
preconditioned Roe scheme or the results of Birken & Meister.

The remedies to the stability issues of the convective scheme which do not
sacrifice accuracy are the same as for the Roe-type schemes. The most common
solution is to use low Mach number preconditioning, which results in a
$\mathcal{O}(M^{0})$ spectral radius on the convective timescale, but destroys
the time accuracy of the scheme [60]. When time accuracy is required, an
implicit scheme can be used, often in conjunction with dual-time stepping
[24], which allows the use of low Mach number preconditioning while
maintaining time accuracy [63, 51, 52]. Using a scheme with the mixed
diffusion scaling would also avoid this issue, but at the price of losing the
robustness against pressure chequerboard instabilities. It should be mentioned
that preconditioning does not entirely remove the need for some cut-off Mach
number because, although the condition number is reduced which improves the
linear long-time convergence, the eigenvectors of the preconditioned system
become increasingly non-normal as $M\to 0$ which impedes the non-linear short-
time convergence [9]. In this case the cut-off Mach number should be applied
only to the preconditioning matrix and $\mu_{11}$ to keep the correct scaling
of $\mu_{22}$.

Turning now to the Zha-Bilgen and Toro-Vasquez splittings, the first thing to
note is that none of the schemes in table 3 use the first forms where
$\delta(pU)=U\delta p$, which gives $\delta\underline{P}=\delta
p\underline{\Pi}$ and maintains the separation of the $\delta p$ diffusion in
the pressure perturbation term and the $\delta U$ diffusion in the velocity
perturbation term. Chen et al. 2018 (AUPM) [4] and Iampietro et al. 2020 [22]
both use the second form of the Zha-Bilgen splitting with $\delta(pU)=p\delta
U+U\delta p$, which gives diffusion equivalent to the Roe scheme on all three
diffusion terms - convective upwinding, and the velocity and pressure
perturbations. Qu et al. 2018 (E-AUSMPWAS) [45] use $\delta(pU)=p\delta U$,
which is equivalent to using the Liou-Steffen form for both the velocity and
pressure perturbation diffusion terms. According to the findings of section 3,
the differences between different forms of the $\delta p$ pressure
perturbation terms are much less important than the form of the velocity
perturbation term, so the difference in behaviour between these two
expressions for $\delta(pU)$ should be small.

In the three Toro-Vasquez schemes in table 3, Sun et al. 2017 (TV-MAS) [56],
Lin et al. 2018 (TV-MAS2) [31], and Chen et al. 2021 (TVAP) [7], the
convection term $U_{1/2}\underline{\varphi}$ adds only convective upwind
diffusion, and not the velocity perturbation $\delta U$ diffusion. Instead,
all of the $\delta U$ and $\delta p$ diffusion terms are included through the
pressure perturbation term $\delta\underline{P}$. $\delta U$ diffusion terms
are added only to the density and energy equations, with the Liou-Steffen form
used in the energy equation. TVAP modifies this slightly, with the static
enthalpy $(\rho h\delta U)$ used instead of the total enthalpy $(\rho E\delta
U)$. None of the TV schemes apply $\delta U$ diffusion to the momentum
equations. In their review of low Mach Roe-type schemes, Li & Gu 2013 [30]
note that inconsistent application of the velocity perturbation diffusion
across the different equations reduces the effectiveness of the convective
diffusion scaling at suppressing pressure chequerboards. Whether the same is
true for the Toro-Vasquez flux with the mixed diffusion scaling has not yet
been addressed.

TVAP adds a $\delta p$ diffusion term only to the momentum equations, so in
effect uses the Liou-Steffen form for both the velocity and pressure
perturbation diffusion terms, albeit with no $\delta U$ diffusion on the
momentum equations. TV-MAS keeps the Toro-Vasquez form of the $\delta p$ terms
in the momentum and energy equations, but uses the mixed scaling for the
$\mu_{22}$ coefficient in the momentum equations and the acoustic scaling in
the energy equations. Lin et al. [31] show that having different $\mu_{22}$
scalings in each equation results in the density variations scaling as
$\nabla\rho\sim\mathcal{O}(M)$, instead of the expected
$\nabla\rho\sim\mathcal{O}(M^{2})$ for isentropic flow, even when the pressure
variations are correct at $\nabla p\sim\mathcal{O}(M^{2})$. They rectify this
in TV-MAS2, and show that two low Mach number Roe-type schemes with the mixed
diffusion scaling also require this correction [58, 15].

In this section we have discussed the diffusion scaling and form of a number
of low Mach number convection-pressure flux-vector splitting schemes from the
literature. The Liou-Steffen splitting is the most popular, although schemes
using the Zha-Bilgen and Toro-Vasquez splittings have also appeared recently.
All of the Liou-Steffen and Zha-Bilgen schemes surveyed match one of the
general forms used here. On the other hand, the Toro-Vasquez schemes do not
match exactly, applying the velocity perturbation only to the density and
energy equations, although the current analysis is presumed to still give a
reasonable indication of their low Mach number behaviour. Of the Zha-Bilgen
and Toro-Vasquez schemes surveyed, every one has $\delta U$ diffusion
equivalent to the Liou-Steffen and Roe forms on the energy equation. The
majority of schemes use the mixed diffusion scaling, especially in the more
recent literature. One notable exception to this is the AUSM+-up scheme which
uses the convective diffusion scaling, so should be treated as asymptotically
equivalent to the preconditioned Roe scheme as $M\to 0$.

##  5 Numerical examples

Three numerical examples are presented to demonstrate the behaviour of the
different splittings with both the convective and mixed diffusion scalings for
convective and acoustic flows. Results are shown for both the first and second
forms of the Zha-Bilgen and Toro-Vasquez forms, with particular attention paid
to the predicted entropy field for each example to show the effect of the
erroneous diffusion terms found in section 3.

The same finite volume code used in [21] is used here, with the only
difference being the interface flux functions. These are constructed using
(12), (13) and (10). $\underline{\varphi}$ is taken as the upwind value based
on the sign of $U_{1/2}$, and the form of $\delta\underline{P}$ is as
described for each splitting in section 3.1. The scheme is first order in
space, and uses first order Euler forward integration in time. Note that this
means that $1/M$ times as many timesteps are required for the convective
scheme due to the increased spectral radius of the diffusion. An ideal gas law
is used for all cases, with $\gamma=1.4$ and $R=287.058$.

###  5.1 One dimensional examples

####  5.1.1 Isolated soundwave

(a)

(b)

Figure 1: Non-dimensional gauge pressure profiles for an isolated soundwave at
$M=0.01$ after a single period using (a) convective diffusion scaling
$\underline{\underline{A}}^{c}$ (b) mixed diffusion scaling
$\underline{\underline{A}}^{m}$.

The performance of the discrete schemes for purely acoustic flow is
demonstrated using the isolated soundwave test case from [21]. The initial
profile, shown in figure 1, is a sinusoidal right travelling isentropic
soundwave with a background flow at $M=0.01$ and 16 points per wavelength.
Pressure profiles after a single period are shown for each of the flux-
splittings with the convective and mixed diffusion scaling in figures 1(a) and
1(b) respectively, with the Roe scheme results shown for comparison.

With the convective diffusion scaling, all schemes except the first Toro-
Vasquez form quickly damp the acoustic variations and closely match the
results with the Roe scheme, as predicted by the limit equations (37) and
(41). The first form of the Toro-Vasquez scheme does not damp the acoustic
wave at all due to the lack of pressure diffusion on the pressure equation.
The results with the mixed diffusion scaling, shown in figure 1(b), also agree
with the findings of section 3. The Liou-Steffen scheme and the second Zha-
Bilgen and Toro-Vasquez forms match the Roe scheme almost exactly, and the
first Zha-Bilgen form has similar resolution but slightly reduced diffusion.
The first Toro-Vasquez form predicts the soundwave amplitude growing, which is
unsurprising given that the scheme asymptotically approaches a central scheme
for acoustic flow, and a first order Euler forward discretisation is used for
the time derivative.

Entropic profiles of $s=p/\rho^{\gamma}$ predicted by the different schemes
are shown in figure 2. The Roe and Liou-Steffen schemes preserve the
isentropic nature of the soundwave for both convective and mixed diffusion
scalings, as predicted by the limit equation analysis. The second Zha-Bilgen
and Toro-Vasquez forms preserve the isentropic profile to a similar degree,
but have fractionally more entropy generation just visible in figure 2(b).
However, the first Zha-Bilgen and Toro-Vasquez forms both generate spurious
entropy variations. For the convective diffusion scaling, this can be seen
from the limit equations (41) and (45) to be due to the erroneous pressure
diffusion term in the acoustic entropy. This term is a factor of
$\gamma/(\gamma-1)$ larger in the Toro-Vasquez splitting and $\nabla^{2}p$ is
larger because the soundwave is undamped, which results in a higher entropy
generation for this scheme. The variations are less than $1\%$ of the
background value, but this is after the soundwave has travelled just a single
wavelength. With the mixed diffusion scaling, these erroneous pressure
diffusion terms are one order of $M$ smaller than the leading order terms in
the entropy equations of these two schemes, so disappear from the limit
equations (43) and (47). Figure 2(b) shows that the entropy variations
predicted with the mixed diffusion scaling are a factor of $M$ smaller than
those predicted by the convective diffusion scaling, consistent with this
scaling of the pressure diffusion terms. The entropy profiles are essentially
unchanged if the test case is rerun with the $\delta p$ terms in the interface
pressure perturbation $\delta\underline{P}$ switched off, which further
reinforces that the entropy generation of the first Zha-Bilgen and Toro-
Vasquez forms is due to the form of the $\delta U$ diffusion terms in the
energy equation.

(a)

(b)

Figure 2: Non-dimensional entropy profiles for an isolated soundwave at
$M=0.01$ after a single period using (a) convective diffusion scaling
$\underline{\underline{A}}^{c}$ (b) mixed diffusion scaling
$\underline{\underline{A}}^{m}$.

###  5.2 Two dimensional examples

####  5.2.1 Circular cylinder

The next test case is inviscid flow at $M=0.01$ around a circular cylinder,
which will demonstrate the performance of the different schemes for purely
convective stationary flow. The setup is identical to that in [21]. A
stretched O-grid is used with the first cell height of $0.036r$, where $r$ is
the cylinder radius, and the farfield boundary is $50r$ from the centre of the
cylinder. The farfield boundary conditions are enforced by setting ghost cell
values using upstream velocity and entropy and downstream static pressure. At
the cylinder wall, curvature corrected boundary conditions are used which
reduce spurious entropy production at the wall [8]. Convergence to steady
state is accelerated using local timestepping and Weiss & Smith’s low Mach
number preconditioner [68].

The profiles of the gauge pressure, the pressure error relative to the
analytical potential flow solution, the entropy, and the density using the Roe
and Liou-Steffen schemes with the convective artificial diffusion scaling are
shown in figures 3 and 4. As for the isolated soundwave, the Liou-Steffen
results are almost indistinguishable from the Roe results. The numerical
results are close to the potential solution, with the pressure variations 3(a)
and 4(a) on the correct order of magnitude. The pressure loss downstream of
the cylinder is due to the fact a first order scheme has been used, but the
smooth error profiles 3(b) and 4(b) indicate that the pressure is free of
chequerboard modes. There is a very small amount of entropy generated at the
wall visible in 3(c) and 4(c) \- which is unavoidable with standard upwind
schemes, although is reduced by the curvature corrected boundary conditions
(see [16] for a recent example of low Mach number entropy conserving/stable
schemes) - and the density profiles 3(d) and 4(d) resemble the pressure
profiles as expected.

The results for the first Zha-Bilgen splitting with the convective diffusion
scaling are shown in figure 5. The pressure field 5(a) is essentially the same
as that of the Roe and Liou-Steffen schemes, although the error plot 5(b)
shows that the pressure loss downstream of the cylinder is marginally less due
to the reduced diffusion in the pressure equation. However, the entropy
profile 5(c) shows large unphysical variations due to the pressure anti-
diffusion term in the entropy equation, which peak at around $30\%$ of the
freestream values close to the wall where the pressure Poisson $\nabla^{2}p$
is largest. The effect of this entropy generation can be seen in the density
profile 5(d), where large variations can also be seen close to the wall. The
results for the second Zha-Bilgen splitting with the convective diffusion
scaling are shown in figure 6. The pressure and pressure error fields are
essentially indistinguishable from the first form, but the entropy and density
fields show clear improvement. The entropy variations have been reduced to
approximately the same level as the Roe and Liou-Steffen schemes, and the
density profile is vastly improved.

The results for the first Toro-Vasquez form with the convective diffusion
scaling are shown in figure 7. For this case, the previous schemes converged
in around $25,000-30,000$ iterations, however the first Toro-Vasquez form did
not converge, with NaNs occuring in the solution after around $42,000$
iterations. Results are shown here after $30,000$ iterations. Significant
radial chequerboard modes can be seen in the pressure and pressure error
profiles 7(a) and 7(b). The density and entropy profiles 7(c) and 7(d) also
have radial chequerboard modes and large errors close to the wall due to the
pressure anti-diffusion in the entropy equation. These errors grow until the
solution eventually diverges and breaks down. The major difference between the
limit equations of the first Zha-Bilgen and Toro-Vasquez forms with the
convective diffusion scaling is the diffusion term in the pressure equation
which, by comparing figures 5 and 7, implies that this term somewhat damps the
errors arising from the spurious entropy generation. The results for the
second Toro-Vasquez form with the convective diffusion scaling are shown in
figure 8. As for the Zha-Bilgen scheme, the entropy and density fields show a
marked improvement on the first form, although the entropy variations are
about twice that seen for the Roe and Liou-Steffen schemes.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 3: Flow around a circular cylinder at $M=0.01$ using the Roe scheme
with convective diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 4: Flow around a circular cylinder at $M=0.01$ using the Liou-Steffen
scheme with convective diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 5: Flow around a circular cylinder at $M=0.01$ using the first Zha-
Bilgen scheme with convective diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 6: Flow around a circular cylinder at $M=0.01$ using the second Zha-
Bilgen scheme with convective diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 7: Flow around a circular cylinder at $M=0.01$ using the first Toro-
Vasquez scheme with convective diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 8: Flow around a circular cylinder at $M=0.01$ using the second Toro-
Vasquez scheme with convective diffusion scaling.

The results for the Roe and Liou-Steffen schemes with the mixed diffusion
scaling are shown in figures 9 and 10 respectively. The pressure and density
plots are very similar to the results found with the convective diffusion
scaling, and the entropy generated at the walls is even lower. Inspection of
figures 9(b) and 10(b) reveals that the error profiles are no longer smooth,
and a slight chequerboard exists in both solutions, as is expected for schemes
with the mixed diffusion scaling.

The results for the first Zha-Bilgen form in figure 11 with the mixed
diffusion scaling also show that the pressure field is generally accurate, but
a chequerboard mode also exists in the solution which is visible in 11(b). The
entropy and density variations 11(c) and 11(d) are approximately a factor of
$M$ smaller than for the convective scaling 5(c) and 5(d), although are still
much larger than those found with the Roe and Liou-Steffen schemes. The
results for the second Zha-Bilgen form are shown in figure 12 and closely
resemble those of the Roe and Liou-Steffen schemes, including the pressure
chequerboard and the level of the entropy variations.

The results for the first Toro-Vasquez form with mixed diffusion scaling in
figure 13 show a vast improvement on the results with the convective scaling
in figure 7. The pressure field much more closely matches the potential
solution, although still with a chequerboard mode visible in 13(b). As for the
Zha-Bilgen scheme, the entropy and density variations are much reduced
compared to the convective scaling results - especially the entropy variations
- while still being larger than the Roe and Liou-Steffen solutions. The
results for the second Toro-Vasquez form with mixed diffusion scaling are
shown in figure 14. As for with the convective diffusion scaling, the results
are a significant improvement on the first form, although the entropy
generation is still around twice that of the Roe and Liou-Steffen schemes.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 9: Flow around a circular cylinder at $M=0.01$ using the Roe scheme
with mixed diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 10: Flow around a circular cylinder at $M=0.01$ using the Liou-Steffen
scheme with mixed diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 11: Flow around a circular cylinder at $M=0.01$ using the first Zha-
Bilgen scheme with mixed diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 12: Flow around a circular cylinder at $M=0.01$ using the second Zha-
Bilgen scheme with mixed diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 13: Flow around a circular cylinder at $M=0.01$ using the first Toro-
Vasquez scheme with mixed diffusion scaling.

(a) Gauge pressure

(b) Pressure error

(c) Entropy

(d) Density

Figure 14: Flow around a circular cylinder at $M=0.01$ using the second Toro-
Vasquez scheme with mixed diffusion scaling.

####  5.2.2 Isentropic vortex

The final test case is an inviscid stationary Gresho vortex [17], which will
be used to further demonstrate the impact of the form of the entropy limit
equations. The maximum velocity of the vortex $u_{m}$ is at a radius $r$ and
Mach number $M=0.01$. The initial conditions can be found in the soundwave-
vortex interaction test case of [21], with the exception that in the present
study the density profile is calculated using $d\rho=dp/a^{2}$ to enforce
isentropic conditions instead of using a constant value. The domain is
$10r\times 10r$, and is periodic in both $x$ and $y$. Results are shown after
a quarter turn of the vortex, i.e. $t=\frac{\pi r}{2u_{m}}$.

The pressure profile on a line through the centre of the vortex predicted by
each of the schemes is shown in figure 15, where it can be seen that all of
the schemes diffuse the vortex to some extent. With the convective diffusion
scaling, the Roe, Liou-Steffen, and second Zha-Bilgen and Toro-Vasquez schemes
match once again, the solution from the first Zha-Bilgen scheme is slightly
less diffusive, and the first Toro-Vasquez scheme is less diffusive still, in
agreement with the limit equation analysis. The limit equations with the mixed
diffusion scaling at the convective limit are the same for all four schemes,
and figure 15(b) shows that the predicted pressure profiles are all almost
identical to each other. The solution with the first Toro-Vasquez scheme with
the mixed diffusion scaling broke down and is not shown.

The entropy profiles found with the convective and mixed diffusion scalings
are shown in figures 16 and 17 respectively. The Roe, Liou-Steffen, and second
Zha-Bilgen schemes generate very little entropy, and the second Toro-Vasquez
generates only a small amount more. For the Roe and the second Zha-Bilgen and
Toro-Vasquez schemes, the entropy generation is predominantly at the four
corners where the flow is the most misaligned with the grid. The Liou-Steffen
splitting generates slightly more entropy than the Roe scheme, which is the
opposite of the trend seen in the circular cylinder test case. The profiles in
16(a) and 17(a) are similar to the entropy production modes of an entropy
stable Roe scheme shown in figure 11 of [16] smeared by a quarter turn, and
the entropy production after a single timestep closely resembles the
production mode. The entropy profiles found with the first Zha-Bilgen and
Toro-Vasquez forms and the convective scaling have much larger variations, and
their shape is consistent with what would be expected from the limit equation
(40) i.e. there appears to be an entropy sink where $\nabla^{2}p>0$ in the
vortex core, and an entropy source where $\nabla^{2}p<0$ at the edge of the
vortex. The entropy profile for the first Zha-Bilgen form with the mixed
diffusion scaling is qualitatively very similar to that with the convective
diffusion scaling, except with the variations reduced by around an order of
$M$. However, the solution for the first Toro-Vasquez for with the mixed
diffusion scaling has entirely broken down, as can be seen from figure 17(e),
with only very high wavenumber oblique modes remaining. As with the previous
examples, the second forms of Zha-Bilgen and Toro-Vasquez significantly
improve the solution compared to the first forms, with the entropy field for
the second Zha-Bilgen form almost the same as with the Roe scheme. The
benefits of the second form are especially clear for the Toro-Vasquez scheme
with mixed diffusion scaling, preventing the total breakdown that was seen
with the first form, with the entropy variations on the same order of
magnitude as for the Roe scheme.

(a)

(b)

Figure 15: Non-dimensional gauge pressure profiles for a Gresho vortex at
$M=0.01$ after $1/4$ rotation using (a) convective diffusion scaling
$\underline{\underline{A}}^{c}$ (b) mixed diffusion scaling
$\underline{\underline{A}}^{m}$. $p_{\text{min}}$ is the minimum pressure at
the centre of the vortex core at $t=0$.

(a) Roe

(b) Liou-Steffen

(c) Zha-Bilgen 1

(d) Zha-Bilgen 2

(e) Toro-Vasquez 1

(f) Toro-Vasquez 2

Figure 16: Non-dimensional entropy profiles for a Gresho vortex at $M=0.01$
after $1/4$ rotation using the convective diffusion scaling
$\underline{\underline{A}}^{c}$.

(a) Roe

(b) Liou-Steffen

(c) Zha-Bilgen 1

(d) Zha-Bilgen 2

(e) Toro-Vasquez 1

(f) Toro-Vasquez 2

Figure 17: Non-dimensional entropy profiles for a Gresho vortex at $M=0.01$
after $1/4$ rotation using the mixed diffusion scaling
$\underline{\underline{A}}^{m}$.

We have shown three numerical examples which demonstrate the behaviour of the
three flux-vector-splittings for both acoustic and convective flow. Most
schemes have the expected property that the mixed diffusion scaling is
suitable for both convective and acoustic flow, whereas the convective
diffusion scaling will only properly resolve convective flow and will quickly
damp out acoustic waves. However, the first Toro-Vasquez form will not damp
out acoustic waves even with the convective scaling. All three schemes suffer
from pressure chequerboard modes with the mixed scaling, but the convective
diffusion scaling resolves this issue for all schemes except the first Toro-
Vasquez scheme. The most significant result is the demonstration that the
erroneous pressure anti-diffusion terms in the entropy variable limit
equations of the Zha-Bilgen and Toro-Vasquez schemes do in fact degrade the
solutions found with these schemes. On the other hand, the Liou-Steffen
results almost exactly match the results found with the Roe-type scheme for
all cases and with both the convective and mixed diffusion scaling, as
predicted by the limit equations analysis.

##  6 Conclusions

In this study. we have analysed the form of the artificial diffusion at low
Mach number for three convection-pressure flux-vector splittings: Liou-
Steffen, Zha-Bilgen, and Toro-Vasquez. The approximate diffusion form of the
Liou-Steffen and Zha-Bilgen splittings bear a close resemblance to that of the
Roe scheme, whereas the Toro-Vasquez form does not. We identified two forms of
the energy equation component of the pressure perturbation term for the Zha-
Bilgen and Toro-Vasquez splittings - the first form results in diffusion terms
with a close correspondence to the original flux-vector splitting, whereas the
second results in diffusion terms closer to the Roe and Liou-Steffen forms.

We then transformed the artificial diffusion to the entropy variables to gain
more insight into the differences between the splittings. The limit equations
of the Liou-Steffen splitting and the second forms of the Zha-Bilgen and Toro-
Vasquez splittings are identical to those of the Roe scheme for both
convective and acoustic flow with both the convective and mixed diffusion
scalings. One the other hand, the first forms of the Zha-Bilgen and Toro-
Vasquez splittings were found to have erroneous pressure anti-diffusion terms
on the entropy equation, which remain in the limit equations with the
convective diffusion scaling at both the convective and acoustic limits. The
first form of the Toro-Vasquez splitting also completely lacks the usual
pressure diffusion in the pressure equation.

Examining low Mach number convection-pressure flux-vector splitting schemes in
the literature, we see that the Liou-Steffen splitting is clearly the most
popular, and that almost all modern schemes of all splittings use the mixed
diffusion scaling. All of the Liou-Steffen and Zha-Bilgen schemes reviewed had
a form matching the general form assumed in this paper. The Toro-Vasquez
schemes, while mostly matching, apply the velocity perturbation diffusion
inconsistently by applying it only to the density and energy equations and not
the momentum equations. All Zha-Bilgen and Toro-Vasquez schemes reviewed used
the second form of the diffusion which more closely resembles the Roe and
Liou-Steffen forms, instead of the first forms that results in the erroneous
entropy generation terms.

Three numerical examples were used to verify the findings of the analysis,
showing excellent agreement between prediction and results for all the
splittings at both the convective and acoustic limits, including the spurious
entropy generation of the first Zha-Bilgen and Toro-Vasquez forms and the
significant improvements when using the second forms.

The differences between the three flux-splitting schemes considered here show
that the diffusion on the energy equation is an important factor in the design
of numerical schemes for low Mach number. In particular, it is crucial to
obtain the correct form of the $\delta U$ diffusion term - the spurious
entropy generation of the first Zha-Bilgen and Toro-Vasquez forms is entirely
due to this term, irrespective of the form of the $\delta p$ diffusion term.
This is also relevant for future studies of low Mach number numerical schemes
because it is not uncommon for the barotropic Euler equations to be used for
such studies, which removes the energy equation by enforcing additional
constraints on the equation of state. The barotropic equations are simpler and
allow for more in-depth mathematical analysis, but clearly care must be taken
when extending the findings of such analysis to the full Euler equations where
the energy equation must also be considered. Future extensions of convection-
pressure flux-vector splittings could include alternative methods for ensuring
correct resolution of the entropy field, the effect of including $\delta U$ in
deciding the upwind direction, and other forms for the $\delta(pU)$ term in
the pressure perturbation.

## Acknowledgements

The authors gratefully acknowledge support from the EPSRC Center for Doctoral
Training in Gas Turbine Aerodynamics and Rolls-Royce plc under grant number
EP/L015943/1.

## References

  * [1] Wasilij Barsukow  “Truly multi-dimensional all-speed schemes for the Euler equations on Cartesian grids”  In _Journal of Computational Physics_ 435, 2021, pp. 110216  DOI: 10.1016/j.jcp.2021.110216
  * [2] P. Birken and A. Meister  “Stability of Preconditioned Finite Volume Schemes at Low Mach Numbers”  In _BIT Numerical Mathematics_ 45.3, 2005, pp. 463–480  DOI: 10.1007/s10543-005-0009-0
  * [3] Pascal Bruel, Simon Delmas, Jonathan Jung and Vincent Perrier  “A low Mach correction able to deal with low Mach acoustics”  In _Journal of Computational Physics_ 378, 2019, pp. 723–759  DOI: 10.1016/j.jcp.2018.11.020
  * [4] Shu-sheng Chen et al.  “A novel flux splitting scheme with robustness and low dissipation for hypersonic heating prediction”  In _International Journal of Heat and Mass Transfer_ 127, 2018, pp. 126–137  DOI: 10.1016/j.ijheatmasstransfer.2018.06.121
  * [5] Shu-sheng Chen et al.  “An improved AUSM-family scheme with robustness and accuracy for all Mach number flows”  In _Applied Mathematical Modelling_ 77, 2020, pp. 1065–1081  DOI: 10.1016/j.apm.2019.09.005
  * [6] Shu-sheng Chen, Chao Yan, Shuai Lou and Bo-xi Lin  “An improved entropy-consistent Euler flux in low Mach number”  In _Journal of Computational Science_ 27, 2018, pp. 271–283  DOI: 10.1016/j.jocs.2018.06.006
  * [7] Shusheng Chen et al.  “A low-diffusion robust flux splitting scheme towards wide-ranging Mach number flows”  In _Chinese Journal of Aeronautics_ 34.5, 2021, pp. 628–641  DOI: 10.1016/j.cja.2020.12.010
  * [8] A. Dadone and B. Grossman  “Surface boundary conditions for the numerical solution of the Euler equations” Publisher: American Institute of Aeronautics and Astronautics _eprint: https://doi.org/10.2514/3.11983  In _AIAA Journal_ 32.2, 1994, pp. 285–293  DOI: 10.2514/3.11983
  * [9] D.. Darmofal and P.. Schmid  “The Importance of Eigenvectors for Local Preconditioners of the Euler Equations”  In _Journal of Computational Physics_ 127.2, 1996, pp. 346–362  DOI: 10.1006/jcph.1996.0180
  * [10] S. Dellacherie, J. Jung, P. Omnes and P.-A. Raviart  “Construction of modified Godunov-type schemes accurate at any Mach number for the compressible Euler system” Publisher: World Scientific Publishing Co.  In _Mathematical Models and Methods in Applied Sciences_ 26.13, 2016, pp. 2525–2615  DOI: 10.1142/S0218202516500603
  * [11] Stéphane Dellacherie  “Analysis of Godunov type schemes applied to the compressible Euler system at low Mach number”  In _Journal of Computational Physics_ 229.4, 2010, pp. 978–1016  DOI: 10.1016/j.jcp.2009.09.044
  * [12] J.. Edwards  “Reflections on the early development of the “AUSM family” of Riemann solvers”  In _Shock Waves_ 29.5, 2019, pp. 601–609  DOI: 10.1007/s00193-018-0863-8
  * [13] Jack Edwards  “Towards unified CFD simulations of real fluid flows”  In _15th AIAA Computational Fluid Dynamics Conference_ , Fluid Dynamics and Co-located Conferences  American Institute of AeronauticsAstronautics, 2001  DOI: 10.2514/6.2001-2524
  * [14] Jack R. Edwards and Meng-Sing Liou  “Low-Diffusion Flux-Splitting Methods for Flows at All Speeds” Publisher: American Institute of Aeronautics and Astronautics  In _AIAA Journal_ 36.9, 1998, pp. 1610–1617  DOI: 10.2514/2.587
  * [15] Philippe Fillion, Augustin Chanoine, Stéphane Dellacherie and Anela Kumbaro  “FLICA-OVAP: A new platform for core thermal–hydraulic studies”  In _Nuclear Engineering and Design_ 241.11, 13th International Topical Meeting on Nuclear Reactor Thermal Hydraulics (NURETH-13), 2011, pp. 4348–4358  DOI: 10.1016/j.nucengdes.2011.04.048
  * [16] Ayoub Gouasmi, Scott M. Murman and Karthik Duraisamy  “Entropy-Stable Schemes in the Low-Mach-Number Regime: Flux-Preconditioning, Entropy Breakdowns, and Entropy Transfers”  In _Journal of Computational Physics_ , 2022, pp. 111036  DOI: 10.1016/j.jcp.2022.111036
  * [17] Philip M. Gresho and Stevens T. Chan  “On the theory of semi-implicit projection methods for viscous incompressible flow and its implementation via a finite element method that also introduces a nearly consistent mass matrix. Part 2: Implementation” _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/fld.1650110510  In _International Journal for Numerical Methods in Fluids_ 11.5, 1990, pp. 621–659  DOI: 10.1002/fld.1650110510
  * [18] H. Guillard and B. Nkonga  “Chapter 8 - On the Behaviour of Upwind Schemes in the Low Mach Number Limit: A Review”  In _Handbook of Numerical Analysis_ 18, Handbook of Numerical Methods for Hyperbolic Problems  Elsevier, 2017, pp. 203–231  DOI: 10.1016/bs.hna.2016.09.002
  * [19] Hervé Guillard and Angelo Murrone  “On the behavior of upwind schemes in the low Mach number limit: II. Godunov type schemes”  In _Computers & Fluids_ 33.4, 2004, pp. 655–675  DOI: 10.1016/j.compfluid.2003.07.001
  * [20] Hervé Guillard and Cécile Viozat  “On the behaviour of upwind schemes in the low Mach number limit”  In _Computers & Fluids_ 28.1, 1999, pp. 63–86  DOI: 10.1016/S0045-7930(98)00017-6
  * [21] Joshua Hope-Collins and Luca Mare  “Artificial diffusion for convective and acoustic low Mach number flows I: Analysis of the modified equations, and application to Roe-type schemes”  In _Journal of Computational Physics_ , 2022, pp. 111858  DOI: 10.1016/j.jcp.2022.111858
  * [22] D. Iampietro, F. Daude and P. Galon  “A low-diffusion self-adaptive flux-vector splitting approach for compressible flows”  In _Computers & Fluids_ 206, 2020, pp. 104586  DOI: 10.1016/j.compfluid.2020.104586
  * [23] Antony Jameson  “Artificial diffusion, upwind biasing, limiters and their effect on accuracy and multigrid convergence in transonic and hypersonic flows”  In _11th Computational Fluid Dynamics Conference_ , Fluid Dynamics and Co-located Conferences  American Institute of AeronauticsAstronautics, 1993  DOI: 10.2514/6.1993-3359
  * [24] Antony Jameson  “Time dependent calculations using multigrid, with applications to unsteady flows past airfoils and wings”  In _10th Computational Fluid Dynamics Conference_ , Fluid Dynamics and Co-located Conferences  American Institute of AeronauticsAstronautics, 1991  DOI: 10.2514/6.1991-1596
  * [25] Keiichi Kitamura and Atsushi Hashimoto  “Reduced dissipation AUSM-family fluxes: HR-SLAU2 and HR-AUSM+-up for high resolution unsteady flow simulations”  In _Computers & Fluids_ 126, 2016, pp. 41–57  DOI: 10.1016/j.compfluid.2015.11.014
  * [26] Keiichi Kitamura, Meng-Sing Liou and Chih-Hao Chang  “Extension and Comparative Study of AUSM-Family Schemes for Compressible Multiphase Flow Simulations” Publisher: Cambridge University Press  In _Communications in Computational Physics_ 16.3, 2014, pp. 632–674  DOI: 10.4208/cicp.020813.190214a
  * [27] Keiichi Kitamura and Eiji Shima  “Towards shock-stable and accurate hypersonic heating computations: A new pressure flux for AUSM-family schemes”  In _Journal of Computational Physics_ 245, 2013, pp. 62–83  DOI: 10.1016/j.jcp.2013.02.046
  * [28] Bram Leer  “Flux-vector splitting for the Euler equations”  In _Eighth International Conference on Numerical Methods in Fluid Dynamics_ , Lecture Notes in Physics  Berlin, Heidelberg: Springer, 1982, pp. 507–512  DOI: 10.1007/3-540-11948-5˙66
  * [29] Xue-song Li and Chun-wei Gu  “An All-Speed Roe-type scheme and its asymptotic analysis of low Mach number behaviour”  In _Journal of Computational Physics_ 227.10, 2008, pp. 5144–5159  DOI: 10.1016/j.jcp.2008.01.037
  * [30] Xue-song Li and Chun-wei Gu  “Mechanism of Roe-type schemes for all-speed flows and its application”  In _Computers & Fluids_ 86, 2013, pp. 56–70  DOI: 10.1016/j.compfluid.2013.07.004
  * [31] Bo-Xi Lin, Chao Yan and Shu-Sheng Chen  “Density enhancement mechanism of upwind schemes for low Mach number flows”  In _Acta Mechanica Sinica_ 34.3, 2018, pp. 431–445  DOI: 10.1007/s10409-017-0737-9
  * [32] Meng-Sing Liou  “A sequel to AUSM, Part II: AUSM+-up for all speeds”  In _Journal of Computational Physics_ 214.1, 2006, pp. 137–170  DOI: 10.1016/j.jcp.2005.09.020
  * [33] Meng-Sing Liou  “A Sequel to AUSM: AUSM+”  In _Journal of Computational Physics_ 129.2, 1996, pp. 364–382  DOI: 10.1006/jcph.1996.0256
  * [34] Meng-Sing Liou  “The Root Cause of the Overheating Problem”  In _23rd AIAA Computational Fluid Dynamics Conference_ , AIAA AVIATION Forum  American Institute of AeronauticsAstronautics, 2017  DOI: 10.2514/6.2017-4107
  * [35] Meng-Sing Liou and J.. Edwards  “AUSM schemes and extensions for low Mach and multiphase flows”  In _30th Computational Fluid Dynamics Lecture Series_ von Karman Institute for Fluid Dynamics, 1999
  * [36] Meng-Sing Liou and Jack Edwards  “Numerical speed of sound and its application to schemes for all speeds” _eprint: https://arc.aiaa.org/doi/pdf/10.2514/6.1999-3268  In _14th Computational Fluid Dynamics Conference_ American Institute of AeronauticsAstronautics, 1999  DOI: 10.2514/6.1999-3268
  * [37] Meng-Sing Liou and Christopher J. Steffen  “A New Flux Splitting Scheme”  In _Journal of Computational Physics_ 107.1, 1993, pp. 23–39  DOI: 10.1006/jcph.1993.1122
  * [38] Yen Liu and Marcel Vinokur  “Upwind algorithms for general thermo-chemical nonequilibrium flows”  In _27th Aerospace Sciences Meeting_ , Aerospace Sciences Meetings  American Institute of AeronauticsAstronautics, 1989  DOI: 10.2514/6.1989-201
  * [39] Shingo Matsuyama  “Performance of all-speed AUSM-family schemes for DNS of low Mach number turbulent channel flow”  In _Computers & Fluids_ 91, 2014, pp. 130–143  DOI: 10.1016/j.compfluid.2013.12.010
  * [40] Yann Moguen et al.  “Pressure–velocity coupling allowing acoustic calculation in low Mach number flow”  In _Journal of Computational Physics_ 231.16, 2012, pp. 5522–5541  DOI: 10.1016/j.jcp.2012.05.001
  * [41] Yann Moguen, Erik Dick, Jan Vierendeels and Pascal Bruel  “Pressure–velocity coupling for unsteady low Mach number flow simulations: An improvement of the AUSM+-up scheme”  In _Journal of Computational and Applied Mathematics_ 246, Fifth International Conference on Advanced COmputational Methods in ENgineering (ACOMEN 2011), 2013, pp. 136–143  DOI: 10.1016/j.cam.2012.10.029
  * [42] B Müller  “Low Mach Number Asymptotics of the Navier-Stokes Equations and Numerical Implications”  In _30th Computational Fluid Dynamics Lecture Series_ von Karman Institute for Fluid Dynamics, 1999
  * [43] H. Paillere, C. Corre and J.. Garcia Cascales  “On the extension of the AUSM+ scheme to compressible two-fluid models”  In _Computers & Fluids_ 32.6, 2003, pp. 891–916  DOI: 10.1016/S0045-7930(02)00021-X
  * [44] Mark Potsdam, Venkateswaran Sankaran and Shishir Pandya  “Unsteady Low Mach Preconditioning with Application to Rotorcraft Flows” _eprint: https://arc.aiaa.org/doi/pdf/10.2514/6.2007-4473  In _18th AIAA Computational Fluid Dynamics Conference_ American Institute of AeronauticsAstronautics, 2007  DOI: 10.2514/6.2007-4473
  * [45] Feng Qu, Di Sun and Chao Yan  “A new flux splitting scheme for the Euler equations II: E-AUSMPWAS for all speeds”  In _Communications in Nonlinear Science and Numerical Simulation_ 57, 2018, pp. 58–79  DOI: 10.1016/j.cnsns.2017.09.002
  * [46] Felix Rieper  “A low-Mach number fix for Roe’s approximate Riemann solver”  In _Journal of Computational Physics_ 230.13, 2011, pp. 5263–5287  DOI: 10.1016/j.jcp.2011.03.025
  * [47] P. Roe  “Approximate Riemann solvers, parameter vectors, and difference schemes”  In _Journal of Computational Physics_ 43.2, 1981, pp. 357–372  DOI: 10.1016/0021-9991(81)90128-5
  * [48] Jai Sachdev, Ashvin Hosangadi and V. Sankaran  “Improved Flux Formulations for Unsteady Low Mach Number Flows”  In _42nd AIAA Fluid Dynamics Conference and Exhibit_ , Fluid Dynamics and Co-located Conferences  American Institute of AeronauticsAstronautics, 2012  DOI: 10.2514/6.2012-3067
  * [49] Lewis B. Schiff and Joseph L. Steger  “Numerical Simulation of Steady Supersonic Viscous Flow” Publisher: American Institute of Aeronautics and Astronautics  In _AIAA Journal_ 18.12, 1980, pp. 1421–1430  DOI: 10.2514/3.50902
  * [50] Eiji Shima  “On the improvement of the all-speed flux scheme for very low Mach number flows” _eprint: https://arc.aiaa.org/doi/pdf/10.2514/6.2013-2696  In _21st AIAA Computational Fluid Dynamics Conference_ American Institute of AeronauticsAstronautics, 2013  DOI: 10.2514/6.2013-2696
  * [51] Eiji Shima and Keiichi Kitamura  “CFD Method for Aero-acoustics using All-speed Numerical Flux and Preconditioned Implicit Time Integration”  In _20th AIAA Computational Fluid Dynamics Conference_ , Fluid Dynamics and Co-located Conferences  American Institute of AeronauticsAstronautics, 2011  DOI: 10.2514/6.2011-3045
  * [52] Eiji Shima and Keiichi Kitamura  “New approaches for computation of low Mach number flows”  In _Computers & Fluids_ 85, International Workshop on Future of CFD and Aerospace Sciences, 2013, pp. 143–152  DOI: 10.1016/j.compfluid.2012.11.017
  * [53] Eiji Shima and Keiichi Kitamura  “On New Simple Low-Dissipation Scheme of AUSM-Family for All Speeds”  In _47th AIAA Aerospace Sciences Meeting including The New Horizons Forum and Aerospace Exposition_ , Aerospace Sciences Meetings  American Institute of AeronauticsAstronautics, 2009  DOI: 10.2514/6.2009-136
  * [54] Eiji Shima and Keiichi Kitamura  “Parameter-Free Simple Low-Dissipation AUSM-Family Scheme for All Speeds” Publisher: American Institute of Aeronautics and Astronautics  In _AIAA Journal_ 49.8, 2011, pp. 1693–1709  DOI: 10.2514/1.J050905
  * [55] Joseph L Steger and R. Warming  “Flux vector splitting of the inviscid gasdynamic equations with application to finite-difference methods”  In _Journal of Computational Physics_ 40.2, 1981, pp. 263–293  DOI: 10.1016/0021-9991(81)90210-2
  * [56] Di Sun, Chao Yan, Feng Qu and Ruofan Du  “A robust flux splitting method with low dissipation for all-speed flows” _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/fld.4337  In _International Journal for Numerical Methods in Fluids_ 84.1, 2017, pp. 3–18  DOI: 10.1002/fld.4337
  * [57] B. Thornber et al.  “An improved reconstruction method for compressible flows with low Mach number features”  In _Journal of Computational Physics_ 227.10, 2008, pp. 4873–4894  DOI: 10.1016/j.jcp.2008.01.036
  * [58] B… Thornber and D. Drikakis  “Numerical dissipation of upwind schemes in low Mach flow” _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/fld.1628  In _International Journal for Numerical Methods in Fluids_ 56.8, 2008, pp. 1535–1541  DOI: https://doi.org/10.1002/fld.1628
  * [59] E.. Toro and M.. Vázquez-Cendón  “Flux splitting schemes for the Euler equations”  In _Computers & Fluids_ 70, 2012, pp. 1–12  DOI: 10.1016/j.compfluid.2012.08.023
  * [60] E. Turkel  “Preconditioning Techniques in Computational Fluid Dynamics” _eprint: https://doi.org/10.1146/annurev.fluid.31.1.385  In _Annual Review of Fluid Mechanics_ 31.1, 1999, pp. 385–416  DOI: 10.1146/annurev.fluid.31.1.385
  * [61] E. Turkel  “Review of preconditioning methods for fluid dynamics”  In _Applied Numerical Mathematics_ 12.1, SPECIAL ISSUE, 1993, pp. 257–284  DOI: 10.1016/0168-9274(93)90122-8
  * [62] Eli Turkel  “Preconditioned methods for solving the incompressible and low speed compressible equations”  In _Journal of Computational Physics_ 72.2, 1987, pp. 277–298  DOI: 10.1016/0021-9991(87)90084-2
  * [63] S Venkateswaran and Charles Merkle  “Dual time-stepping and preconditioning for unsteady computations”  In _33rd Aerospace Sciences Meeting and Exhibit_ , Aerospace Sciences Meetings  American Institute of AeronauticsAstronautics, 1995  DOI: 10.2514/6.1995-78
  * [64] S. Venkateswaran and C. Merkle  “Efficiency and accuracy issues in contemporary CFD algorithms” _eprint: https://arc.aiaa.org/doi/pdf/10.2514/6.2000-2251  In _Fluids 2000 Conference and Exhibit_ American Institute of AeronauticsAstronautics, 2000  DOI: 10.2514/6.2000-2251
  * [65] Sankaran Venkateswaran and Charles Merkle  “Artficial Dissipation Control for Viscous and Unsteady Computations” _eprint: https://arc.aiaa.org/doi/pdf/10.2514/6.2003-3695  In _16th AIAA Computational Fluid Dynamics Conference_ American Institute of AeronauticsAstronautics, 2003  DOI: 10.2514/6.2003-3695
  * [66] Yasuhiro Wada and Meng-Sing Liou  “A flux splitting scheme with high-resolution and robustness for discontinuities”  In _32nd Aerospace Sciences Meeting and Exhibit_ , Aerospace Sciences Meetings  American Institute of AeronauticsAstronautics, 1994  DOI: 10.2514/6.1994-83
  * [67] R. Warming and B. Hyett  “The modified equation approach to the stability and accuracy analysis of finite-difference methods”  In _Journal of Computational Physics_ 14.2, 1974, pp. 159–179  DOI: 10.1016/0021-9991(74)90011-4
  * [68] Jonathan M. Weiss and Wayne A. Smith  “Preconditioning applied to variable and constant density flows” Publisher: American Institute of Aeronautics and Astronautics  In _AIAA Journal_ 33.11, 1995, pp. 2050–2057  DOI: 10.2514/3.12946
  * [69] G-C. Zha and E. Bilgen  “Numerical solutions of Euler equations by using a new flux vector splitting scheme” _eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/fld.1650170203  In _International Journal for Numerical Methods in Fluids_ 17.2, 1993, pp. 115–144  DOI: 10.1002/fld.1650170203