CN111241728A - 一种欧拉方程的间断伽辽金有限元数值求解方法 - Google Patents
一种欧拉方程的间断伽辽金有限元数值求解方法 Download PDFInfo
- Publication number
- CN111241728A CN111241728A CN202010005524.5A CN202010005524A CN111241728A CN 111241728 A CN111241728 A CN 111241728A CN 202010005524 A CN202010005524 A CN 202010005524A CN 111241728 A CN111241728 A CN 111241728A
- Authority
- CN
- China
- Prior art keywords
- matrix
- pod
- error
- finite element
- dimension
- 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
Images
Classifications
-
- 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
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Complex Calculations (AREA)
Abstract
本发明属于计算流体力学数值求解领域,涉及一种欧拉方程的间断伽辽金有限元数值求解方法,基于正交分解模型降阶。本发明将间断伽辽金有限元(DG)方法与正交分解(POD)方法结合,通过构造POD基来对求解变量做近似,以此来降低模型维数,利用POD模型降阶把多维的物理过程进行低维的近似描述,降低模型维数以使得DG算法的计算过程中,在维护可接受的精度的同时使得计算模型的自由度减少以降低计算成本,减少时间和内存消耗。最终本发明既保持了有限元和有限体积的优点又克服了它们的不足。
Description
技术领域
本发明属于计算流体力学数值求解领域,涉及一种欧拉方程的间断伽辽金有限元数值求解方法,基于正交分解模型降阶。
背景技术
计算流体力学(Computational fluid dynamics,CFD)是现代流体力学的一个重要学科分支,它是计算机科学、流体力学和计算数学以及计算物理相结合的产物。在用CFD解决实际问题的过程中,首先要将实际问题描述成相应的数学模型,一般情况下是给定初值及边界条件的偏微分方程系统,常解的模型为Euler/Navier-Stokes(N-S)方程。在CFD领域,数值方法已经是其求解计算最重要的方式之一,包括有限元法、有限差分法、有限体积法或间断有限元方法等。
间断有限元方法,也称间断伽辽金(Discontinuous Galerkin,DG)方法,是利用完全间断的分片多项式空间作为近似解和试验函数空间的有限元方法,目前该方法已经成熟应用于计算流体力学中。DG方法不仅具有局部守恒、稳定等特性,并且可以通过提高插值函数的阶数来提高精度,所以容易处理复杂的几何形状以及有悬挂节点的不规则网格,并且在不同单元里具有不同度的逼近多项式。这样间断有限元方法既保持了有限元和有限体积的优点又克服了它们的不足。
尽管DG方法有诸多优点,但其有个主要的缺点,尤其对于稳定问题特别敏感,即DG方法在每个单元需要求解的变量数增加且随着精度的提高变量数是非线性的增长,对系统做数值模拟时会导出一个很大的偏微分方程组,导致计算时的自由度太多,需要大量计算时间和内存容量,使得对真实时间估计不可控制,因此极大的限制了其实际应用。而目前在计算流体力学数值求解领域并没有关于优化DG方法缺陷的技术手段。
发明内容
针对上述存在问题或不足,为解决DG方法在求解Euler方程过程中自由度太多而造成的大量计算时间和内存消耗问题,提高计算效率,本发明提供了一种欧拉方程的间断伽辽金有限元数值求解方法,即基于正交分解(Proper Orthogonal Decomposition,POD)模型降阶的间断伽辽金(DG)有限元数值方法。
具体技术方案包括以下步骤:
A、将目标Euler方程用间断Galerkin有限元算法进行离散,建立一个广义系统;
B、根据步骤A建立的广义系统求解出瞬时解的集合U,作为样本进行选取来构造瞬像矩阵A;
C、对步骤B中得到的瞬像矩阵A进行奇异值分解得到互相关矩阵的特征系统C=ATA;
D、求出步骤C中所得特征系统C=ATA的特征值和对应的特征向量,根据误差界从中选择POD基向量;
F、对步骤E带入近似的状态向量后的广义系统做Galerkin投影,得到降阶模型并求解。
本发明将DG方法与POD方法结合,通过构造POD基来对求解变量做近似,以此来降低模型维数,利用POD模型降阶把多维的物理过程进行低维的近似描述,降低模型维数以使得DG算法的计算过程中,在维护可接受的精度的同时使得计算模型的自由度减少以降低计算成本,减少时间和内存消耗。最终本发明既保持了有限元和有限体积的优点又克服了它们的不足。
附图说明
图1是本发明的流程图。
图2是实施例步骤B-D求POD基的流程示意图。
图3为实施例三维球DG与POD降维后的压力场图。
图4为实施例三维球降维前后表面场压力系数变化对比。
图5为实施例三维球降维前后表面场相对误差变化。
具体实施方式
下面结合附图来详细说明本发明的技术方案。
本发明属于计算流体力学数值求解领域,涉及一种欧拉方程的间断伽辽金有限元数值求解方法,基于正交分解模型降阶,参照图1,包括以下步骤:
A.将目标Euler方程用间断Galerkin有限元算法进行离散,建立一个广义系统;
三维非定常可压缩欧拉/纳维-斯托克斯(Euler/Navier-Stokes)方程组的守恒形式为:
当ivis=0时为Euler方程,ivis=1时为N-S方程。
式中ρ,p分别为密度和压力;x,y,z分别为笛卡尔坐标系下的坐标分量;u,v,w分别为笛卡尔坐标系下的速度分量;e,h分别为内能和焓;E,H分别为总能和总焓。
对于完全气体,有:
其中R表示气体常数,对于空气,一般取287.053焦耳/千克开。cp,cv,γ分别为定压比热、定容比热和比热比。
用间断Galerkin有限元算法进行离散,化简后得如下的广义系统:
M为质量矩阵,其元素为mij,只与单元坐标及单元类型相关;ui为单元自由度,通过时间推进求解,R(u)为右端项,u=(ρ,ρu,ρv,ρw,ρE)T,t∈[0,TF],TF为总的求解时间。
B.根据步骤A建立的广义系统求解出瞬时解的集合U,作为样本进行选取来构造瞬像矩阵A;
瞬像矩阵的构造对于求解POD基很关键,一般在一个很小的时间区间[0,T0](T0远小于TF)上进行选取,将这个区间不同时刻的瞬时解当做样本。现将时间区间[0,T0]上分成N个相等子区间,即时间步长则可得到样本矩阵:
其中n为单元个数,ui(tnl)(i=1,2,...,4n,l=1,2,...,N)为第i个点上第nl时刻的解。接下来从样本矩阵里选取L个时间点的解构造瞬像矩阵如下:
对于L的确定,我们可以通过样本数据的后验误差来进行判断:
L确定的具体步骤如下:
②重复①中的步骤求出L0+1对应的误差error(1),若error(1)≤error(0),则循环步骤②直至error(k)>error(k-1)时停止计算,L0+k-1即为所求的L。
C.对步骤B中得到的瞬像矩阵A进行奇异值分解得到互相关矩阵的特征系统C=ATA,具体如下:
对瞬像矩阵做奇异值(SVD)分解
r为瞬像矩阵的秩,Σr×r=diag(σ1,σ2,...,σr),其中σi(i=1,2,...,r)为A按递减顺序排列的奇异值,且有σ1≥σ2≥...≥σr>0。令U=(φ1,φ2,...,φN),分别为N×N和L×L的酉矩阵,其中φi(i=1,2,...,N)为AAT的特征向量;同样的为ATA的特征向量。此时POD基可由降维的特征系统C=ATA得到。
D.求出步骤C中所得特征系统C=ATA的特征值和对应的特征向量,根据误差界从中选择POD基向量;
或者定义为
其中(an)i为特征向量里的元素,Un为A的列向量,式(7)与式(8)完全等价。式(7)与式(8)此时的维数为特征系统C=ATA的秩,但并非最佳POD基的维数,为估计POD基的维数d,定义如下误差界:
通过(9)式,得到POD基的维数d的估计值,则取前d项构成POD基。
在由步骤D获得POD基向量之后,我们扩展任意几何形状的流动解。为了说明情况,三维流场的自由度将被表示为:
统一上述守恒变量系数,作如下近似:
记自由度为d=dρ+dρu+dρv+dρw+dρE,其中ψ=(ψρ,ψρu,ψρv,ψρw,ψρE)为4n×d的矩阵,α(t)为d×5的矩阵:
F.对步骤E带入近似的状态向量后广义系统做Galerkin投影,得到降阶模型并求解;最后,对步骤E中的广义系统,即(11)式做Galerkin投影,则可以得到如下的降维模型:
由于d<<4n,求解系统的迭代矩阵从原来4n×4n的M矩阵降维到维数只有d×d的Mn矩阵,右端项从原来4n×5的矩阵降维到维数只有d×5的Rn(u)矩阵,所需求解的自由度,也由原来的4n×5个降到d个,从而大大降低了求解系统的维数,减少计算时间,提高计算效率。
实施例:
以一个半径为0.5的球体以及[-20,20]×[-20,20]×[-20,20]的空气盒子为模型验证POD模型降阶的有效性,考虑一个亚声速绕流问题,其马赫数为M∞=0.2,在单元数为7893的简单球体绕流模型上进行计算。其场图的结果对比如图3所示:
图3-a为原DG的压力场图,图3-b为POD降维后的压力场图,由图可知两者得到的场图结果一致;图4为POD降维前后表面场压力系数变化的对比,图5为POD降维前后表面场的相对误差变化,由图可知降维前后相对误差控制在0.2%左右,并不影响计算精度。
表1则给出了计算结果
表1三维球降维前后计算结果对比
由表中可见,降维前算法自由度为31572个,POD降维后5个守恒变量的自由度在每个点上都相等,分别为1,2,17,17,1,总自由度仅38个,降维前RKDG算法时间步长为10-5,POD降维后为5*10-5,时间步长扩大5倍;而整体计算时间提高了2.5倍。采用该实施例,进一步说明,本发明的方法减少了计算时间,提高了计算效率。
Claims (4)
1.一种欧拉方程的间断伽辽金有限元数值求解方法,包括以下步骤:
A、将目标欧拉Euler方程用间断伽辽金Galerkin有限元算法进行离散,建立一个广义系统;
M为质量矩阵,其元素为mij,只与单元坐标及单元类型相关;ui为单元自由度,通过时间推进求解,R(u)为右端项,u=(ρ,ρu,ρv,ρw,ρE)T,t∈[0,TF],TF为总的求解时间;
B、根据步骤A建立的广义系统求解出瞬时解的集合U,作为样本进行选取来构造瞬像矩阵A;
其中ui(tnl)为第i个点上第nl时刻的解,n为单元个数,i=1,2,...,4n;l=1,2,...,N;接下来从样本矩阵里选取L个时间点的解构造瞬像矩阵如下:
对于L的确定,通过样本数据的后验误差来进行判断:
L确定的具体步骤如下:
②重复①中的步骤求出L0+1对应的误差error(1),若error(1)≤error(0),则循环步骤②直至error(k)>error(k-1)时停止计算,L0+k-1即为所求的L;
C、对步骤B中得到的瞬像矩阵A进行奇异值分解得到互相关矩阵的特征系统C=ATA;
D、求出步骤C中所得特征系统C=ATA的特征值和对应的特征向量,根据误差界从中选择POD基向量;
F、对步骤E带入近似状态向量后的广义系统做Galerkin投影,得到降阶模型并求解。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010005524.5A CN111241728B (zh) | 2020-01-03 | 2020-01-03 | 一种欧拉方程的间断伽辽金有限元数值求解方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202010005524.5A CN111241728B (zh) | 2020-01-03 | 2020-01-03 | 一种欧拉方程的间断伽辽金有限元数值求解方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN111241728A true CN111241728A (zh) | 2020-06-05 |
CN111241728B CN111241728B (zh) | 2023-05-05 |
Family
ID=70868666
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202010005524.5A Active CN111241728B (zh) | 2020-01-03 | 2020-01-03 | 一种欧拉方程的间断伽辽金有限元数值求解方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN111241728B (zh) |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116504341A (zh) * | 2022-05-20 | 2023-07-28 | 大连理工大学 | 数据驱动识别偏微分方程的序列奇异值过滤方法 |
CN117194859A (zh) * | 2023-08-23 | 2023-12-08 | 哈尔滨工程大学 | 基于间断伽辽金法的非结构网格自适应细分与高效并行高精度算法框架的构建方法及系统 |
CN117236104A (zh) * | 2023-08-22 | 2023-12-15 | 西南科技大学 | 一种基于环境模拟时序仿真的特种阀高阶模型降阶方法 |
Citations (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105389839A (zh) * | 2015-11-06 | 2016-03-09 | 北京航空航天大学 | 基于流体分析的流体参数估计方法 |
GB201703586D0 (en) * | 2017-03-07 | 2017-04-19 | Kompetenzzentrum-Das Virtuelle Fahrzeug | Method to reduce complexity in flow noise calculation |
CN107944141A (zh) * | 2017-11-24 | 2018-04-20 | 电子科技大学 | 基于杂交时域间断伽辽金法的时域计算电磁学数值方法 |
CN108108332A (zh) * | 2018-01-26 | 2018-06-01 | 浙江大学 | 一种基于广义流体力学非线性本构方程的耦合求解方法 |
CN108197358A (zh) * | 2017-12-20 | 2018-06-22 | 北京石油化工学院 | 一种高效快速模拟水力压裂的方法 |
CN108536954A (zh) * | 2018-04-08 | 2018-09-14 | 南京航空航天大学 | 一种基于交点间断伽辽金的高精度格子波尔兹曼方法 |
CN109190169A (zh) * | 2018-08-02 | 2019-01-11 | 电子科技大学 | 一种三维时域电磁学杂交时域间断伽辽金数值方法 |
CN109858158A (zh) * | 2019-02-01 | 2019-06-07 | 中国人民解放军军事科学院国防科技创新研究院 | 一种计算流体力学模拟的参数配置方法及系统 |
CN110516316A (zh) * | 2019-08-03 | 2019-11-29 | 电子科技大学 | 一种间断伽辽金法求解欧拉方程的gpu加速方法 |
-
2020
- 2020-01-03 CN CN202010005524.5A patent/CN111241728B/zh active Active
Patent Citations (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN105389839A (zh) * | 2015-11-06 | 2016-03-09 | 北京航空航天大学 | 基于流体分析的流体参数估计方法 |
GB201703586D0 (en) * | 2017-03-07 | 2017-04-19 | Kompetenzzentrum-Das Virtuelle Fahrzeug | Method to reduce complexity in flow noise calculation |
CN107944141A (zh) * | 2017-11-24 | 2018-04-20 | 电子科技大学 | 基于杂交时域间断伽辽金法的时域计算电磁学数值方法 |
CN108197358A (zh) * | 2017-12-20 | 2018-06-22 | 北京石油化工学院 | 一种高效快速模拟水力压裂的方法 |
CN108108332A (zh) * | 2018-01-26 | 2018-06-01 | 浙江大学 | 一种基于广义流体力学非线性本构方程的耦合求解方法 |
CN108536954A (zh) * | 2018-04-08 | 2018-09-14 | 南京航空航天大学 | 一种基于交点间断伽辽金的高精度格子波尔兹曼方法 |
CN109190169A (zh) * | 2018-08-02 | 2019-01-11 | 电子科技大学 | 一种三维时域电磁学杂交时域间断伽辽金数值方法 |
CN109858158A (zh) * | 2019-02-01 | 2019-06-07 | 中国人民解放军军事科学院国防科技创新研究院 | 一种计算流体力学模拟的参数配置方法及系统 |
CN110516316A (zh) * | 2019-08-03 | 2019-11-29 | 电子科技大学 | 一种间断伽辽金法求解欧拉方程的gpu加速方法 |
Non-Patent Citations (3)
Title |
---|
FAHEEM AHMED等: "Numerical solution of 2D Euler equations for the transonic flow past over NACA0012 and RAE2822 airfoils using high order accurate Runge-Kutta Discontinuous Galerkin method" * |
朱兰: "基于HDG的混合网格法的应用研究" * |
梁俊凯: "欧拉方程间断伽辽金有限元解法的研究" * |
Cited By (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116504341A (zh) * | 2022-05-20 | 2023-07-28 | 大连理工大学 | 数据驱动识别偏微分方程的序列奇异值过滤方法 |
CN116504341B (zh) * | 2022-05-20 | 2023-11-07 | 大连理工大学 | 数据驱动识别偏微分方程的序列奇异值过滤方法 |
CN117236104A (zh) * | 2023-08-22 | 2023-12-15 | 西南科技大学 | 一种基于环境模拟时序仿真的特种阀高阶模型降阶方法 |
CN117236104B (zh) * | 2023-08-22 | 2024-04-02 | 西南科技大学 | 一种基于环境模拟时序仿真的特种阀高阶模型降阶方法 |
CN117194859A (zh) * | 2023-08-23 | 2023-12-08 | 哈尔滨工程大学 | 基于间断伽辽金法的非结构网格自适应细分与高效并行高精度算法框架的构建方法及系统 |
CN117194859B (zh) * | 2023-08-23 | 2024-09-17 | 哈尔滨工程大学 | 基于间断伽辽金法的非结构网格自适应细分与高效并行高精度算法框架的构建方法及系统 |
Also Published As
Publication number | Publication date |
---|---|
CN111241728B (zh) | 2023-05-05 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN111241728A (zh) | 一种欧拉方程的间断伽辽金有限元数值求解方法 | |
CN106767780B (zh) | 基于Chebyshev多项式插值逼近的扩展椭球集员滤波方法 | |
CN108345741B (zh) | 基于无网格rkpm各向异性材料二维热变形和热应力分析方法 | |
CN110781626A (zh) | 有限差分多重分辨三角函数weno格式的模拟方法 | |
CN110336057A (zh) | 一种交叉流电堆二维温度分布观测器构建方法及其应用 | |
CN116205153A (zh) | 一种基于流场物理信息的网格密度优化方法 | |
Alekseenko et al. | Deterministic solution of the Boltzmann equation using a discontinuous Galerkin velocity discretization | |
CN118070621A (zh) | 固壁边界的处理方法、装置、终端设备及存储介质 | |
Du et al. | Computing critical eigenvalues of power systems using inexact two-sided Jacobi-Davidson | |
CN105160092B (zh) | 一种适用于热防护系统瞬态温度场计算的热环境插值方法 | |
Lefebvre et al. | Least squares finite element solution of compressible and incompressible flows | |
CN114818550B (zh) | 一种飞机振动试验中时变气动载荷地面等效模拟方法 | |
Onur et al. | Effects of the Jacobian evaluation on Newton's solution of the Euler equations | |
CN114036806B (zh) | 基于热导率各向异性介质的三维地温场数值模拟方法 | |
Li et al. | The adaptive GRP scheme for compressible fluid flows over unstructured meshes | |
Winter et al. | Nonlinear reduced-order modeling of unsteady aerodynamic loads based on dynamic local linear neuro-fuzzy models | |
Balan et al. | A stable spectral difference method for triangles | |
CN112307684B (zh) | 一种多重分辨率weno格式结合ilw边界处理的定点快速扫描方法 | |
CN115995277B (zh) | 一种材料动力学特性评估方法、装置、设备及介质 | |
CN113792461B (zh) | 一种工程结构在极端荷载下动态响应的复合时域分析方法 | |
CN118378554B (zh) | 基于应力约束的多自由度柔顺机构稳健性拓扑优化方法 | |
DOLEJŠÍ et al. | Anisotropic mesh adaptation for transonic and supersonic flow simulation | |
CN116305547A (zh) | 一种改善高超声速数值激波稳定性的计算方法和装置 | |
Dao et al. | Comparison of upwind and centered schemes for low Mach number flows | |
Li et al. | Unified Solution of Conjugate Fluid and Solid Heat Transfer–Part I. Solid Heat Conduction |
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 |