4.3.1 Central Scheme with Artificial Dissipation 含人工耗散的中心格式[cfd-4-3-1]

The central scheme with artificial dissipation is very simple compared to other discretisation methods. It is easy to implement with either the cell-centred scheme or with both cell-vertex schemes. For these reasons the scheme became very wide-spread.

The basic idea of the central scheme is to compute the convective fluxes at a face of the control volume from the arithmetic average of the conservative variables on both sides of the face. Since this would allow for odd-even decoupling of the solution (generation of two independent solutions of the discretised equations) and for overshoots at shocks, artificial dissipation (which is similar to the viscous fluxes) has to be added for stability. The scheme was first implemented for the Euler equations by Jameson et al. [6]. Because of the names of the authors, it is also abbreviated as the JST scheme.

The central scheme is generally less accurate in the resolution of discontinuities and boundary layers than, let say, the upwind schemes. However, it is computationally considerably cheaper. Therefore, attempts were made to improve the accuracy of the scheme, while still keeping the numerical effort low. For example, improvements were devised to reduce the amount of the artificial dissipation in boundary layers [31], [32], or to enhance the shock resolution [33]-[35]. Another idea, followed in [35], is to utilise the Jacobian matrix of the convective fluxes, in order to scale the dissipation independently for each conservation equation. This successful approach is known as the matrix dissipation scheme. It can be viewed as a compromise between the original scalar scheme and the upwind schemes. The basic scheme also employs a single, pressure-based sensor to switch from second- to first-order accuracy at discontinuities to prevent non-physical oscillations of the flow variables. In [36], Jameson developed a concept called the SLIP (Symmetric Limited Positive) scheme, where a limiter is applied separately for each conservation equation. A brief description of the previous approaches and comparisons for inviscid and viscous 2-D flows was presented in [37].

The convective fluxes (Eq. (2.21)) through a face of the control volume are approximated using the average of variables, according to the Equations (4.16), (4.21), (4.36), or (4.41), respectively. Artificial dissipation is then added to the central fluxes for stability [6], [38]. Thus, the total convective flux at face \((I+1/2,J,K)\) reads

where the flow variables are averaged as (see also Fig. 4.8)

In the case of the cell-vertex scheme with dual control volumes, node indices \((i,j,k)\) would be used instead. For simplicity, \((I+1/2,J,K)\) will be abbreviated as \((I+1/2)\) hereafter. The artificial dissipation flux consists of a blend of adaptive second- and fourth-order differences, which result from the sum of first- and third-order difference operators

From Eq. (4.50) we can see that the scheme possesses a compact 9-point stencil in two dimensions and a 13-point stencil in three dimensions. The dissipation is scaled by the sum of the spectral radii of the convective flux Jacobians in all coordinate directions

The spectral radius at the cell face \((I+1/2)\), e.g., in I-direction (represented by the superscript \(I\)), results from the average

It is evaluated using the formula

where \(V\) stands for the contravariant velocity (2.22) and \(c\) for the speed of sound, respectively. A pressure-based sensor is used to switch off the fourth-order differences at shocks, where they would lead to strong oscillation of the solution. The sensor also switches off the second-order differences in smooth parts of the flow field, in order to reduce the dissipation to the lowest possible level. Herewith, the coefficients \(\epsilon^{(2)}\) and \(\epsilon^{(4)}\) in Eq. (4.50) are defined as

with the pressure sensor given by

Typical values of the parameters are \(k^{(2)} = 1/2\) and \(1/128 \le k^{(4)} \le 1/64\). In order to reduce the amount of artificial dissipation across a viscous shear layer, we can re-define the scaling factors in Eq. (4.50) as follows [31], [32]

These are then employed instead of Eq. (4.51). The directionally dependent coefficients \(\phi\) are given by the relations

The parameter \(\sigma\) is usually set equal to 1/2 or 2/3. This formulation decreases the scaling of the dissipation terms in the direction along the shorter side of a control volume, whose longer side is aligned with the flow.

In order to improve the accuracy by reducing the numerical dissipation, the preceding JST scheme can be modified to become more like an upwind scheme. The idea is to use a matrix - the convective flux Jacobian - instead of the scalar value \(\hat{\Lambda}^{S}\) to scale the dissipation terms [35]. In this way, each equation is scaled properly by the corresponding eigenvalue. Hence, the Eq. (4.50) becomes

The scaling matrix corresponds to the convective flux Jacobian \((\bar{A}_c = \partial\vec{F}_c/\partial\vec{W})\) diagonalised with absolute values of the eigenvalues

The matrices of right \((\bar{T})\) and left \((\bar{T}^{-1})\) eigenvectors as well as the diagonal matrix of the eigenvalues \(\bar{\Lambda}_c\) can be found in the Appendix A.11. The eigenvalues must be limited at stagnation points and sonic lines to prevent the dissipation from becoming zero. An efficient way of computing the product of \(|\bar{A}_c|\) with \(\vec{W}\) is provided in [35]. It should be noted that by setting \(\epsilon^{(2)} = 1/2\) and \(\epsilon^{(4)} = 0\), we obtain a first-order accurate, fully upwind scheme.

As it was already stated before, the idea here was to develop a scheme which accuracy is close to that of upwind schemes, but which is still computationally only slightly more expensive (about 15-20%) than the scalar dissipation approach. Results of comparisons with flux-vector splitting schemes (CUSP and AUSM) were recently reported in [37].