8.3 远场

8.3 远场

对翼型、机翼、汽车及其他构型周围的外部流动进行数值模拟，必须采用有限区域，因此需要人工远场边界条件。远场边界条件的数值实现必须满足两个基本要求：首先，与无限区域相比，截断区域不应对流动解产生明显影响；其次，任何向外传播的扰动均不得反射回流场 $ [9] $。亚声速和跨声速问题具有椭圆型特性，对远场边界条件尤其敏感。不恰当的实现会显著减慢向定常解的收敛，并可能降低解的精度。人们已发展出多种能在人工边界吸收外传波动的方法 $ [10-15] $。不同无反射边界条件的综述见文献 $ [16] $。

下面两节讨论 Whitfield 和 Janus $ [12] $ 所述的特征变量概念，并介绍针对升力物体的远场边界条件扩展。

8.3.1 特征变量的概念

根据对流通量 Jacobian 特征值的正负（A.11 节，式（A.84）或（A.88）），信息沿特征线传出或传入计算域。例如，亚声速入流时，有四条特征线进入计算域（三维情况），一条离开计算域（式（A.88）中的 $ \Lambda_{5} $）。亚声速出流时则相反。根据 Kreiss 的一维理论 [17]，应从外部施加于边界的条件数必须等于进入计算域的特征线数，其余条件应由区域内部的解确定。

Whitfield 和 Janus 的方法 $ [12] $ 基于沿边界法向的 Euler 方程特征形式（2.44）。该方法在结构网格、非结构网格以及多种流动问题中均表现很好，不仅适用于远场边界，也适用于无黏固体壁面（8.2.1 节）。

图 8.6 示出了远场边界上的两种基本情况：流动可以进入或离开计算域。因此，根据局部马赫数，需要处理四种不同的远场边界条件：

• 超声速入流；

• 超声速出流；

• 亚声速入流；

• 亚声速出流。

超声速入流

对于超声速入流，所有特征值同号。由于流动进入物理区域，边界上的守恒变量（图 8.6 的点 b）仅由自由来流值确定，因此

$$ \vec{W}_{b}=\vec{W}_{a}.\tag{8.20} $$

值 $ \vec{W}_{a} $ 依据给定的马赫数 $ M_{\infty} $ 和两个流动角度（迎角与侧滑角）确定。

超声速出流

这种情况下，所有特征值也同号，但流动离开物理区域，因此边界上的全部守恒变量必须由区域内部的解确定。可以直接设置

$$ \vec{W}_{b}=\vec{W}_{d},\tag{8.21} $$

这相当于常值外推。也可以采用线性外推，但可能降低求解格式的稳定性。

亚声速入流

此时，四条特征线进入物理区域，一条离开。因此，根据自由来流值指定四个特征变量，另一个特征变量则从物理区域内部外推，得到以下边界条件 $ [12] $：

$$ \begin{aligned}p_b&=\frac{1}{2}\big\{p_{a}+p_{d}-\rho_{0}c_{0}\left[n_{x}(u_{a}-u_{d})+n_{y}( v_{a}- v_{d})+n_{z}(w_{a}-w_{d})\right]\big\}\\ \rho_{b}&=\rho_{a}+(p_{b}-p_{a})/c_{0}^{2}\\ u_{b}&=u_{a}-n_{x}(p_{a}-p_{b})/(\rho_{0}c_{0})\\ v_{b}&= v_{a}-n_{y}(p_{a}-p_{b})/(\rho_{0}c_{0})\\ w_{b}&=w_{a}-n_{z}(p_{a}-p_{b})/(\rho_{0}c_{0}),\end{aligned}\tag{8.22} $$

其中，$ \rho_{0} $ 和 $ c_{0} $ 表示参考状态，通常取为内部点（图 8.6 的点 d）的状态。点 a 的值由自由来流状态确定。

亚声速出流

对于亚声速出流，四个流动变量（密度及三个速度分量）必须从物理区域内部外推，剩余的第五个变量（压力）必须从外部指定。远场边界的原始变量按文献 $ [12] $ 的关系计算：

$$ \begin{aligned}&p_{b}=p_{a}\\&\rho_{b}=\rho_{d}+(p_{b}-p_{d})/c_{0}^{2}\\&u_{b}=u_{d}+n_{x}(p_{d}-p_{b})/(\rho_{0}c_{0})\\& v_{b}= v_{d}+n_{y}(p_{d}-p_{b})/(\rho_{0}c_{0})\\&w_{b}=w_{d}+n_{z}(p_{d}-p_{b})/(\rho_{0}c_{0})\\ \end{aligned}\tag{8.23} $$

