基于映射函数的中心型三阶格式

刘朋欣, 李 沁, 张涵信

(1. 中国空气动力研究与发展中心 空气动力学国家重点实验室, 四川 绵阳 621000;2. 中国空气动力研究与发展中心 计算空气动力研究所, 四川 绵阳 621000)

基于映射函数的中心型三阶格式

刘朋欣1,2,*, 李 沁1,2, 张涵信1,2

(1. 中国空气动力研究与发展中心 空气动力学国家重点实验室, 四川 绵阳 621000;2. 中国空气动力研究与发展中心 计算空气动力研究所, 四川 绵阳 621000)

高阶有限差分格式具有更小的数值耗散和色散误差，可以更好地捕捉流场的精细结构。在工程计算中，追求计算格式高分辨率的同时，也要兼顾其稳定性。Henrick首先提出映射函数的方法，作用于WENO5格式的非线性权，使其满足收敛阶数要求并提高了格式的性能。不同于Henrick的构造思路，李沁等提出了一种分段多项式形式的映射函数，改进后的WENO5格式在保证稳定性的前提下，拥有更好的分辨率。本文采用映射函数的思想，将一种五次分段多项式映射函数(Piecewise Polynomial Mapping Method, PPM)运用到中心型三阶格式上，构造出SWENO3-PPM5格式；将Henrick等提出的映射函数运用到中心型三阶格式上，构造出SWENO3-M格式，并与NND、WENO3格式的性能进行比较。近似色散关系分析表明，SWENO3-PPM5和SWENO3-M的色散和耗散关系相比于WENO3得到明显的提高。典型一维算例(Shu-osher问题)、二维算例(激波/层流边界层干扰问题、激波/旋涡干扰问题、绕25°/55°尖双锥流动问题)的计算结果表明，新格式具有比NND、WENO3更高的分辨能力，可以更好地捕捉到流场中的激波、接触间断面等复杂的流动现象，且由于分段多项式映射函数相比于Henrick提出的映射函数，可以使最终的非线性权值更加接近于格式的线性权，因此SWENO3-PPM5的分辨率要稍微高于SWENO3-M。

有限差分；映射函数；中心型三阶格式；分辨率；稳定性

0 引 言

随着所研究的流动问题越来越复杂，研究和应用对差分格式的计算精度也提出了更高的要求。相比于传统的二阶精度格式(如TVD格式)，高阶精度格式具有更小的数值耗散和色散误差，可以更好地刻画流场的精细结构。

高阶有限差分格式是高阶精度格式的重要组成部分。国外自从Harten提出ENO[1]和Liu、Osher和Chan等提出WENO[2]格式以来，众多学者从不同的思路出发，提出各种形式的高阶格式。其中，Jiang&Shu[3]提出了一种有效地实现WENO的方式，并将其用于有限差分，其五阶形式(WENO5)成为最流行的高阶格式之一，其非线性权求解过程及光滑度量因子(Indicator of Smoothness, IS)形式也逐渐成为后续研究的参考标准。尽管拥有很多的优势，但WENO5还存在一定的缺陷。Henrick等[4]首先指出，WENO5在极值点不能保持五阶精度，他们分析了保持精度的充分必要条件，并提出了一种精心设计的映射函数对格式的非线性权进行修正，相应的格式称为WENO-M。WENO-M可以保持精度一致，且分辨率有所提高。对非线性权进行改进并提出新格式的还有Borges[5]等提出的WENO-Z格式，Ha[6]等提出的WENO-NS格式，Ping Fan[7]等提出的WENO-η格式，Hui Feng[8]等提出的WENO-PM格式等。最近，Li等[9]从另一种构造映射函数的思路出发，独立地提出了一种分段多项式映射函数，对WENO5格式的非线性权进行修正，提出了WENO-PPM格式。新格式在保持计算稳定性的同时，进一步提高了格式的分辨率。国内张涵信[10]较早地提出了构造差分格式应该满足的四个原则，提出了一种无振荡、无自由参数耗散算法，即二阶精度的NND格式。由于NND形式简单、物理概念清楚、易于编程和计算量小，在国内工程领域得到了广泛的应用。随后，张涵信及其团队进一步发展了一系列高精度计算格式，包括ENN[11]、WNND[12]、WENN[13]、中心型三阶格式[14]等。

相比于基础和机理研究所采用的高阶格式，工程应用的需求通常具有以下特点：1) 实际问题中网格结构复杂，希望计算格式使用较小的网格模板；2) 由于实际问题流动的复杂性，要求格式具有类TVD格式的强鲁棒性；3) 格式具有较好的耗散特性，能够提升对旋涡等流动结构的描述能力。已有计算显示(如直升机旋翼的翼梢涡、涡流发生器产生流向涡的数值模拟等)，当采用TVD格式对这些结构进行数值模拟时，结构往往很快被耗散掉[14]。因此有必要进一步发展符合工程需求的三阶格式，在保证计算稳定性的前提下，改善计算的精度和分辨能力，且能够方便地推广到有限体积程序。

