6.2.2 Evaluation of the Flux Jacobian 通量雅可比的求值[cfd-6-2-2]

视底层空间离散化格式的类型而定ï¼解析求出方程(6.28)中的通量雅可比 $\partial\vec{R}/\partial\vec{W}$ 可能非常复杂ï¼甚至不可能。为使概念更清晰ï¼我们先对无黏流动导出通量雅可比ï¼然后再讨论向Navier-Stokes方程的推广。en

Depending on the type of the underlying spatial discretisation scheme, an analytical evaluation of the flux Jacobian $\partial\vec{R}/\partial\vec{W}$ in Eq. (6.28) may become very complex if not impossible. In order to make the concepts more clear, we shall derive the flux Jacobian for inviscid flows and then discuss the extension to the Navier-Stokes equations.

Central Scheme 中心格式

在中心空间离散化的情形下ï¼通量雅可比最容易构造。如图6.1的例子所示, $\partial\vec{R}/\partial\vec{W}$ 由对流通量雅可比组成(参见方程(6.34))ï¼它们可以解析导出(见A.9节)。人工黏性通常以简化形式纳入ï¼即不含非线性压力传感器(方程(4.55))。这一点将在6.2.3小节再讨论。en

The flux Jacobian is most easily formulated in the case of the central spatial discretisation. As we already saw for the example in Fig. 6.1, $\partial\vec{R}/\partial\vec{W}$ consists of the convective flux Jacobians (cf. Eq. (6.34)), which can be derived analytically (see Section A.9). Artificial viscosity is usually included in a simplified form, without the non-linear pressure sensor (Eq. (4.55)). We shall return to this point below in Subsection 6.2.3.

Flux-Vector Splitting Scheme 通量向量分裂格式

当以通量向量分裂格式(4.3.2小节)之一为基础来导出通量雅可比时ï¼其求值变得更加复杂。为说明起见ï¼考察Steger和Warming[30]的格式。已有的研究[31]表明ï¼对隐式算子的各种上风离散化,Steger-Warming分裂比例如Van Leer的通量向量分裂格式(方程(4.60))更可取。en

The evaluation of the flux Jacobian becomes more involved when one of the flux-vector splitting schemes (Subsection 4.3.2) is used as the basis for its derivation. Let us, for illustration, consider the scheme due to Steger and Warming [30]. Previous investigations [31] revealed that the Steger-Warming splitting is preferable over, e.g., the Van Leer's flux-vector splitting scheme (Eq. (4.60)) for various upwind discretisations of the implicit operator.

Steger-Warming通量向量分裂格式的基本思想是把对流通量分成正、负两部分ï¼即en

The basic idea of the Steger-Warming flux-vector splitting scheme is to divide the convective fluxes into a positive and a negative part, i.e.,

\vec{F}_c = \vec{F}_c^{+} + \vec{F}_c^{-} \tag{6.36}

其中通量定义为en

with the fluxes defined as

\vec{F}_c^{\pm} = \bar{A}_{SW}^{\pm}\vec{W} = \left(\bar{T}\bar{\Lambda}^{\pm}\bar{T}^{-1}\right)\vec{W}. \tag{6.37}

在方程(6.37)中, $\bar{A}_{SW}^{\pm}$ 表示正/负Steger-Warming通量分裂雅可比; $\bar{T}$ 表示右特征向量矩阵, $\bar{T}^{-1}$ 为左特征向量矩阵, $\bar{\Lambda}^{\pm}$ 代表正/负特征值构成的对角矩阵(参见A.11节)。特征值矩阵定义为[30]en

In Eq. (6.37), $\bar{A}_{SW}^{\pm}$ denotes the positive/negative Steger-Warming flux-splitting Jacobian. Furthermore, $\bar{T}$ represents the matrix of right eigenvectors, $\bar{T}^{-1}$ the matrix of left eigenvectors, and $\bar{\Lambda}^{\pm}$ stands for the diagonal matrix of positive/negative eigenvalues, respectively (cf. Section A.11). The eigenvalue matrices are defined as [30]

\bar{\Lambda}^{\pm} = \frac{1}{2}\left(\bar{\Lambda}_c \pm \left|\bar{\Lambda}_c\right|\right), \tag{6.38}

其中 $\bar{\Lambda}_c$ 由方程(A.84)给出。利用方程(6.36)定义的分裂ï¼对通量雅可比与方程(6.29)中更新量 $\Delta\vec{W}^n$ 的乘积得到en