其中，$ p_{a} $ 为指定的静压。

虚拟单元或节点中的物理量可由状态 b 和 d 作线性外推得到。

8.3.2 针对升力物体的修正

上述特征远场边界条件假设环量为零，这对于亚声速或跨声速流动中的升力物体并不正确。因此，远场边界必须位于距物体很远的位置，否则流动解会不准确。若在自由来流中计及单个涡的影响（三维时为马蹄涡），则可显著缩短远场距离，约减少一个数量级。假设该涡位于升力物体的中心，其强度与物体产生的升力成正比。下面介绍二维及三维涡修正的实现。

二维涡修正

这里介绍 Usab 和 Murman $ [18] $ 提出的方法。假设流动可压缩，修正后的自由来流速度分量为

$$ \begin{gathered}u_{\infty}^{*}=u_{\infty}+\left(\frac{\Gamma\sqrt{1-M_{\infty}^{2}}}{2\pi d}\right)\frac{1}{1-M_{\infty}^{2}\sin^{2}(\theta-\alpha)}\sin\theta,\\v_{\infty}^{*}= v_{\infty}-\left(\frac{\Gamma\sqrt{1-M_{\infty}^{2}}}{2\pi d}\right)\frac{1}{1-M_{\infty}^{2}\sin^{2}(\theta-\alpha)}\cos\theta\end{gathered}\tag{8.24} $$

其中，$ \Gamma $ 为环量，$ (d,\theta) $ 为远场点的极坐标，$ \alpha $ 为迎角，$ M_{\infty} $ 表示自由来流马赫数。环量由下式得到：

$$ \Gamma=\frac{1}{2}\|\vec{ v}_{\infty}\|_{2}a C_{\mathrm{L}}\tag{8.25} $$

这里采用 Kutta–Joukowsky 定理。式（8.25）中，a 表示翼型弦长，$ C_{L} $ 为通过表面压力积分求得的升力系数。式（8.24）的极坐标计算为

$$ \begin{gathered}d=\sqrt{(x-x_{\mathrm{ref}})^{2}+( y- y_{\mathrm{ref}})^{2}}\\\theta=\tan\left(\frac{ y- y_{\mathrm{ref}}}{x-x_{\mathrm{ref}}}\right),\end{gathered}\tag{8.26} $$

校注：原书极角公式写为正切函数 $\tan$。由笛卡尔坐标求极角应使用反正切，实际计算宜采用 $\operatorname{atan2}(y-y_{\mathrm{ref}},x-x_{\mathrm{ref}})$ 以区分象限。译文保留原式。

其中，$ x_{ref} $、$ y_{\mathrm{ref}} $ 为参考点的坐标，即涡的位置，例如四分之一弦长处。

修正后的自由来流压力 $ p_{\infty}^{*} $ 为

$$ p_{\infty}^{*}=\left[p_{\infty}^{(\gamma-1)/\gamma}+\left(\frac{\gamma-1}{\gamma}\right)\frac{\rho_{\infty}\left(\|\vec{ v}_{\infty}\|_{2}^{2}-\|\vec{ v}_{\infty}^{*}\|_{2}^{2}\right)}{2p_{\infty}^{1/\gamma}}\right]^{\gamma/(\gamma-1)}\tag{8.27} $$

其中，$ \|\vec{ v}_{\infty}^{*}\|_{2}^{2}=(u_{\infty}^{*})^{2}+( v_{\infty}^{*})^{2} $。修正后的自由来流密度由状态方程得到：

$$ \rho_{\infty}^{*}=\rho_{\infty}\left(\frac{p_{\infty}^{*}}{p_{\infty}}\right)^{1/\gamma}.\tag{8.28} $$

将修正后的量 $ u_{\infty}^{*} $、$ v_{\infty}^{*} $、$ p_{\infty}^{*} $ 和 $ \rho_{\infty}^{*} $ 代入式（8.22）或（8.23），替代 $ u_{a} $、$ v_{a} $、$ p_{a} $ 和 $ \rho_{a} $。