本文采用映射函数的思想，将五次分段多项式映射函数应用到中心型三阶格式上，构造了新型三阶格式SWENO3-PPM5；并与Henrick等提出的映射函数作对比(SWENO3-M)。通过对典型算例的计算，验证比较了新格式的性能。

1 数值计算方法

1.1 中心型三阶格式

李沁等[14]通过分析格式的耗散关系，提出了一类耗散可调的中心型三阶格式，提升了三阶格式的性能。本文基于此格式，做进一步的改进。首先简单介绍中心型三阶格式。

以一维双曲守恒律方程为例：

考虑守恒型差分格式，式(1)的数值求解可以写成如下半离散化形式：

数值通量表示为：

本文计算中取p=2，ε=10-16。其中，ISk代表第k个模板的光滑度量因子。参考WENO系列格式非线性加权的基本框架ISk具体形式为：

由于中心型三阶格式基于NND格式增加了一个下风模板，格式的线性权系数中就引入了一个自由参数，记为n3rd。线性权系数见式(8)。在约化波数为π处，中心型三阶线性格式的耗散为WENO3线性格式的1/n3rd。通过改变n3rd的值就可以调整中心型三阶格式的耗散。具体分析可参考文献[14]。本文的计算中取n3rd=10。

为了增强格式的稳定性，在使用(7)式计算得到各个模板的光滑因子IS后，要对最下风模板进行一次修正：

1.2 分段多项式映射函数

通过使用映射函数来修正原来的非线性权ωk，可以使WENO5格式的精度保持五阶。这个函数定义域为[0,1]，且满足：

gk(0)=0,gk(1)=1

此函数将使用式(6)计算得到的非线性权ωk进行一次映射，再将得到的值进行归一化作为最终的非线性权。图1给出了Henrick型映射权函数与分段多项式型映射函数的分布。可以看出，相比于Henrick型映射权函数，采用五次分段多项式映射函数方法使非线性权在ω=Ck处的分布更加平坦，即更加接近于格式的理想线性权值，这将有利于提高格式分辨精细结构的能力。

本文将分段多项式映射函数(式(14))运用于中心型三阶格式，相应的格式记为SWENO3-PPM5。并与采用Henrick等提出的映射函数的三阶格式作比较，记为SWENO3-M格式。

1.3 不同格式的近似色散关系(ADR)对比

图2给出了三种非线性格式的近似色散关系[15](Approximate Dispersion Relations, ADR)对比，可以看到：采用映射函数方法的中心型三阶格式的色散和耗散关系得到明显的提高，且SWENO3-PPM5的表现略优于SWENO3-M。

2 数值计算结果与分析

本节计算了典型的一维和二维算例来检测新型三阶格式的性能，对比了SWENO3-PPM5、SWENO3-M、WENO3、NND四种格式。所有的算例均使用Steger-Warming分裂方式，并在计算中使用特征投影。时间推进方面，对于非定常问题使用具有TVD性质的三阶Runge-Kutta格式，对于定常问题则使用LU-SGS格式。

2.1 Shu-Osher问题

该算例描述了一维情形下马赫数为3.0的右行激波和密度波发生干扰作用，其目的在于分辨主激波波后的系列小尺度结构，通常用于验证格式的分辨率。初始条件为：

(ρ,u,p)=

计算采用801个网格点的均匀网格，结果如图3和图4所示。图中把使用WENO5在1601个网格点上的计算结果作为准确解。从计算结果来看，由于中心型三阶格式的耗散低于WENO3格式，且使用映射函数可以使格式的权值更加接近于线性权，因此SWENO3-PPM5和SWENO3-M的分辨率要高于WENO3格式和NND格式；且由于PPM5映射函数相比于Henrick提出的函数，可以使非线性权在ω=Ck处的分布更加平坦，即更加接近于线性权，因此SWENO3-PPM5的分辨率也要略高于SWENO3-M。

2.2 激波/层流边界层干扰问题

激波/层流边界层干扰流动包含丰富的流动结构，常作为检验格式性能的算例之一。来流状态为M∞=2.0，T∞=293K，Re∞=2.96×105，Pr=0.72，β=32.585°，β为激波入射角；壁面为黏性绝热条件；进出口皆为超声速边界。计算域为：0.0≤x≤2.02，0.0≤y≤1.30，网格为nx×ny=102×121，x方向均匀分布，y方向在壁面附近进行网格加密。

为了研究这几种格式分辨膨胀波系的能力，本文取y=0.3045处的压力分布(即图5中红线位置)，并