where $\bar{\Lambda}_c$ is given by Eq. (A.84). Using the splitting defined in Eq. (6.36), we obtain for the product of the flux Jacobian with the update $\Delta\vec{W}^n$ in Eq. (6.29)

\frac{\partial\vec{R}_I}{\partial\vec{W}}\Delta\vec{W}^n = \sum_{m=1}^{N_F}\left[\frac{\partial\left(\vec{F}_c^{+}\Delta S\right)_m}{\partial\vec{W}_{L,m}}\Delta\vec{W}_{L,m}^n + \frac{\partial\left(\vec{F}_c^{-}\Delta S\right)_m}{\partial\vec{W}_{R,m}}\Delta\vec{W}_{R,m}^n\right]. \tag{6.39}

在上面的方程(6.39)中, $\Delta\vec{W}_{L,m}^n$ 和 $\Delta\vec{W}_{R,m}^n$ 分别表示面 $$m$$ 处左状态和右状态的更新量。在结构网格上ï¼左、右状态可以用MUSCL方法(方程(4.46))求值;在非结构网格上ï¼可以采用5.3.3小节讨论的重构方法。然而ï¼模板随精度提高而变宽ï¼导致系统矩阵的带宽增大。因此ï¼也为了降低数值复杂性ï¼方程(6.39)中通常只采用一阶精度近似。作为一种折中ï¼可以用更高精度重构左、右状态ï¼但在求导数时保留一阶格式的模板[32]。en

In the above Eq. (6.39), $\Delta\vec{W}_{L,m}^n$ and $\Delta\vec{W}_{R,m}^n$ denote the updates of the left and right state at the face $$m$$ , respectively. On structured grids, the left and right state can be evaluated by the MUSCL approach (Eq. (4.46)). On unstructured grids, the reconstruction methods discussed in Subsection 5.3.3 can be applied. However, the stencil becomes wider with increasing accuracy, which leads to larger bandwidth of the system matrix. Therefore, and in order to reduce the numerical complexity, only first-order accurate approximation is usually employed in Eq. (6.39). As a compromise, we could reconstruct the left and right state with higher accuracy but retain the stencil of the first-order scheme for the evaluation of the derivatives [32].

为继续讨论方程(6.39)中导数 $\partial\vec{F}_c^{\pm}/\partial\vec{W}$ 的求值ï¼以面 $$m$$ 处的正通量为例来考察en

To proceed with the discussion on the evaluation of the derivatives $\partial\vec{F}_c^{\pm}/\partial\vec{W}$ in Eq. (6.39), let us consider, e.g., the positive flux at face $$m$$

\frac{\partial\left(\vec{F}_c^{+}\Delta S\right)_m}{\partial\vec{W}_{L,m}} = \frac{\partial}{\partial\vec{W}_{L,m}}\left[\left(\bar{A}_{SW}^{+}\vec{W}\right)_{L,m}\Delta S_m\right] \tag{5}

利用方程(6.37)和(6.38)ï¼它变为en

which becomes with Eqs. (6.37), (6.38)

\begin{aligned}\frac{\partial\left(\vec{F}_c^{+}\Delta S\right)_m}{\partial\vec{W}_{L,m}} = \frac{\Delta S_m}{2}\left[\left(\bar{A}_c\right)_{L,m} + \left|\left(\bar{A}_c\right)_{L,m}\right|\right] \\+ \frac{\Delta S_m}{2}\left[\frac{\partial\left(\bar{A}_c\right)_{L,m}}{\partial\vec{W}_{L,m}} + \frac{\partial\left|\left(\bar{A}_c\right)_{L,m}\right|}{\partial\vec{W}_{L,m}}\right]\vec{W}_{L,m}^n.\end{aligned} \tag{6.40}

对负通量也可以得到与方程(6.40)类似的表达式。可以看到ï¼方程(6.40)的第一项由对流通量雅可比组成(见A.9节)ï¼因而没有困难;但第二项涉及矩阵元素的导数。虽然可以通过手工推导或使用符号代数软件解析地得到这些导数ï¼但这会产生庞大而计算效率低下的代码[33]。另一种可能的做法是假定矩阵 $\bar{A}_c$ 局部为常值ï¼从而可以忽略方程(6.40)中的第二项。然而ï¼视隐式格式的类型而定ï¼这可能会严重限制CFL数[34]。en

