CN109726465B - 基于非结构曲边网格的三维无粘低速绕流的数值模拟方法 - Google Patents
基于非结构曲边网格的三维无粘低速绕流的数值模拟方法 Download PDFInfo
- Publication number
- CN109726465B CN109726465B CN201811600701.3A CN201811600701A CN109726465B CN 109726465 B CN109726465 B CN 109726465B CN 201811600701 A CN201811600701 A CN 201811600701A CN 109726465 B CN109726465 B CN 109726465B
- Authority
- CN
- China
- Prior art keywords
- coordinate system
- tetrahedral mesh
- curved edge
- equation
- reference coordinate
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 39
- 238000004088 simulation Methods 0.000 title claims abstract description 19
- 239000000853 adhesive Substances 0.000 title claims abstract description 13
- 230000001070 adhesive effect Effects 0.000 title claims abstract description 13
- 239000012530 fluid Substances 0.000 claims abstract description 17
- 238000004364 calculation method Methods 0.000 claims description 16
- 239000011159 matrix material Substances 0.000 claims description 15
- 238000009826 distribution Methods 0.000 claims description 14
- 239000006185 dispersion Substances 0.000 claims description 9
- 230000004069 differentiation Effects 0.000 claims description 6
- 230000004907 flux Effects 0.000 claims description 3
- 230000010354 integration Effects 0.000 claims description 3
- 230000009466 transformation Effects 0.000 claims description 3
- 230000001131 transforming effect Effects 0.000 claims description 3
- 238000004422 calculation algorithm Methods 0.000 description 6
- 239000007787 solid Substances 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000007613 environmental effect Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 230000008569 process Effects 0.000 description 1
Images
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Complex Calculations (AREA)
Abstract
本发明属于三维流体力学数值求解技术领域,涉及一种基于非结构曲边网格的三维无粘低速绕流的数值模拟方法。本发明通过非结构曲边网格来精确拟合曲面物面,然后把非结构曲边网格投影到参考坐标系获得直网格,在此基础上,计算固壁边界条件;同时开发相应的曲边网格数值模拟方法代替传统的直网格数值模拟方法,提高数值模拟方法的精度。
Description
技术领域
本发明属于三维流体力学数值求解技术领域,涉及一种基于非结构曲边网格的三维无粘低速绕流的数值模拟方法。
背景技术
计算流体力学(简称CFD)一直以来在汽车制造、土木工程、环境工程、船舶工业以及航空工业等领域有着广泛的应用,其有助于解释、理解理论和实验的结果,是流体力学分析不可缺少的方法。前期由于计算机水平的限制,CFD的实际求解仅限于二维流动,而大部分的真实流动世界主要是三维的。随着计算机水平的发展,今天的CFD已经可以大量求解三维流场,虽然仍需要大量的人力和计算机资源,但是这种求解方法在工业设备中已得到广泛使用。
随着CFD的发展,相应的数值算法在CFD中的应用也随之发展起来,如有限差分、有限体积和有限元法等。随着工业技术的进步,流体动力学对数值算法的精度提出了更高的要求,因此需要高精度的数值模拟方法。间断Galerkin有限元法由于其易于实现高阶精度、灵活处理间断问题、适用于非结构网格、有利于实现并行算法,因此有很好的应用前景和工程实用价值。然而间断Galerkin有限元法的高阶精度实现依赖于边界的精度,这包括算法的精度和模型空间离散的精度,因此需要研究物面边界的空间离散,特别是复杂结构的边界拟合。虽然非结构网格的应用能够更好的拟合复杂结构边界,但是非结构网格也难以精确拟合曲面物面。模型空间离散的误差限制了边界的精度,从而制约了高精度算法的应用。
发明内容
针对上述存在问题或不足,为解决空间离散误差导致数值算法精度低的问题;本发明提供了一种基于非结构曲边网格的三维无粘低速绕流的数值模拟方法,以非结构曲边网格来精确拟合曲面物面,在此基础上开发相应的高精度数值模拟方法。
一种基于非结构曲边网格的三维无粘低速绕流的数值模拟方法,包括以下步骤:
A.将目标结构进行建模,然后建立流体计算域;
B.对步骤A所建流体计算域采用曲边四面体网格进行剖分,转化为离散空间模型;
C.利用二阶拉格朗日节点基函数将步骤B获得的曲边四面体网格变换成参考坐标系下的直四面体网格,获得雅克比矩阵;
如图4所示,曲边四面体网格是在真实的直角坐标系(x,y,z)下的表示,其中点1、2、3、4为曲边四面体网格的四个顶点,点5、6、7、8、9、10为曲边四面体网格每条边的中点,这些点的坐标信息xi,yi,zi(i=1,2,…,10)可以通过一般网格生成程序获得,这里不再详细阐述。
如图5所示,直四面体网格是在参考坐标系(ξ,η,ζ)下的表示的标准四面体网格,顶点1的坐标为(0,0,0),顶点2的坐标为(1,0,0),顶点3的坐标为(0,1,0),顶点4的坐标为(0,0,1),点5、6、7、8、9、10为对应每条边的中点。在参考坐标系下我们定义如下二阶拉格朗日节点基函数:
其中ξ,η,ζ为参考坐标系分量。真实坐标系的分量x,y,z可以表示为:
最后直角坐标系和参考坐标系之间的变换雅克比矩阵为:
将(1)式、(2)式代入(3)式可最终求得雅克比矩阵的显示表达式,雅克比矩阵的逆表示为:
D.在参考坐标系下计算目标结构固壁边界网格的固壁边界条件;
在直角坐标系下,曲边四面体网格的边界面上的值有如下固壁边界条件:
其中u,v,w分别为直角坐标系下的速度分量;下标l,r分别表示本网格和相邻网格;表示直角坐标系的基矢量;表示曲边四面体网格中边界面上外单位法向量;ρ,p分别为密度和压力。由于在边界上,曲边四面体网格的边界面没有相邻网格,因此我们通过(5)式来构造相邻网格在边界面上的值,即为固壁边界条件。曲边四面体网格中,由于其边界面为曲面因此无法直接计算(5)式,需要在参考坐标系下的直四面体网格中计算。
基于参考坐标系,有:
其中:
将(6)式,(8)式代入(5)式。最后,参考坐标系下直四面体网格中的边界面上的值有如下固壁边界条件:
E.利用间断Galerkin有限元法,将三维无粘低速绕流的控制方程在步骤C所得每一个直四面体网格上进行空间离散,获得一个关于时间微分的有限元方程;
对于无粘绕流问题,在参考坐标系下求解如下三维守恒形式的欧拉方程:
式中Q为守恒变量,D=[Fc,Gc,Hc]为无粘通量张量,其具体形式如下:
其中u,v,w分别为直角坐标系下的速度分量;ρ,p分别为密度和压力;e为总能;U,V,W为逆变速度。(11)式是5个方程的组合,为了表示方便我们用下标h(h=1,2,3,4,5)表示(11)式中第h几个方程以及Q和D中的第h个分量Qh,Dh。
在参考坐标系下,对于间断Galerkin有限元法,变量在直四面体网格Ω内的分布采用下面的多项式近似表达:
φj(ξ,η,ζ)表示插值基函数,在这里我们选择正交基函数,N表示插值基函数的个数。(11)式两端乘以测试函数φi(ξ,η,ζ)(i=0,…N),然后在Ω内积分并且把(14)式代入,可以得到(11)式中第h个方程的Galerkin方法弱形式:
上述积分项采用数值积分计算,(15)式最后简化为:
其中q=[q1 q2 q3 q4 q5]T,RHS=[RHS1 RHS2 RHS3 RHS4 RHS5]T
F.对步骤E所得有限元方程进行时间离散,获得迭代方程;
G.对步骤F获得的迭代方程,给定每一个由曲边四面体网格剖分后得到的四面体单元初值,进行循环迭代,直至满足迭代终止条件,获得整个计算域的场分布。
本发明通过非结构曲边网格来精确拟合曲面物面,然后把非结构曲边网格投影到参考坐标系获得直网格,在此基础上,计算目标结构固壁边界条件;同时开发相应的曲边网格数值模拟方法代替传统的直网格数值模拟方法,提高数值模拟方法的精度。
附图说明
图1是本发明流程图;
图2是实施例的流体计算域模型剖面图;
图3是实施例球表面网格示意图;
图4是实施例直角坐标系下曲边四面体网格示意图;
图5是实施例参考坐标系下直四面体网格示意图;
图6是实施例场分布剖视图;
图7是现有非结构直网格三维无粘低速绕流数值模拟方法场分布剖视图。
具体实施方式
下面结合附图和实施例来详细说明本发明的技术方案。
参照附图1,一种基于非结构曲边网格的三维无粘低速绕流的数值模拟方法,包括以下步骤:
A.建立球型结构的几何模型,然后建立流体计算域,其结构剖视图如图2所示。
B.对步骤A所建流体计算域采用曲边四面体网格进行剖分,转化为离散空间模型;
采用曲边四面体网格剖分步骤A所建流体计算域,剖分后的计算域被人为分割为多个三维曲边四面体网格,从而将连续的几何空间转化为离散的网格空间,球的表面网格如图3所示。
C.利用二阶拉格朗日节点基函数将步骤B获得的曲边四面体网格变换成参考坐标系下的直四面体网格,获得雅克比矩阵;
如图4所示,曲边四面体网格是在真实的直角坐标系(x,y,z)下的表示,其中点1、2、3、4为曲边四面体网格的四个顶点,点5、6、7、8、9、10为曲边四面体网格每条边的中点,这些点的坐标信息xi,yi,zi(i=1,2,…,10)可以通过一般网格生成程序获得,这里不再详细阐述。如图5所示,直四面体网格是在参考坐标系(ξ,η,ζ)下的表示的标准四面体网格,顶点1的坐标为(0,0,0),顶点2的坐标为(1,0,0),顶点3的坐标为(0,1,0),顶点4的坐标为(0,0,1),点5、6、7、8、9、10为对应每条边的中点。在参考坐标系下我们定义如下二阶拉格朗日节点基函数:
其中ξ,η,ζ为参考坐标系分量。真实坐标系的分量x,y,z可以表示为:
最后直角坐标系和参考坐标系之间的变换雅克比矩阵为:
将(1)式、(2)式代入(3)式可最终求得雅克比矩阵的显示表达式,雅克比矩阵的逆表示为:
D.在参考坐标系下计算球形结构固壁边界网格的固壁边界条件;
在直角坐标系下,曲边四面体网格的边界面上的值有如下固壁边界条件:
其中u,v,w分别为直角坐标系下的速度分量;下标l,r分别表示本网格和相邻网格;表示直角坐标系的基矢量;表示曲边四面体网格中边界面上外单位法向量;ρ,p分别为密度和压力。由于在边界上,曲边四面体网格的边界面没有相邻网格,因此我们通过(5)式来构造相邻网格在边界面上的值,即为固壁边界条件。曲边四面体网格中,由于其边界面为曲面因此无法直接计算(5)式,需要在参考坐标系下的直四面体网格中计算。
基于参考坐标系,有:
其中:
将(6)式,(8)式代入(5)式。最后,参考坐标系下直四面体网格中的边界面上的值有如下固壁边界条件:
E.利用间断Galerkin有限元法,将三维无粘低速绕流的控制方程在步骤C所得每一个直四面体网格上进行空间离散,获得一个关于时间微分的有限元方程;
对于无粘绕流问题,在参考坐标系下求解如下三维守恒形式的欧拉方程:
式中Q为守恒变量,D=[Fc,Gc,Hc]为无粘通量张量,其具体形式如下:
其中u,v,w分别为直角坐标系下的速度分量;ρ,p分别为密度和压力;e为总能;U,V,W为逆变速度。(11)式是5个方程的组合,为了表示方便我们用下标h(h=1,2,3,4,5)表示(11)式中第h几个方程以及Q和D中的第h个分量Qh,Dh。
在参考坐标系下,对于间断Galerkin有限元法,变量在直四面体网格Ω内的分布采用下面的多项式近似表达:
φj(ξ,η,ζ)表示插值基函数,在这里我们选择正交基函数,N表示插值基函数的个数。(11)式两端乘以测试函数φi(ξ,η,ζ)(i=0,…N),然后在Ω内积分并且把(14)式代入,可以得到(11)式中第h个方程的Galerkin方法弱形式:
上述积分项采用数值积分计算,(15)式最后简化为:
其中q=[q1 q2 q3 q4 q5]T,RHS=[RHS1 RHS2 RHS3 RHS4 RHS5]T
F.对步骤E所得有限元方程进行时间离散,获得迭代方程;
时间离散上采用二阶龙格库塔方法,显示二阶龙格库塔方法如下:
其中k表示时间步。上述方程是一个随时间步迭代的方程,可以通过前一个时刻k的值计算下一个时刻k+1的值。
G.对步骤F获得的迭代方程,给定每一个由曲边四面体网格剖分后得到的四面体单元初值,进行循环迭代,直至满足迭代终止条件,获得整个计算域的场分布。
根据实际问题给定所有四面体单元的初值,按照(21)式计算所有四面体单元当前时刻的值,然后再把当前时刻的值当成初值,继续计算下一时刻的值,以此反复循环迭代,直至计算结果收敛。最后根据(14)式计算每个单元的场分布,给出整个计算域的场分布。
图2示出了实施例的流体计算域模型剖视图;图6示出了实施例场分布剖视图;图7示出了现有非结构直网格三维无粘低速绕流数值模拟方法场分布剖视图。对比图6、图7的场分布可以看出实施例的场分布更加具有对称性,这更加接近真实情况,从而说明本发明相比现有技术提高了数值模拟方法的精度。
Claims (1)
1.一种基于非结构曲边网格的三维无粘低速绕流的数值模拟方法,包括以下步骤:
A.将目标结构进行建模,然后建立流体计算域;
B.对步骤A所建流体计算域采用曲边四面体网格进行剖分,转化为离散空间模型;
C.利用二阶拉格朗日节点基函数将步骤B获得的曲边四面体网格变换成参考坐标系下的直四面体网格,获得雅克比矩阵;
曲边四面体网格是在真实的直角坐标系(x,y,z)下的表示,其中点1、2、3、4为曲边四面体网格的四个顶点,点5、6、7、8、9、10为曲边四面体网格每条边的中点,并获得坐标信息xi,yi,zi(i=1,2,…,10);
直四面体网格是在参考坐标系(ξ,η,ζ)下的表示的标准四面体网格,顶点1的坐标为(0,0,0),顶点2的坐标为(1,0,0),顶点3的坐标为(0,1,0),顶点4的坐标为(0,0,1),点5、6、7、8、9、10为对应每条边的中点;在参考坐标系下定义如下二阶拉格朗日节点基函数:
其中ξ,η,ζ为参考坐标系分量,真实坐标系的分量x,y,z可以表示为:
最后直角坐标系和参考坐标系之间的变换雅克比矩阵为:
将(1)式、(2)式代入(3)式可最终求得雅克比矩阵的显示表达式,雅克比矩阵的逆表示为:
D.在参考坐标系下计算目标结构固壁边界网格的固壁边界条件;
在直角坐标系下,曲边四面体网格的边界面上的值有如下固壁边界条件:
通过(5)式来构造相邻网格在边界面上的值,即为固壁边界条件;在参考坐标系下的直四面体网格中计算,基于参考坐标系,有:
其中:
将(6)式,(8)式代入(5)式;参考坐标系下直四面体网格中的边界面上的值有如下固壁边界条件:
E.利用间断Galerkin有限元法,将三维无粘低速绕流的控制方程在步骤C所得每一个直四面体网格上进行空间离散,获得一个关于时间微分的有限元方程;
对于无粘绕流问题,在参考坐标系下求解如下三维守恒形式的欧拉方程:
式中Q为守恒变量,D=[Fc,Gc,Hc]为无粘通量张量,其具体形式如下:
其中u,v,w分别为直角坐标系下的速度分量;ρ,p分别为密度和压力;e为总能;U,V,W为逆变速度;(11)式是5个方程的组合,用下标h(h=1,2,3,4,5)表示(11)式中第h几个方程以及Q和D中的第h个分量Qh,Dh;
在参考坐标系下,对于间断Galerkin有限元法,变量在直四面体网格Ω内的分布采用下面的多项式近似表达:
φj(ξ,η,ζ)表示插值基函数,在这里我们选择正交基函数,N表示插值基函数的个数;(11)式两端乘以测试函数φi(ξ,η,ζ)(i=0,…N),然后在Ω内积分并且把(14)式代入,可以得到(11)式中第h个方程的Galerkin方法弱形式:
上述积分项采用数值积分计算,(15)式最后简化为:
其中q=[q1 q2 q3 q4 q5]T,RHS=[RHS1 RHS2 RHS3 RHS4 RHS5]T
F.对步骤E所得有限元方程进行时间离散,获得迭代方程;
G.对步骤F获得的迭代方程,给定每一个由曲边四面体网格剖分后得到的四面体单元初值,进行循环迭代,直至满足迭代终止条件,获得整个计算域的场分布。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811600701.3A CN109726465B (zh) | 2018-12-26 | 2018-12-26 | 基于非结构曲边网格的三维无粘低速绕流的数值模拟方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811600701.3A CN109726465B (zh) | 2018-12-26 | 2018-12-26 | 基于非结构曲边网格的三维无粘低速绕流的数值模拟方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN109726465A CN109726465A (zh) | 2019-05-07 |
CN109726465B true CN109726465B (zh) | 2022-07-29 |
Family
ID=66297164
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811600701.3A Active CN109726465B (zh) | 2018-12-26 | 2018-12-26 | 基于非结构曲边网格的三维无粘低速绕流的数值模拟方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN109726465B (zh) |
Families Citing this family (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN110909511B (zh) * | 2019-11-22 | 2022-10-14 | 电子科技大学 | 一种无曲面体积分的无粘低速绕流数值模拟方法 |
CN111783343A (zh) * | 2020-07-08 | 2020-10-16 | 中国人民解放军国防科技大学 | 非结构网格方向模板选择方法与流场求解方法 |
CN111859825B (zh) * | 2020-07-28 | 2022-11-22 | 北京航空航天大学 | 包含任意流固界面非定常不可压流动数值模拟方法及设备 |
CN113178011B (zh) * | 2021-04-28 | 2022-08-02 | 哈尔滨工程大学 | 一种用于求解vof对流方程的切割体网格thinc方法 |
CN113505544A (zh) * | 2021-06-18 | 2021-10-15 | 清华大学 | 基于有限体积法自行车运动虚拟数值风洞*** |
CN116151084B (zh) * | 2023-04-21 | 2023-07-14 | 中国空气动力研究与发展中心计算空气动力研究所 | 基于结构网格的模拟方法、装置、终端设备及存储介质 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106055827A (zh) * | 2016-06-15 | 2016-10-26 | 中国石油化工股份有限公司胜利油田分公司勘探开发研究院 | 一种油藏数值模拟参数敏感性分析装置及方法 |
CN107577857A (zh) * | 2017-08-28 | 2018-01-12 | 电子科技大学 | 一种基于热辐射边界条件的三维有限元模拟方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
FI20002296A (fi) * | 2000-10-17 | 2002-04-18 | Lumeo Software Oy | Mekaanisen alisysteemin ja hydraulisen alisysteemin omaavan systeemin simulointi |
-
2018
- 2018-12-26 CN CN201811600701.3A patent/CN109726465B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106055827A (zh) * | 2016-06-15 | 2016-10-26 | 中国石油化工股份有限公司胜利油田分公司勘探开发研究院 | 一种油藏数值模拟参数敏感性分析装置及方法 |
CN107577857A (zh) * | 2017-08-28 | 2018-01-12 | 电子科技大学 | 一种基于热辐射边界条件的三维有限元模拟方法 |
Non-Patent Citations (2)
Title |
---|
尹俊辉 等.微波管结构动力学分析的快速有限元方法及CAD技术研究.《 中国电子学会真空电子学分会第二十一届学术年会》.2018,619-623. * |
张帆.非结构网格有限体积法的空间离散算法研究.《中国优秀博硕士学位论文全文数据库(博士)基础科学辑》.2018,(第07期),A004-10. * |
Also Published As
Publication number | Publication date |
---|---|
CN109726465A (zh) | 2019-05-07 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109726465B (zh) | 基于非结构曲边网格的三维无粘低速绕流的数值模拟方法 | |
CN108763683B (zh) | 一种三角函数框架下新weno格式构造方法 | |
CN112016167B (zh) | 基于仿真和优化耦合的飞行器气动外形设计方法及*** | |
CN109726433B (zh) | 基于曲面边界条件的三维无粘低速绕流的数值模拟方法 | |
CN106682262B (zh) | 一种获取飞行器流场的数值模拟方法 | |
CN111709171A (zh) | 一种热流强耦合问题的等几何求解及散热拓扑生成方法 | |
CN104268322B (zh) | Weno差分方法的一种边界处理技术 | |
CN111859766A (zh) | 可变计算域的拉格朗日积分点有限元数值仿真***及方法 | |
CN104036095A (zh) | 基于区域分解的耦合高精度复杂外形流场快速算法 | |
CN116629079B (zh) | 混合有限元空间构造及求解线弹性力学问题的方法及装置 | |
Gerace et al. | A model-integrated localized collocation meshless method for large scale three-dimensional heat transfer problems | |
Yang et al. | High-order gas-kinetic scheme on three-dimensional unstructured meshes for compressible flows | |
CN111079326A (zh) | 二维各向异性网格单元度量张量场光滑化方法 | |
CN108763777B (zh) | 基于泊松方程显式解的vlsi全局布局模型建立方法 | |
CN116757026B (zh) | 一种基于等几何分析的板架结构极限强度分析实现方法 | |
CN111159956B (zh) | 一种基于特征的流场间断捕捉方法 | |
CN110909511B (zh) | 一种无曲面体积分的无粘低速绕流数值模拟方法 | |
CN105160130B (zh) | 一种基于三维图像的有限差分法预测材料热导率的方法 | |
Yu et al. | Shock capturing for the correction procedure via reconstruction method using artificial viscosity and diffusivity | |
Luo | A finite volume method based on weno reconstruction for compressible flows on hybrid grids | |
CN104992046A (zh) | 流体力学计算***及方法 | |
CN111008492B (zh) | 一种基于无雅克比矩阵的高阶单元欧拉方程数值模拟方法 | |
Huang et al. | A new grid deformation technology with high quality and robustness based on quaternion | |
Duan et al. | High order FR/CPR method for overset meshes | |
CN118052167B (zh) | 一种多维相关响应的流场模型构建方法 |
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 |