与WENO5的计算结果进行比较，如图6所示。可以看出，在x=1.4附近，NND格式完全不能分辨该区域；WENO3格式亦不能有效分辨该膨胀波系；虽然与模板较宽的WENO5格式的结果有一定差距，但是SWENO3-PPM5和SWENO3-M格式也能够较好分辨到此区域，分辨率明显优于NND格式和WENO3格式。

2.3 激波/旋涡干扰问题

图7给出了四种格式计算得到的密度等值线。可以看到，SWENO3-PPM5和SWENO3-M格式分辨出的流场结构明显比WENO3和NND更加丰富。由于初始流场含有激波间断，为了适应Euler方程的解算器，在x=-3.6处会形成垂直方向的二次干扰波。NND和WENO3格式由于耗散较大，并没有较好地捕捉到此干扰波；而经过映射后的中心型三阶格式则较为清晰地分辨出此处的干扰波。同样还有图中方框区域内的流场结构，可以看出，SWENO3-PPM5分辨出的结构较SWENO3-M更加细致。

2.4 绕25°/55°尖双锥流动问题

绕尖双锥流动是一类典型的高超声速绕流。本文对绕25°/55°尖双锥流动进行模拟，来流条件为M∞=9.59，T∞=185.6K，Tw=293.3K，Re=139436m-1。网格为nstream×nnormal=262×148。计算结果与Gnoffo[17]给出的实验和数值模拟结果进行对比。Gnoffo数值模拟时使用的网格是512×256。

图8给出了几种格式计算得到的密度梯度等值线和流线分布。从密度梯度等值线可以看出，几种激波捕捉格式基本上都捕捉到了流场中的激波、接触间断面、超声速射流区、分离区等复杂的流动现象，但是SWENO3-PPM5和SWENO3-M分辨出的结构要更加清晰一些，尤其是在分离区及超声速射流区处；且在分离区，SWENO3-PPM5分辨出的旋涡数最多，NND格式最少。

图9和图10比较了几种格式计算的壁面压力和热流分布。可以看出，SWENO3-PPM5的计算出的分离涡较大，且壁面压力和热流分布更加接近于Gnoffo的计算和实验结果，显示了SWENO3-PPM5较好的热流预测能力。

3 结 论

本文采用映射函数的思想，将一种五次分段多项式映射函数运用到中心型三阶格式上，构造了SWENO3-PPM5格式，同时与Henrick等提出的映射函数作对比(SWENO3-M格式)。Shu-Osher算例、激波/层流边界层干扰、激波/旋涡干扰以及绕25°/55°尖双锥流动的计算表明，使用映射函数的新三阶格式具有比NND、WENO3更高的分辨能力，可以有效地捕捉到流场中的激波、接触间断面等复杂的流动现象。且由于较Henrick的映射函数，分段多项式映射函数可以使最终的非线性权更加接近于格式的线性权，因此SWENO3-PPM5的分辨能力也要稍微高于SWENO3-M。

[1]Harten A, Enquist B, Osher S, et al. Uniformly high order essentially non-oscillatory schemes[J]. Journal of Computational Physics, 1987, 71: 231-303.

[2]Liu X, Osher S, Chan T. Weighted essentially non-oscillatory schemes[J]. Journal of Computational Physics, 1994, 115: 200-212.

[3]Jiang G, Shu Qiwang. Efficient implementation of weighted ENO schemes[J]. Journal of Computational Physics, 1996, 126: 202-228.

[4]Henrick A K, Aslam T D, Powers J M. Mapped weighted essentially non-oscillatory schemes: achieving optimal order near critical points[J]. Journal of Computational Physics, 2005, 207: 542-567.

[5]Borges R, Carmona M, Costa B, et al. An improved weighted essentially non-oscillatory scheme for hyperbolilc conservation laws[J]. Journal of Computational Physics, 2008, 227: 3191-3211.

[6]Ha Y, Kim C H, Lee Y J, et al. An improved weighted essentially non-oscillatory scheme with a new smoothness indicator[J]. Journal of Computational Physics, 2013, 232: 68-86.

[7]Ping Fan, Yiqing Shen, Naolin Tian, et al. A new smoothness indicator for improving the weighted essentially non-oscillatory scheme[J]. Journal of Computational Physics, 2014, 269: 329-354.

[8]Hui Feng, Fuxing Hu, Rong Wang. A new mapped weighted

essentially non-oscillatory scheme[J]. Journal of Scientific Computing, 2012, 51: 449-473.

[9]Li Q, Liu P X, Zhang H X. Piecewise polynomial mapping method and corresponding WENO scheme with improved resolution[J]. Communication in Computational Physics, 2015, 18(5): 1417-1444.

[10]Zhang H X. Non-oscillatory and non-free-parameter dissipation difference scheme[J]. Acta Aerodynamica Sinica, 1988, 6(2): 143-165.(in Chinese)张涵信. 无波动、无自由参数的耗散差分格式[J]. 空气动力学学报, 1988, 6(2): 143-165.