An expression similar to Eq. (6.40) can also be found for the negative flux. As we can see, the first term in Eq. (6.40) consists of convective flux Jacobians (see Section A.9) and thus presents no difficulty. However, the second term involves derivatives of matrix elements. Although it is possible to obtain the derivatives analytically either by hand calculation or by using a symbolic algebra package, this will produce a large, computationally inefficient code [33]. Alternatively, it is possible to assume the matrix $\bar{A}_c$ is locally constant so that the second term in Eq. (6.40) can be neglected. However, depending on the type of the implicit scheme, this may severely restrict the CFL number [34].

计算方程(6.39)中导数 $\partial\vec{F}_c^{\pm}/\partial\vec{W}$ 的其他可行方法还有源代码的自动微分(例如使用ADIFOR[35])或有限差分法(参见例如[26]、[33])。这样ï¼向量 $\vec{F}$ 的第 $$i$$ 个分量对因变量 $\vec{X}$ 第 $$j$$ 个分量的导数可以近似为en

Other approaches that we could use to compute the derivatives $\partial\vec{F}_c^{\pm}/\partial\vec{W}$ in Eq. (6.39) would be the automatic differentiation of the source code (e.g., using ADIFOR [35]) or the finite-difference method (see, e.g., [26], [33]). Herewith, the derivative of the $$i$$ -th component of a vector $\vec{F}$ with respect to the $$j$$ -th component of a dependent variable $\vec{X}$ can be approximated as

\frac{\partial f_i}{\partial x_j} \approx \frac{f_i\left(\vec{X} + h_j\vec{e}^{\,j}\right) - f_i\left(\vec{X}\right)}{h_j}, \tag{6.41}

其中 $\vec{e}^{\,j}$ 表示第 $$j$$ 个标准基向量。Dennis和Schnabel[36]建议步长 $$h_j$$ 取如下形式en

where $\vec{e}^{\,j}$ denotes the $$j$$ -th standard basis vector. Dennis and Schnabel [36] suggested a stepsize $$h_j$$ of the form

h_j = \sqrt{\epsilon}\,\max\left\{\left|x_j\right|, \text{typ}\,x_j\right\}\,\text{sign}(x_j) \tag{6.42}

其中 $\epsilon$ 为机器精度, $\text{typ}\,x_j$ 为 $$x_j$$ 的典型大小。关于雅可比矩阵的高效数值求值ï¼还可参阅[37]和[38]中的提示。en

with $\epsilon$ being the machine accuracy and $\text{typ}\,x_j$ a typical size of $$x_j$$ . The reader is also referred to [37] and [38] for hints on efficient numerical evaluation of Jacobian matrices.

Flux-Difference Splitting Scheme 通量差分分裂格式

对Roe的通量差分分裂格式(4.3.3小节ï¼方程(4.91))ï¼通量雅可比与方程(6.29)中更新量的乘积可以写为en

In the case of the flux-difference splitting scheme due to Roe (Subsection 4.3.3, Eq. (4.91)), we can write the product of the flux Jacobian with the update in Eq. (6.29) as

\begin{aligned}\frac{\partial\vec{R}_I}{\partial\vec{W}}\Delta\vec{W}^n = \sum_{m=1}^{N_F}\frac{\Delta S_m}{2}\Big\{ &\left(\bar{A}_c\right)_{L,m}\Delta\vec{W}_{L,m}^n + \left(\bar{A}_c\right)_{R,m}\Delta\vec{W}_{R,m}^n \\&- \frac{\partial}{\partial\vec{W}_{L,m}}\left[\left|\bar{A}_{Roe}\right|_m\left(\vec{W}_{R,m}^n - \vec{W}_{L,m}^n\right)\right]\Delta\vec{W}_{L,m}^n \\&- \frac{\partial}{\partial\vec{W}_{R,m}}\left[\left|\bar{A}_{Roe}\right|_m\left(\vec{W}_{R,m}^n - \vec{W}_{L,m}^n\right)\right]\Delta\vec{W}_{R,m}^n\Big\}.\end{aligned} \tag{6.43}

与通量向量分裂类似ï¼表达式(6.43)既包含对流通量雅可比ï¼也包含Roe矩阵 $\bar{A}_{Roe}$ 的导数。这些导数已在[34]中给出。由于相应的公式非常复杂(另见文献[33])ï¼更好的做法是像上文讨论的那样ï¼用数值方法求值 $\partial\vec{F}_c/\partial\vec{W}$ 项。en