上述涡修正（式（8.24）—（8.28））严格来说仅适用于亚声速流动，但修正自由来流条件对跨声速流动也有帮助。图 8.7 研究了升力系数对远场边界距离的依赖性，远场半径分别取 5、20、50 和 99 倍弦长。可以看出，不采用涡修正时，模拟结果对远场距离高度敏感；采用涡修正时，远场距离缩短至约 20 倍弦长仍能保持足够精度，从而显著减少网格单元或点数。文献 [19] 表明，采用涡修正的高阶项后，远场边界可以仅距物体约 5 倍弦长，而没有明显的精度损失。

三维涡修正

机翼对远场边界的影响可以用马蹄涡近似。对于可压缩流动，修正后的自由来流速度分量由以下关系得到 $ [20, 21] $：

$$ \begin{aligned}u_{\infty}^{*}&=u_{\infty}+\frac{\Gamma\beta^{2}}{2\pi}\mathcal{A},\\ v_{\infty}^{*}&= v_{\infty}-\frac{\Gamma}{2\pi}\left[\frac{z+l}{(z+l)^{2}+ y^{2}}\mathcal{B}-\frac{z-l}{(z-l)^{2}+ y^{2}}\mathcal{C}+\frac{x\beta^{2}}{x^{2}+ y^{2}\beta^{2}}\mathcal{A}\right],\\w_{\infty}^{*}&=w_{\infty}+\frac{\Gamma}{2\pi}\left[\frac{ y}{(z+l)^{2}+ y^{2}}\mathcal{B}-\frac{ y}{(z-l)^{2}+ y^{2}}\mathcal{C}\right],\end{aligned}\tag{8.29} $$

其中，$ \Gamma $ 表示环量，$ (x, y, z) $ 是远场点的笛卡尔坐标，l 为半翼展。此外，式（8.29）假设流动沿 x 轴正方向，机翼沿 z 轴布置。式（8.29）中的 A、B、C 为

$$ \begin{array}{l}\mathcal{A}=\frac{z+l}{\sqrt{\psi_{+}}}-\frac{z-l}{\sqrt{\psi_{-}}},\\\mathcal{B}=1+\frac{x}{\sqrt{\psi_{+}}},\\\mathcal{C}=1+\frac{x}{\sqrt{\psi_{-}}}.\end{array}\tag{8.30} $$

各简记量为

$$ \begin{aligned}\psi_{+}&=x^{2}+\beta^{2}(z+l)^{2}+ y^{2}\beta^{2},\\\psi_{-}&=x^{2}+\beta^{2}(z-l)^{2}+ y^{2}\beta^{2},\\\beta&=\sqrt{1-M_{\infty}^{2}}\end{aligned}\tag{8.31} $$

其中，$ M_{\infty} $ 为自由来流马赫数。环量 $ \Gamma $ 利用式（8.25）计算，此时 a 表示平均弦长。修正压力 $ p_{\infty}^{*} $ 和密度 $ \rho_{\infty}^{*} $ 分别由式（8.27）、（8.28）得到。式（8.22）或（8.23）中的量 $ u_{a} $、$ v_{a} $、$ w_{a} $、$ p_{a} $ 和 $ \rho_{a} $ 分别替换为其修正值 $ u_{\infty}^{*} $、$ v_{\infty}^{*} $、$ w_{\infty}^{*} $、$ p_{\infty}^{*} $ 和 $ \rho_{\infty}^{*} $。

式（8.29）中修正速度分量 $ v_{\infty}^{*} $ 和 $ w_{\infty}^{*} $ 的表达式，在涡线穿过出口边界的位置趋于无穷大。这些点为

$$ \begin{aligned}z&=+l,\quad y=0,\\z&=-l,\quad y=0,\end{aligned} $$

以及 $ x = x_{farf} $。为避免数值奇点，文献 $ ^{[21]} $ 建议将

$$ \begin{array}{l}(z+l)^{2}\ + y^{2},\\ (z-l)^{2}\ + y^{2}\end{array} $$

在式（8.29）中的值限制为四分之一翼展，即 l/2。这一简单措施减小了涡线 z=l 和 z=-l 周围 l/2 距离以内的速度 $ v_{\infty} $ 和 $ w_{\infty} $ 修正量。

文献 $ [21] $ 的数值结果表明，采用式（8.29）的涡修正可以降低升力和阻力系数对远场距离的敏感性。研究发现，远场距离为 $ 7 \cdot l $ 时即可获得准确结果。