[11]Zhang H X, He G H, Zhang L. Some important problems for high order accurate difference scheme solving gas dynamic equations[J]. Acta Aerodynamica Sinica, 1993, 11(4): 341-356.(in Chinese)张涵信, 贺国宏, 张雷, 等. 高精度求解气动方程的几个问题[J]. 空气动力学学报, 1993, 11(4): 341-356.

[12]Liu W. Nonlinear dynamics analysis for mechanism of slender wing rock and study of numerical simulation method[D]. National University of Defense Technology, 2004.(in Chinese)刘伟. 细长机翼摇滚机理的非线性动力学分析及数值模拟方法研究[D]. 湖南: 国防科学技术大学, 2004.

[13]Li Q, Zhang H X, Gao S C. Numerical simulations on supersonic shear layer flow[J]. Acta Aerodynamica Sinica, 2000, 18(S): 67-77. (in Chinese)李沁, 张涵信, 高树椿. 关于超声速剪切流动的数值模拟[J]. 空气动力学学报, 2000, 18(增刊): 67-77.

[14]Li Q, Sun D, Zheng Y K, et al. On a class of center-type third order difference scheme orienting to engineering utilizations[J]. Acta Aerodynamica Sinica, 2013, 31(4): 466-472. (in Chinese)李沁, 孙东, 郑永康, 等. 一类中心型三阶格式及应用[J]. 空气动力学学报, 2013, 31(4):466-472.

[15]Sergio Pirozzoli. On the spectral properties of shock-capturing schemes[J]. Journal of Computational Physics, 2006, 219: 489-497.

[16]Inoue O, Hattori Y. Sound generation by shock vortex interactions[J]. Journal of Fluid Mechanics, 1990, 380: 81-116.

[17]Gnoffo P. CFD validation studies for hypersonic flow prediction[R]. AIAA 2001-1025.

A kind of third order central scheme based on mapping functions

Liu Pengxin1,2,*, Li Qin1,2, Zhang Hanxin1,2

(1.StateKeyLaboratoryofAerodynamicsofChinaAerodynamicsResearchandDevelopmentCenter,Mianyang621000,China;2.ComputationalAerodynamicsInstituteofChinaAerodynamicsResearchandDevelopmentCenter,Mianyang621000,China)

High-order finite difference schemes have low numerical dissipation and dispersion error. They can capture flow structures in detail. In engineering applications, robustness should be considered, while high resolution is pursued. The method of mapping function was first proposed by Henrick to adjust nonlinear weights in [0,1] for the fifth-order WENO scheme. With the help of this method, the requirement of convergence order is satisfied, and the performance of the scheme is improved. Different from Henrick’s method, a concept of piecewise polynomial function is proposed by authors in this paper, and the corresponding improved WENO schemes have better resolution with the preservation of robustness. In this study, a fifth-order piecewise polynomial function (PPM5) is applied to a third order central scheme, which is referred as the SWENO3-PPM5. Henrick's mapping function is also applied to the third order central scheme, which is referred as the SWENO3-M. The performance of the SWENO3-PPM5 and the SWENO3-M are compared with that of the NND and the WENO3. The analyses of approximate dispersion relations (ADR) indicate that the dissipation and dispersion relations of SWENO3-PPM5 and SWENO3-M have an improvement compared with those of the WENO3. The performance of the new schemes is evaluated by testing typical problems, i.e. one-dimensional Shu-osher problem, two-dimensional ones including shock-boundary layer interaction problem, shock/vortex interaction problem, and hypersonic flows of the sharp double cone. The results indicate that the new schemes have better resolution than the WENO3 and the NND and can present a better description on complex flow structures, such as shock waves and contact discontinuity. Because the PPM5 can make the final weight closer to the idea weight of the scheme, the SWENO3-PPM5 possesses a little better resolution than the SWENO3-M.

finite difference scheme; mapping functions; third order central scheme; resolution; robustness

0258-1825(2017)01-0071-07

2015-09-01;

2016-04-13

国家重点基础研究发展计划(2014CB744100)，国家自然科学基金重大研究计划(91541105)

刘朋欣*(1990-)，男，河南商丘人，硕士，主要从事计算流体力学及爆震燃烧研究. E-mail：liupengxin2008@163.com

刘朋欣, 李沁, 张涵信. 基于映射函数的中心型三阶格式[J]. 空气动力学学报, 2017, 35(1): 71-77.

10.7638/kqdlxxb-2015.0172 Liu P X, Li Q, Zhang H X. A kind of third order central scheme based on mapping functions[J]. Acta Aerodynamica Sinica, 2017, 35(1): 71-77.

V211.3

A doi: 10.7638/kqdlxxb-2015.0172