Similar to flux-vector splitting, the expression (6.43) contains convective flux Jacobians as well as derivatives of the Roe matrix $\bar{A}_{Roe}$ . The derivatives were presented in [34]. Since the corresponding formulae are very complex (see also Ref. [33]), it is a better idea to evaluate the term $\partial\vec{F}_c/\partial\vec{W}$ numerically, as discussed above.

不过ï¼我们也可以假设Roe矩阵局部为常数[39]ï¼从而简化式(6.43):en

However, we can also simplify Eq. (6.43) by assuming locally constant Roe matrices [39]

\frac{\partial\vec{R}_I}{\partial\vec{W}}\Delta\vec{W}^n \approx \sum_{m=1}^{N_F}\frac{\Delta S_m}{2}\Big\{ \left(\bar{A}_c\right)_{L,m}\Delta\vec{W}_{L,m}^n + \left(\bar{A}_c\right)_{R,m}\Delta\vec{W}_{R,m}^n - \left|\bar{A}_{Roe}\right|_m\left(\Delta\vec{W}_{R,m}^n - \Delta\vec{W}_{L,m}^n\right)\Big\}. \tag{6.44}

与Steger-Warming通量向量分裂格式不同ï¼上述近似线性化(6.44)对隐式格式性能的降低非常有限[34]。en

In contrast to the Steger-Warming flux-vector splitting scheme, the above approximate linearisation (6.44) degrades the performance of the implicit scheme only slightly [34].

Viscous Flows 黏性流动

对于Navier-Stokes方程ï¼还必须在隐式算子中计入黏性通量。导数 $\partial\vec{F}_v/\partial\vec{W}$ ï¼即式(6.30)中的黏性通量雅可比ï¼一般不容易直接得到。额外的复杂性来自黏性通量向量本身包含流动变量的导数这一事实。因此ï¼要么用有限差分(式(6.41))求黏性通量雅可比ï¼要么采用简化的形式。en

For the Navier-Stokes equations, we have to account also for the viscous fluxes in the implicit operator. The derivative $\partial\vec{F}_v/\partial\vec{W}$ , i.e., the viscous flux Jacobian in Eq. (6.30) is in general not straightforward to obtain. Additional complexity arises due to the fact that the viscous flux vector contains derivatives of flow variables. For this reason, we have either to evaluate the viscous flux Jacobian by finite differences (Eq. (6.41)), or we have to use a simplified formulation.

在Navier-Stokes方程的TSL近似下(参见2.4.3小节与附录A.6节)ï¼通过假设动力黏度与热传导系数局部为常数ï¼可以解析地求出黏性通量雅可比。于是ï¼根据附录A.10中的讨论ï¼式(6.29)中与黏性通量相关的项变为en

In the case of the TSL approximation of the Navier-Stokes equations (cf. Subsection 2.4.3 and Section A.6), it is possible to find the viscous flux Jacobian analytically by assuming locally constant dynamic viscosity and thermal conductivity coefficients. Then, according to the discussion in Appendix A.10, the term related to the viscous fluxes in Eq. (6.29) becomes

\frac{\partial(\vec{F}_v\Delta S)_m}{\partial\vec{W}}\Delta\vec{W}^n \approx \left[\left(\bar{A}^{*}_v\right)_{R,m}\Delta\vec{W}_{R,m}^n - \left(\bar{A}^{*}_v\right)_{L,m}\Delta\vec{W}_{L,m}^n\right]\Delta S_m. \tag{6.45}

在上式(6.45)中, $\bar{A}^{*}_v$ 表示由式(A.71)或式(A.75)给出的黏性通量雅可比ï¼但不含空间算子 $\partial_y(\cdot)$ (参见式(A.74)与(A.79))。除动力黏度取算术平均外ï¼这些雅可比对其余所有变量要么用左状态、要么用右状态求值。应当指出ï¼假如左右状态以一阶精度计算ï¼式(6.45)在隐式算子中给出的是二阶中心差分近似。en

In above Eq. (6.45), $\bar{A}^{*}_v$ stands for the viscous flux Jacobian given by Eq. (A.71) or Eq. (A.75) but without the spatial operators $\partial_y(\cdot)$ (cf. Eq. (A.74) and (A.79)). The Jacobians are evaluated using either the left or the right state for all variables except for the dynamic viscosity, which is determined from an arithmetic average. It should be mentioned that Eq. (6.45) leads to a second-order central difference approximation in the implicit operator, supposed the left and right state are computed with first-order accuracy.