一种基于流场密度阶跃的高精度间断迦辽金人工粘性激波捕捉方法与流程

本文涉及计算流体力学技术流域，具体涉及一种高精度间断迦辽金人工粘性激波捕捉方法。

背景技术：

由于高精度间断迦辽金(DG)方法在数值耗散和色散方面的优良特性，受到越来越多CFD学者的高度重视并得到了蓬勃的发展。根据Godunov原理，高于一阶精度的线性格式都不会有单调解。高精度DG方法在激波间断附近会产生Gibbs现象，引起非物理解的产生，导致计算中断。如何有效地进行激波捕捉，保证格式的单调性和鲁棒性，已经成为DG方法的重要发展方向之一。DG方法中抑制间断附近数值震荡的方法一般有限制器和重构方法。

限制器的思想最早源于总变差减少(不增)格式(TVD，Total Variation Diminishing)，但是在限制过程中，往往只用到了临近单元的平均值信息，高精度情况下的表现不理想。在间断附近会自动降为一阶精度，对间断的分辨率不够高，激波的抹平效应明显，无法有效捕捉。采用重构进行激波捕捉，在非结构网格上实现时，无法保持DG格式的紧致特性，在工程实用化方面仍有困难。

限制器和重构方法的核心思想是对流场当前步的解进行后处理操作，但是在进行操作之前，都要对流场中的激波区域进行探测，称之为激波探测，然后在探测到激波存在的区域进行限制器或者重构的操作。这种操作方式不仅形式复杂，而且激波探测的准确性对后续的处理影响很大，而以往的激波探测过程都存在探测不准确的情况。

技术实现要素：

本发明的目的是提出一种基于流场密度阶跃的高精度间断迦辽金人工粘性激波捕捉方法，该方法针对以往的高精度DG格式中激波捕捉需要提前进行激波探测、过程复杂等问题，基于流场中物理解的特性，在保证量纲一致的情况下，在保证鲁棒性和稳定性的情况下，省去了激波探测过程，简化了激波捕捉所需的构造变量。

为了实现上述目的，本发明采用如下技术方案：

一种基于流场密度阶跃的高精度间断迦辽金人工粘性激波捕捉方法，包括以下步骤：

步骤一：采用Euler方程作为控制方程，建立DG高精度框架，将整个计算区域进行网格剖分，对每个网格单元的高斯积分点进行计算和存储；

步骤二：选择流场中的密度作为阶跃变量，同时选择网格单元的尺度和速度作为变量，与阶跃变量进行组合，构造交界面的人工粘性系数；

步骤三：将交界面处的人工粘性系数划归到单元内，然后带入到控制方程中，通过迭代计算，得到仿真的气动结果和流场。

本发明建立DG高精度框架，控制方程采用Euler方程，计算区域采用非结构网格进行剖分。方程的对流项采用高阶HLLC格式离散求解，同时在主控方程的基础上添加人工粘性项，在保证鲁棒性、计算精度的基础上捕捉激波。人工粘性项的添加基于网格单元面上的阶跃，无需激波探测技术，采用归一化的方式保证量纲一致并将阶跃量划归到网格单元内。

与现有的手段比较，本发明的特点在于：

1、选择密度作为阶跃变量来构造人工粘性，不需要额外的激波探测等过程，构造算法所需的变量少，计算过程得到简化。

2、相比之前的激波捕捉方法，该方法的的激波捕捉效率得到提高，有效捕捉激波所需的迭代步数明显减少，所用的CPU时间明显减少。

附图说明

本文将通过例子并参照附图的方式说明，其中：

图1为本发明方法的流程图；

图2a、图2b为本发明方法在高精度情况下对NACA0012翼型算例的激波捕捉问题的测试。从图中可以看出，

具体实施方式

本说明书中公开的所有特征，或公开的所有方法或过程中的步骤，除了互相排斥的特征和/或步骤以外，均可以以任何方式组合。

本发明的一种基于流场密度阶跃的高精度间断迦辽金人工粘性激波捕捉方法，如图1所示包括三个部分，其中：

第一部分：采用Euler方程，建立DG高精度框架。包括以下步骤：

步骤101、构建微分形式下的Euler方程

其中U代表流场中的守恒量,U是一个矢量，U＝(ρ,ρu,ρv,ρw,ρE)^T，代表守恒量相对于时间t的偏导数，F^c代表守恒通量，代表守恒通量的散度。

步骤102、对流场区域网格剖分，并计算得到网格单元内的Gauss积分点。

步骤103、选择相应的基函数φi，基函数的选取没有限制，可以采用Taylor基。对公式(1)在单元Ω内进行积分并乘以基函数，利用格林高斯公式，得到弱形式下的Euler方程。

其中为单元边界面，为单元边界面的外法线方向，代表基函数的梯度。

步骤104、将守恒量U用基函数的线性组合代替，并带入到方程(2)中，最终得到高精度DG框架下的求解方程。

其中，M＝∫Ωφiφj代表质量矩阵，u代表守恒变量用基函数线性组合时采用的系数又称为自由度，代表自由度对时间的导数。

第二部分：选择流场中的密度作为阶跃变量，与单元内的速度和单元尺度组合，构造人工粘性。包括以下步骤：

步骤201、选择拉普拉斯人工粘性模型，带入公式(1)并重复步骤103和104得到如下形式

其中，ε为人工粘性系数，决定单元内人工粘性的数值，代表自由度的梯度，代表自由度对时间的导数。

步骤202、将单元交界面处密度的阶跃定义为[ρ]＝|ρ⁺-ρ^-|，交界面处密度的平均定义为ρ⁺和ρ^-分别代表单元交界面左右两侧的值，单元交界面处的人工粘性系数εjump形式如下

其中，V代表单元内的速度，c代表当地声速，代表当地网格单元的特征尺度。Cε为经验参数，一般取0.1。

第三部分：将人工粘性系数带入到控制方程中，迭代计算，得到结果。具体步骤如下

步骤301、步骤202中得到的人工粘性系数分配到左右两侧的单元中，为了简化处理，认为两侧单元内的人工粘性系数ε＝εjump。

步骤302、将单元内的人工粘性系数ε带入到方程(4)中，并对方程(4)迭代求解，得到流场的解。

如图2a所示，为本发明方法在高精度情况下对NACA0012翼型算例的激波捕捉问题的测试，从图2b可以看出，在高精度下，该方法能有效的捕捉激波，同时收敛所需的迭代步数明显减少，所用的CPU时间明显减少。

本发明并不局限于前述的具体实施方式。本发明扩展到任何在本说明书中披露的新特征或任何新的组合，以及披露的任一新的方法或过程的步骤或任何新的组合。