CN114638173A - 一种高阶非线性激波捕捉空间离散方法 - Google Patents
一种高阶非线性激波捕捉空间离散方法 Download PDFInfo
- Publication number
- CN114638173A CN114638173A CN202210086558.0A CN202210086558A CN114638173A CN 114638173 A CN114638173 A CN 114638173A CN 202210086558 A CN202210086558 A CN 202210086558A CN 114638173 A CN114638173 A CN 114638173A
- Authority
- CN
- China
- Prior art keywords
- flux
- format
- grid point
- grid
- order
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 230000035939 shock Effects 0.000 title claims abstract description 45
- 238000000034 method Methods 0.000 title claims abstract description 41
- 239000006185 dispersion Substances 0.000 title claims abstract description 16
- 230000004907 flux Effects 0.000 claims abstract description 119
- 238000004364 calculation method Methods 0.000 claims abstract description 46
- 238000009499 grossing Methods 0.000 claims description 6
- 230000008901 benefit Effects 0.000 abstract description 3
- 238000011438 discrete method Methods 0.000 abstract description 2
- 238000011160 research Methods 0.000 description 7
- 230000008878 coupling Effects 0.000 description 6
- 239000011159 matrix material Substances 0.000 description 4
- 230000008569 process Effects 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 230000010355 oscillation Effects 0.000 description 3
- 239000000853 adhesive Substances 0.000 description 2
- 238000010168 coupling process Methods 0.000 description 2
- 238000005859 coupling reaction Methods 0.000 description 2
- 239000012530 fluid Substances 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 230000007246 mechanism Effects 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 230000003068 static effect Effects 0.000 description 2
- 238000012795 verification Methods 0.000 description 2
- 241000287196 Asthenes Species 0.000 description 1
- 241001463014 Chazara briseis Species 0.000 description 1
- 230000001070 adhesive effect Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000004422 calculation algorithm Methods 0.000 description 1
- 230000001413 cellular effect Effects 0.000 description 1
- 230000002860 competitive effect Effects 0.000 description 1
- 238000012937 correction Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000000737 periodic effect Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000010183 spectrum analysis Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T3/00—Geometric image transformations in the plane of the image
- G06T3/40—Scaling of whole images or parts thereof, e.g. expanding or contracting
- G06T3/4007—Scaling of whole images or parts thereof, e.g. expanding or contracting based on interpolation, e.g. bilinear interpolation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/08—Fluids
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- Fluid Mechanics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Computing Systems (AREA)
- Pure & Applied Mathematics (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- Geometry (AREA)
- General Engineering & Computer Science (AREA)
- Algebra (AREA)
- Complex Calculations (AREA)
Abstract
本发明适用于流体力学计算领域,提供了一种高阶非线性激波捕捉空间离散方法,包括如下步骤:读取流场数据,得到各网格点上的物理通量和特征通量;根据所述特征通量判断各网格点的间断性,所述间断性包括光滑和间断;当所述网格点的间断性为光滑时,通过k次样条插值格式得到半网格点上通量;当所述网格点的间断性为间断时,通过单调保持格式得到半网格点上通量;根据所述半网格点上通量,并结合各网格点上的物理通量,得到整体的空间离散格式;将时间推进至预设结束时间,得到预设结束时间的流场数据。本发明提供的一种高阶非线性激波捕捉空间离散方法具有低耗散、高精度、高分辨率和间断的精确捕捉等优势。
Description
技术领域
本发明涉及流体力学计算领域,尤其是涉及一种高阶非线性激波捕捉空间离散方法。
背景技术
随着航空技术飞速发展,先进飞行器的气动噪声问题越来越受到广泛重视。航空界把飞行器的声舒适性和对环境噪声影响作为研究的出发点,尽量降低飞行器的气动噪声,提高商用飞机的竞争能力。求解包含激波的气动声学问题是一个重要的研究方向,数值模拟这类问题对计算格式提出了更高的要求:求解声场要求数值格式具有高阶精度、高分辨率和低耗散的特性,而模拟激波运动还要求数值格式必须能够捕捉激波。目前,还不存在比较理想的模拟这类激波噪声问题的数值格式,但激波噪声却是航天飞行器噪声的一个重要来源,具有重要的工程应用背景和理论研究价值。因此非常有必要研究适合于模拟激波噪声问题的数值格式,从而为我们开展激波噪声问题的研究提供帮助,并用于激波噪声产生、传播过程的数值模拟,促进对物理本质的认识,从而对激波噪声的降低与控制提供新的手段。
其中,控制方程中对流项的数值离散会直接影响无粘流区域的激波的分辨,并间接影响激波噪声的预测。由于中心紧致差分格式没有数值耗散并且色散误差小,所以对于计算多尺度问题有着天然的优势。但是直接把高阶中心紧致格式应用到非线性双曲守恒律也不现实,因为差分离散格式不能分辨的高频波将会累积误差,最后导致数值不稳定。
传统的激波捕捉格式WENO格式虽然能够达到较好的分辨率,但其计算模板上节点数据较多,为边界处理带来了一定程度的困难,且增加了计算工作量,耗散较高,在光滑解极值点并没有达到目标精度,且不能够很好的计算剪切层等精细结构,并且WENO格式对噪声的识别能力不够准确。
直接采用线性中心紧致差分格式应用到非线性双曲守恒律,则差分离散格式不能分辨的高频波将会累积误差,最后导致数值不稳定。为了克服中心紧致格式的数值不稳定缺陷,加入人工耗散来消去非物理的波动,这样的经验性很强,不好控制。
发明内容
本发明的目的是解决当流场存在激波间断时出现非物理振荡而导致的数值不稳定性的技术问题,提供了一种高阶非线性激波捕捉空间离散方法,具有低耗散、高精度、高分辨率和间断的精确捕捉等优势。
本发明提供了一种高阶非线性激波捕捉空间离散方法,包括如下步骤:
步骤S10:读取流场数据,得到各网格点上的物理通量和特征通量;
步骤S20:根据所述特征通量判断各网格点的间断性,所述间断性包括光滑和间断;
步骤S30:当所述网格点的间断性为光滑时,通过k次样条插值格式得到半网格点上通量,其中,k为样条插值格式的最高次数项,k为正整数且k≥3;当所述网格点的间断性为间断时,通过单调保持格式得到半网格点上通量;
步骤S40:根据所述半网格点上通量,并结合各网格点上的物理通量,得到整体的空间离散格式;
步骤S50:将时间推进至预设结束时间,得到预设结束时间的流场数据。
进一步的,所述初始单调保持通量uMP的计算公式为uMP=uj+midmod[uj+1-uj,4(uj-uj-1)],其中,uj+1为第j+1个网格点的特征通量,uj-1为第j-1个网格点的特征通量。
进一步的,所述半网格点的光滑因子的计算公式为其中,Mj+1为第j+1个网格点上k次样条插值函数的k-1次导数,Mj为第j个网格点上k次样条插值函数的 k-1次导数,k次样条插值函数为S(j),a≤j≤b,b-a小于预设区间,h为相邻网格点之间的距离。
进一步的,步骤S30中,当所述网格点的间断性为光滑时,通过三次样条插值格式得到半网格点上通量。
进一步的,步骤S10中所述特征通量的计算方法为:根据控制方程计算各网格点上的初始通量,通过对所述初始通量分裂得到正负通量,对所述正负通量进行特征投影,得到特征通量。
进一步的,步骤S40中所述空间导数通过线性紧致差分格式得到,所述紧致差分格式为其中,a为第一系数,b为第二系数,α为第三系数,f′j-1为第j-1个网格点上物理通量的导数的差分近似式,f′j为第j个网格点上物理通量的导数的差分近似式,f′j+1为第j+1个网格点上物理通量的导数的差分近似式,为第j个网格点的右半网格点上通量的反投影,为第j个网格点的左半网格点上通量的反投影,fj+1为第j+1 个网格点上物理通量,fj-1为第j-1个网格点上物理通量;
进一步的,步骤S50中采用三阶龙格库塔方法对时间项进行离散,将时间推进至预设结束时间。
综上所述,本发明至少能够实现如下技术效果:
1.对于光滑流场,多次样条插值耦合中心差分格式可以很精确地计算流场。本发明在当网格点的间断性为光滑时,通过k次样条插值格式得到半网格点上通量,实现了流场单元界面处通量的的快速求解,同时提高了光滑流场的计算精度;
2.当流场存在间断时,非物理振荡可能导致数值不稳定性,为了能够计算包含激波等的复杂多尺度问题,将多次样条插值后的半网格点通量采用单调保持格式来修正,获得对激波等间断问题的精确捕捉,实现对超声速问题高精度、高分辨率的计算;
3.本发明还提出了一种新的网格点间断性的判断方法,能够直接由空间格式自身中计算得到,具有节省计算步骤,计算效率高的特点;
4.本发明通过单调保持非线性样条插值的办法来求解半网格点处的函数值,从而得到了高阶非线性紧致格式,所述高阶非线性紧致格式保持了原有中心紧致差分格式的高阶精度、高分辨率和低耗散特性,同时由于其单调保持非线性插值的机制,使得高阶非线性紧致格式能够捕捉强激波,一系列数值算例表明所述样条插值非线性紧致格式是高阶精度、高分辨率和低耗散数值计算方法,具有广泛的科研以及应用前景。
附图说明
为了更清楚地说明本发明实施例的技术方案,下面将对本发明实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面所描述的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1是本发明中高阶非线性激波捕捉空间离散方法的示意图;
图2是Exact、MPSI4、CNCS4和WENO5的分辨率曲线;
图3是本发明中MPSI4格式计算Lax问题;
图4是本发明中MPSI4格式计算激波双马赫反射;
图5是本发明中MPSI4格式计算前台阶问题;
图6是本发明中MPSI4格式计算气动噪声传播问题;
图7是本发明中MPSI4格式计算气动噪声传播的噪声脉动分布。
具体实施方式
以下的说明提供了许多不同的实施例、或是例子,用来实施本发明的不同特征。以下特定例子所描述的元件和排列方式,仅用来精简的表达本发明,其仅作为例子,而并非用以限制本发明。
下面将结合本发明实施例中附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。通常在此处附图中描述和示出的本发明实施例的组件可以以各种不同的配置来布置和设计。因此,以下对在附图中提供的本发明的实施例的详细描述并非旨在限制要求保护的本发明的范围,而是仅仅表示本发明的选定实施例。基于本发明的实施例,本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例,都属于本发明保护的范围。
应注意到:相似的标号和字母在下面的附图中表示类似项,因此,一旦某一项在一个附图中被定义,则在随后的附图中不需要对其进行进一步定义和解释。同时,在本发明的描述中,术语“第一”、“第二”等仅用于区分描述,而不能理解为指示或暗示相对重要性。
实施例一:
如图1所示,本发明实施例一提供了一种高阶非线性激波捕捉空间离散方法,包括如下步骤:
步骤S10:读取流场数据,得到各网格点上的物理通量和特征通量;
步骤S20:根据所述特征通量判断各网格点的间断性,所述间断性包括光滑和间断;
步骤S30:当所述网格点的间断性为光滑时,通过k次样条插值格式得到半网格点上通量,其中,k为样条插值格式的最高次数项,k为正整数且k≥3;当所述网格点的间断性为间断时,通过单调保持格式得到半网格点上通量;
步骤S40:根据所述半网格点上通量,并结合各网格点上的物理通量,得到整体的空间离散格式;
步骤S50:将时间推进至预设结束时间,得到预设结束时间的流场数据。
对于光滑流场,多次样条插值耦合中心差分格式可以很精确地计算流场。所述k次样条插值格式的最高次数项k仅需要满足k为正整数且k≥3即可,具体的最高次数项k的数值越大对间断的捕捉精度越高,得到的分辨率也会越高,如五次样条插值格式的计算精度>四次样条插值格式的计算精度>三次样条差值格式的计算精度。
本实施例在当网格点的间断性为光滑时,通过k次样条插值格式得到半网格点上通量,实现了流场单元界面处通量的的快速求解,同时提高了光滑流场的计算精度。
当流场存在间断时,非物理振荡可能导致数值不稳定性,为了能够计算包含激波等的复杂多尺度问题,将多次样条插值后的半网格点通量采用单调保持格式来修正,获得对激波等间断问题的精确捕捉,实现对超声速问题高精度、高分辨率的计算。
所述网格点间断性的判断方法能够直接由空间格式自身中计算得到,具有节省计算步骤,计算效率高的特点。
进一步的,所述初始单调保持通量uMP的计算公式为 uMP=uj+mid mod[uj+1-uj,4(uj-uj-1)],其中,uj+1为第j+1个网格点的特征通量,uj-1为第j-1个网格点的特征通量。
进一步的,所述半网格点的光滑因子的计算公式为其中,Mj+1为第j+1个网格点上k次样条插值函数的k-1次导数,Mj为第j个网格点上k次样条插值函数的 k-1次导数,k次样条插值函数为S(j),a≤j≤b,b-a小于预设区间,h为相邻网格点之间的距离。
进一步的,步骤S30中,当所述网格点的间断性为光滑时,通过三次样条插值格式得到半网格点上通量。
对于k次样条插值格式的最高次数项k,虽然最高次数项k的数值越大对间断的捕捉精度越高,得到的分辨率也会越高,但最高次数项k的数值越大计算过程也将更加复杂,运算效率也会降低。因而,k次样条插值格式的最高次数项k优选为3,此时,既能够实现高精度、高分辨率,也可确保计算过程的简单化和运算的高效率。
以半网格点j+1/2为例,本实施例中以三次样条插值为例,构造半网格点上通量的高阶混合计算方法计算步骤如下:
将三次样条插值定义为函数S(x)∈C2[a,b],在每个小区间上是三次多项式,其中a=x0<x1<…<xn-1<xn=b是给定的网格节点。若函数S(x)在节点上满足以下公式,则称S(x)是三次样条函数:
S0(x0)=f(x0),Sn-1(xn)=f(xn)
Sj-1(xj)=f(xj)=Sj(xj),j=1,…,n-1
S′j-1(xj)=S′j(xj)
S″j-1(xj)=S″j(xj)
上述方程共有4n-2个,为求解所有Sj(xj)的系数,还需两个方程,通常根据物理边界条件来得到,具体如下:
假设S″(xj)=Mj,则得到:
其中hj=xj+1-xj。对上式积分两次,则可以得到如下两个函数:
利用一阶导数的接触面条件S′j-1(xj)=S′j(xj)可得:
由上式可得三对角方程组:μjMj-1+2Mj+λjMj+1=dj,j=1,…,n-1
上述方程组联合边界条件可以求解出S″(xj)=Mj,从而可以得到:
当计算网格是均匀网格时,可以简化为:
dj=uj-1-2uj+uj+1;
min mod(x1,…,xk)=s·min(|x1|,|x2|,…,xk∣)
median(x,y,z)=x+min mod(y-x,z-x)
进一步的,步骤S10中所述特征通量的计算方法为:根据控制方程计算各网格点上的初始通量,通过对所述初始通量分裂得到正负通量,对所述正负通量进行特征投影,得到特征通量。
控制方程为在不考虑质量力和源项的情况下,无量纲积分形式的三维Navier-Stokes控制方程:
在牛顿流体中,粘性应力张量各个分量分别为:
Λk为特征矩阵的特征值矩阵,特征值代表了各类波的特征速度,其具体形式如下:
Lk和Rk分别为左、右特征矩阵,具体形式为:
在数值计算中,根据特征传播方向,通常将无粘通量分为正、负两个部分,上标“+”表示这部分的无粘矢通量的特征速度指向计算坐标轴的正方向,上标“-”则表示该部分的特征速度指向计算坐标轴的负方向传播。
进一步的,步骤S40中所述空间导数通过线性紧致差分格式得到,所述紧致差分格式为其中,a为第一系数,b为第二系数,α为第三系数,f′j-1为第j-1个网格点上物理通量的导数的差分近似式,f′j为第j个网格点上物理通量的导数的差分近似式,f′j+1为第j+1 个网格点上物理通量的导数的差分近似式,为第j个网格点的右半网格点上通量的反投影,为第j个网格点的左半网格点上通量的反投影,fj+1为第j+1个网格点上物理通量,fj-1为第j-1个网格点上物理通量;
中心紧致格式中的网格中心处的值采用三次样条插值可以得到四阶紧致格式,进而得到整个的空间离散格式,将本实施例中的空间离散格式记为四阶 MPSI(Monotonicitypreserving Spline Interpolation)格式,简称MPSI4。
当k次样条插值格式的最高次数项k大于等于3时,则得到k+1阶紧致格式,此时空间离散格式记为k阶MPSI(Monotonicity preserving Spline Interpolation)格式,简称MPSIk。
进一步的,步骤S50中采用三阶龙格库塔方法对时间项进行离散,将时间推进至预设结束时间。
采用显式Runge-Kutta时间推进方法来离散时间导数项。
优化的三阶TVD Runge-Kutta时间离散格式为: U(1)=Un+ΔtL(Un)
其中,n为第n时刻步的值。至此,完成了对欧拉或者纳维斯托克斯方程的高阶时空离散,从而将时间从开始时间tini一直推进至预设结束时间tend。
实施例二:
本发明实施例二采用标准算例来验证四阶MPSI格式的精度,令ut+ux=0,其中-1≤x≤1。设定初始条件为u(x,t=0)=sin(πx),精确解为u(x,t)=sin(π(x-t)),给定周期性的边界条件,计算t=1。
求出数值解和精确解在不同网格密度下的L1误差,L1数值精度阶数,L∞误差,L∞数值精度阶数。时间离散采用三阶TVD Runge-Kutta格式,数值精度的验证如表1所示,其中,N为网格数,从表1中可以看到MPSI4格式确实达到了理论的收敛精度。
表1四阶MPSI格式和三阶TVD Runge-Kutta格式的数值精度验证
N | L<sub>1</sub>误差 | L<sub>1</sub>数值精度阶数 | L<sub>∞</sub>误差 | L<sub>∞</sub>数值精度阶数 |
20 | 0.109E-03 | -- | 0.172E-03 | -- |
40 | 0.678E-05 | 4.01 | 0.107E-04 | 4.01 |
60 | 0.134E-05 | 4.00 | 0.210E-05 | 4.00 |
80 | 0.423E-06 | 4.00 | 0.665E-06 | 4.00 |
100 | 0.173E-06 | 4.00 | 0.272E-06 | 4.00 |
同时,采用非线性格式的频谱分析方法来研究MPSI4格式的分辨率。如图 2所示为Exact、MPSI4、CNCS4和WENO5的分辨率曲线,表明MPSI4格式的分辨率远优于传统的紧致差分格式和激波捕捉WENO格式。
(1)Lax激波管问题
如图3所示为采用MPSI4格式计算Lax问题(N=100),该问题的初始条件为:
计算区域取为-5≤x≤5,计算网格采用100个点,初始间断位于x=0处,最终计算时刻取为t=1.3。
采用MPSI4格式计算一维激波管Lax问题,计算结果在激波、膨胀波及接触间断附近都能达到基本无波动。
(2)激波双马赫反射问题
如图4所示为采用MPSI4格式计算激波双马赫反射(N=960×240),双马赫反射问题描述的是Mach数为10的强运动斜激波以与x轴方向呈60度角的方向入射,入射点在(1/6,0),计算区域取[0,4]×[0,1]。激波波前静止空气的密度为1.4,压力为1,波后按激波关系给定初始条件。下边界在[1/6,4]的区域给定壁面边界条件,其它边界区域按照激波运动所在的位置,分别给定波前或波后的值。初始物理参数为:
本发明实施例中采用960×240的均匀网格计算到无量纲时间t=0.2的时刻。计算结果在激波附近都能达到基本无波动,对剪切层的分辨率也非常高。
(3)前台阶问题
如图5所示为采用MPSI4格式计算前台阶问题(N=960×320)。前台阶问题的初始条件为Mach数为3的强运动激波沿着带有台阶的风洞向右传播。台阶高度为0.2,前缘位置在0.6,计算域为[0,3]×[0,1]。初始物理参数为:
(ρ,u,v,p)=(1.4,3.0,0.0,1.0)
左、右边界为开边界,其余边界取壁面反射条件。采用960×320的均匀网格计算截止到无量纲时间t=4的时刻。
计算结果在激波附近都能达到基本无波动,对剪切层的分辨率也非常高。
(4)二维气动噪声传播问题
计算区域取为x∈[-100,100],y∈[-100,100],x方向的马赫数Mx=0.5,y方向的马赫数My=0,初始时刻的物理量为:
采用MPSI4格式计算二维气动噪声传播问题,如图6所示为采用MPSI4 格式计算气动噪声传播问题;如图7所示为采用MPSI4格式计算气动噪声传播的噪声脉动分布。从计算结果与精确解的比较可以看出:MPSI4格式耗散性非常小,长时间气动噪声模拟能力很强。
本发明提供的一种高阶非线性激波捕捉空间离散方法,通过单调保持非线性样条插值的办法来求解半网格点处的函数值,从而得到了高阶非线性紧致格式MPSIk,所述高阶非线性紧致格式MPSIk保持了原有中心紧致差分格式的高阶精度、高分辨率和低耗散特性,同时由于其单调保持非线性插值的机制,使得高阶非线性紧致格式能够捕捉强激波,一系列数值算例表明所述样条插值非线性紧致格式是高阶精度、高分辨率和低耗散数值计算方法,具有广泛的科研以及应用前景。
以上所述仅为本发明的较佳实施例而已,并不用以限制本发明,凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明的保护范围之内。
Claims (10)
1.一种高阶非线性激波捕捉空间离散方法,其特征在于,包括如下步骤:
步骤S10:读取流场数据,得到各网格点上的物理通量和特征通量;
步骤S20:根据所述特征通量判断各网格点的间断性,所述间断性包括光滑和间断;
步骤S30:当所述网格点的间断性为光滑时,通过k次样条插值格式得到半网格点上通量,其中,k为样条插值格式的最高次数项,k为正整数且k≥3;当所述网格点的间断性为间断时,通过单调保持格式得到半网格点上通量;
步骤S40:根据所述半网格点上通量,并结合各网格点上的物理通量,得到整体的空间离散格式;
步骤S50:将时间推进至预设结束时间,得到预设结束时间的流场数据。
3.如权利要求2所述的一种高阶非线性激波捕捉空间离散方法,其特征在于,所述初始单调保持通量uMP的计算公式为uMP=uj+mid mod[uj+1-uj,4(uj-uj-1)],其中,uj+1为第j+1个网格点的特征通量,uj-1为第j-1个网格点的特征通量。
5.如权利要求1-4之一所述的一种高阶非线性激波捕捉空间离散方法,其特征在于,步骤S30中,当所述网格点的间断性为光滑时,通过三次样条插值格式得到半网格点上通量。
8.如权利要求1所述的一种高阶非线性激波捕捉空间离散方法,其特征在于,步骤S10中所述特征通量的计算方法为:根据控制方程计算各网格点上的初始通量,通过对所述初始通量分裂得到正负通量,对所述正负通量进行特征投影,得到特征通量。
9.如权利要求1所述的一种高阶非线性激波捕捉空间离散方法,其特征在于,步骤S40中所述空间导数通过线性紧致差分格式得到,所述紧致差分格式为其中,a为第一系数,b为第二系数,α为第三系数,fj′-1为第j-1个网格点上物理通量的导数的差分近似式,fj′为第j个网格点上物理通量的导数的差分近似式,fj′+1为第j+1个网格点上物理通量的导数的差分近似式,为第j个网格点的右半网格点上通量的反投影,为第j个网格点的左半网格点上通量的反投影,fj+1为第j+1个网格点上物理通量,fj-1为第j-1个网格点上物理通量;
10.如权利要求1所述的一种高阶非线性激波捕捉空间离散方法,其特征在于,步骤S50中采用三阶龙格库塔方法对时间项进行离散,将时间推进至预设结束时间。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210086558.0A CN114638173B (zh) | 2022-01-25 | 2022-01-25 | 一种高阶非线性激波捕捉空间离散方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210086558.0A CN114638173B (zh) | 2022-01-25 | 2022-01-25 | 一种高阶非线性激波捕捉空间离散方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114638173A true CN114638173A (zh) | 2022-06-17 |
CN114638173B CN114638173B (zh) | 2023-06-02 |
Family
ID=81946138
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210086558.0A Active CN114638173B (zh) | 2022-01-25 | 2022-01-25 | 一种高阶非线性激波捕捉空间离散方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114638173B (zh) |
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116384288A (zh) * | 2023-06-05 | 2023-07-04 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种可压缩流动高分辨率数值模拟方法、介质及设备 |
CN116542184A (zh) * | 2023-07-05 | 2023-08-04 | 中国空气动力研究与发展中心计算空气动力研究所 | 粘性通量的计算方法、装置、终端设备和存储介质 |
CN116611369A (zh) * | 2023-07-17 | 2023-08-18 | 中国空气动力研究与发展中心计算空气动力研究所 | 基于光滑度量量级与候选模板点数的插值方法及装置 |
CN118313313A (zh) * | 2024-06-04 | 2024-07-09 | 中国人民解放军国防科技大学 | 基于非线性谱差分的高阶精度流场计算方法 |
CN118313315A (zh) * | 2024-06-05 | 2024-07-09 | 中国人民解放军国防科技大学 | 基于通量重构与非线性谱差分混合的流场计算方法 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20060089803A1 (en) * | 2002-12-27 | 2006-04-27 | Riken | Method and device for numberical analysis of flow field of non-compressive viscous fluid, directly using v-cad data |
CN107391885A (zh) * | 2017-08-29 | 2017-11-24 | 西北工业大学 | 基于有限体积法的剪切滑移动网格方法 |
CN108197367A (zh) * | 2017-12-27 | 2018-06-22 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种基于流场通量阶跃的高精度间断Galerkin人工粘性激波捕捉方法 |
CN110457806A (zh) * | 2019-08-02 | 2019-11-15 | 南京航空航天大学 | 基于交错网格的中心五阶weno格式的全流场模拟方法 |
-
2022
- 2022-01-25 CN CN202210086558.0A patent/CN114638173B/zh active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20060089803A1 (en) * | 2002-12-27 | 2006-04-27 | Riken | Method and device for numberical analysis of flow field of non-compressive viscous fluid, directly using v-cad data |
CN107391885A (zh) * | 2017-08-29 | 2017-11-24 | 西北工业大学 | 基于有限体积法的剪切滑移动网格方法 |
CN108197367A (zh) * | 2017-12-27 | 2018-06-22 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种基于流场通量阶跃的高精度间断Galerkin人工粘性激波捕捉方法 |
CN110457806A (zh) * | 2019-08-02 | 2019-11-15 | 南京航空航天大学 | 基于交错网格的中心五阶weno格式的全流场模拟方法 |
Cited By (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116384288A (zh) * | 2023-06-05 | 2023-07-04 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种可压缩流动高分辨率数值模拟方法、介质及设备 |
CN116384288B (zh) * | 2023-06-05 | 2023-08-25 | 中国空气动力研究与发展中心计算空气动力研究所 | 一种可压缩流动高分辨率数值模拟方法、介质及设备 |
CN116542184A (zh) * | 2023-07-05 | 2023-08-04 | 中国空气动力研究与发展中心计算空气动力研究所 | 粘性通量的计算方法、装置、终端设备和存储介质 |
CN116542184B (zh) * | 2023-07-05 | 2023-09-19 | 中国空气动力研究与发展中心计算空气动力研究所 | 粘性通量的计算方法、装置、终端设备和存储介质 |
CN116611369A (zh) * | 2023-07-17 | 2023-08-18 | 中国空气动力研究与发展中心计算空气动力研究所 | 基于光滑度量量级与候选模板点数的插值方法及装置 |
CN116611369B (zh) * | 2023-07-17 | 2023-09-29 | 中国空气动力研究与发展中心计算空气动力研究所 | 基于光滑度量量级与候选模板点数的插值方法及装置 |
CN118313313A (zh) * | 2024-06-04 | 2024-07-09 | 中国人民解放军国防科技大学 | 基于非线性谱差分的高阶精度流场计算方法 |
CN118313315A (zh) * | 2024-06-05 | 2024-07-09 | 中国人民解放军国防科技大学 | 基于通量重构与非线性谱差分混合的流场计算方法 |
CN118313315B (zh) * | 2024-06-05 | 2024-08-06 | 中国人民解放军国防科技大学 | 基于通量重构与非线性谱差分混合的流场计算方法 |
Also Published As
Publication number | Publication date |
---|---|
CN114638173B (zh) | 2023-06-02 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN114638173A (zh) | 一种高阶非线性激波捕捉空间离散方法 | |
Fidkowski et al. | Review of output-based error estimation and mesh adaptation in computational fluid dynamics | |
Ceze et al. | Drag prediction using adaptive discontinuous finite elements | |
Biancolini et al. | Static aeroelastic analysis of an aircraft wind-tunnel model by means of modal RBF mesh updating | |
CN112214869B (zh) | 一种求解欧拉方程的改进型高阶非线性空间离散方法 | |
Pandya et al. | Assessment of USM3D Hierarchical Adaptive Nonlinear Method Preconditioners for Three-Dimensional Cases | |
CN113158338B (zh) | 一种基于粗网格快速湍流壁面函数气动力预测方法 | |
CN115238397B (zh) | 高超声速飞行器热环境计算方法、装置和计算机设备 | |
CN113723027A (zh) | 一种针对弹性飞机的静气动弹性计算方法 | |
Luo et al. | A fast, matrix-free implicit method for computing low Mach number flows on unstructured grids | |
Diskin et al. | Grid convergence for turbulent flows | |
Diskin et al. | Grid-convergence of Reynolds-averaged Navier–Stokes solutions for benchmark flows in two dimensions | |
CN112580241A (zh) | 一种基于结构降阶模型的非线性气动弹性动响应分析方法 | |
CN114611438A (zh) | 一种湍流流动中目标物的受力状态模拟方法及装置 | |
Fernández | Entropy-stable hybridized discontinuous Galerkin methods for large-eddy simulation of transitional and turbulent flows | |
Huang et al. | Piston problem for the isentropic Euler equations for a modified Chaplygin gas | |
CN113361173B (zh) | 一种完全基于当地流场参数对转捩模型可压缩修正的方法 | |
Diskin et al. | Verification test suite for spalart-allmaras QCR2000 turbulence model | |
Golovatov et al. | Optimization of technological parameters of impregnation of load-bearing rod elements of reflector made of polymer composite materials by transfer molding method | |
Kaya et al. | Discrete adjoint-based aerodynamic shape optimization framework for natural laminar flows | |
Du et al. | Super resolution generative adversarial networks for multi-fidelity pressure distribution prediction | |
Eliasson et al. | Application of a line-implicit scheme on stretched unstructured grids | |
Dulikravich et al. | Multidisciplinary inverse problems | |
Das et al. | 3D topology optimization of aircraft wings with conventional and non-conventional layouts: A comparative study | |
Diskin et al. | Reference solutions for benchmark three dimensional turbulent flows |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |