公开/公告号CN109948202A
专利类型发明专利
公开/公告日2019-06-28
原文格式PDF
申请/专利权人 浙江远算云计算有限公司;成都远算智能科技有限公司;
申请/专利号CN201910162080.3
发明设计人 吴健明;
申请日2019-03-04
分类号
代理机构浙江翔隆专利事务所(普通合伙);
代理人许守金
地址 311799 浙江省杭州市淳安县千岛湖镇阳光路688号301室
入库时间 2024-02-19 11:32:30
法律状态公告日
法律状态信息
法律状态
2022-10-11
专利权的转移 IPC(主分类):G06F30/28 专利号:ZL2019101620803 登记生效日:20220926 变更事项:专利权人 变更前权利人:浙江远算科技有限公司 变更后权利人:浙江远算科技有限公司 变更事项:地址 变更前权利人:311799 浙江省杭州市淳安县千岛湖镇阳光路688号301室 变更后权利人:311799 浙江省杭州市淳安县千岛湖镇阳光路688号301室 变更事项:专利权人 变更前权利人:成都远算智能科技有限公司 变更后权利人:
专利申请权、专利权的转移
2020-08-14
授权
授权
2019-07-23
实质审查的生效 IPC(主分类):G06F17/50 申请日:20190304
实质审查的生效
2019-06-28
公开
公开
技术领域
本发明涉及基于线化矩阵混合求解模式的高阶CFD隐式时间推进方法,属于计算流体力学技术领域。
背景技术
由于在色散和耗散方面的巨大优势,高阶间断迦辽金方法在计算流体力学(英文缩写:CFD)领域中得到了广泛的研究和应用。考虑到稳定性条件等因素的影响,隐式时间推进方法在高阶间断迦辽金方法的应用最为普遍。而在隐式时间推进的过程中,最核心的内容就是求解每个未知量的线化矩阵。如何能够高效准确的求解线化矩阵,是隐式计算方法中的重要研究内容。
以往求解线化矩阵主要有两种方法,手动解析方法和自动微分方法。
手工解析求出线化矩阵,这种方法是最准确的,然而对于Navier-Stokes方程中的很多部分,如Riemann通量、边界条件、提升算子等部分,很难通过解析的方式求出其关于自由度的导数,而且该方法计算复杂,公式推导十分繁琐,适用性也一般,效率较低,无法快速有效的得到流体仿真的计算结果和流场信息。
自动微分方法依托相应的PDE求解库,针对不同的求导问题,编写相应的代码计算导数项。该方法计算量小,但是代码通用性非常差,效率低。针对控制方程中的不同项,计算代码都需要重新编写,不利于编程实现,从而导致流体仿真的工作量大、效率低,同样无法快速有效的得到流体仿真的计算结果和流场信息。
发明内容
针对现有技术的缺陷,本发明的目的在于提供一种通过全微分和链式法则,将依赖项对自由度的导数直接转换为依赖项对变量的导数与变量对自由度导数的乘积,或者转化为依赖项对变量梯度的导数与变量梯度对自由度的导数乘积;能够大大减少计算量,提高计算效率,利于编程实现,同时兼顾线化矩阵的准确度,能够快速得到流体仿真的结果和流场信息的基于线化矩阵混合求解模式的高阶CFD隐式时间推进方法。
为实现上述目的,本发明的技术方案为:
基于线化矩阵混合求解模式的高阶CFD隐式时间推进方法,其包括如下步骤:
步骤一:建立基于非结构混合网格的高精度间断迦辽金计算框架;用于流体动力学的数值仿真计算,其包括Navier-Stokes控制方程、基函数、空间离散,控制方程采用Navier-Stokes方程;
步骤二:建立以Newton/GMRES为基础的隐式迭代方法;
步骤三:通过全微分和链式法则,将依赖项对自由度的导数直接转换为依赖项对变量的导数与变量对自由度导数的乘积,或者转化为依赖项对变量梯度的导数与变量梯度对自由度的导数乘积;采用手动解析和差分扰动混合的方式求解线化矩阵,带入到控制方程,求解得到仿真计算结果。
为了能够简化线化矩阵的求解,提高高精度间断迦辽金方法隐式迭代的效率,本发明提出了基于线化矩阵混合求解模式的高阶CFD隐式时间推进方法,该方法通过全微分和链式法则,将依赖项对自由度的导数直接转换为依赖项对变量的导数与变量对自由度导数的乘积,或者转化为依赖项对变量梯度的导数与变量梯度对自由度的导数乘积。变量(导数)对其自由度的导数,能够通过解析的方式得到。对于Riemann通量或者边界积分,能够采用差分微小扰动的方式进行求解,相比传统方法,本发明能够大大减少计算量,提高计算效率,也更具有通用性,利于编程实现,同时兼顾线化矩阵的准确度,保证隐式迭代有高效的收敛性,能够快速得到流体仿真的结果和流场信息。
作为优选技术措施,
所述步骤一,其具体包括以下步骤:
步骤101、构建微分形式下的Navier-Stokes方程
上式中,t代表时间,矢量U代表方程求解的守恒变量,
步骤102、将流体的计算区域剖分成若干不重合的离散单元κ,在不同的离散单元内,流场中的守恒变量U能够用基函数的线性组合来表示,选择单项式等级基作为基函数,单元内的积分点根据网格单元的类型来确定;
步骤103、在方程两边同时乘以基函数φi并进行积分,得到弱形式下的控制方程(2);
κ和
其中,
作为优选技术措施,
所述步骤二,其具体包括以下步骤:
步骤201、采用Euler后插方法对时间导数进行离散,从而得到控制方程的全离散形式;
上式中,Un代表守恒变量U在第n迭代步下的值,Δt代表迭代步之间的时间间隔;
步骤202、式(3)等效为关于守恒变量U的非线性系统R(U)=0,采用Newton迭代方法对非线性系统进行迭代求解;
通过Newton迭代,非线性系统转化为一系列线性系统的逼近,式(4)中,R(U)为公式(3)中的左侧部分,即
包括时间积分、体积分、面积分和提升算子积分项,
作为优选技术措施,
所述步骤三,其具体包括以下步骤:
步骤301、用基函数的线性组合代替守恒变量U,其具体形式为:
其中,φj代表基函数,j为基函数的下标,uj(t)代表对应基函数的系数,又被称作自由度,是一个关于时间的函数,后续的公式中用uj简化表示uj(t);N代表总的基函数的个数,是一个与计算阶数直接相关的量;
步骤302、将线化矩阵中的守恒变量用(6)式来代替,分别能够得到时间积分、体积分、面积分和提升算子积分项的线化矩阵;
步骤303、计算时间积分项的线化矩阵;
步骤304、利用链式法则和全微分,对体积分项进行分解;
上式中,通量对守恒变量或者守恒变量的梯度的导数能够通过差分微小扰动的方式求出,守恒变量对自由度的偏导数能够转换成基函数,而守恒变量的梯度对自由度的偏导数能够转化为基函数的梯度;
步骤305、同样利用链式法则和全微分,对面积分项和提升算子项进行分解,并计算对应项的线化矩阵;
步骤305、将得到的线化矩阵带入步骤202,更写迭代,得到最终计算结果。
本发明解除了线化矩阵对自由度的依赖,不用对每个变量的自由度进行扰动,只需要对守恒变量或者守恒变量的梯度进行相应的扰动就能求解线化矩阵,减少计算量,同时兼顾线化矩阵的准确度,保证隐式迭代有高效的收敛性,能够快速得到流体仿真的结果和流场信息。相比其他方法,采用本发明方法所需的计算时间有效减少,计算效率得到了提升。
与现有技术相比,本发明具有以下有益效果:
1、本发明将依赖项对自由度的导数转换为对变量的导数与变量对自由度导数的乘积,以及对变量梯度的导数与变量梯度对自由度的导数乘积,解决了Riemann通量积分和边界积分项的导数计算问题,具备良好的通用性。
2、本发明解除了线化矩阵对自由度的依赖,不用对每个变量的自由度进行扰动,只需要对守恒变量或者守恒变量的梯度进行相应的扰动就能求解线化矩阵,减少计算量,提高计算效率,也更具有通用性,利于编程实现,同时兼顾线化矩阵的准确度,保证隐式迭代有高效的收敛性,能够快速得到流体仿真的结果和流场信息。
附图说明
图1为本发明方法的流程图。
图2为本发明采用Newton/GMRES隐式时间推进的流程图。
图3为采用本发明方法计算得到的不同精度下二维圆柱绕流问题的压力分布等值线示图。
图4为计算二维圆柱绕流问题时本发明与其他方法在计算量、CPU计算时间上的对比示图。
图5为采用本发明计算NACA0012二维翼型问题得到的Mach数和流场分布情况示图。
图6为计算NACA0012二维翼型问题时本发明的收敛曲线与传统方法的对比示图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合附图及实施例,对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅仅用以解释本发明,并不用于限定本发明。
相反,本发明涵盖任何由权利要求定义的在本发明的精髓和范围上做的替代、修改、等效方法以及方案。进一步,为了使公众对本发明有更好的了解,在下文对本发明的细节描述中,详尽描述了一些特定的细节部分。对本领域技术人员来说没有这些细节部分的描述也可以完全理解本发明。
如图1-2所示,基于线化矩阵混合求解模式的高阶CFD隐式时间推进方法,包括以下步骤:
步骤1、建立基于非结构混合网格求解Navier-Stokes方程的高精度间断迦辽金框架,主要用来进行空气动力学的数值仿真计算,包括Navier-Stokes控制方程、基函数、空间离散。
其具体包括以下步骤:
步骤101、构建微分形式下的Navier-Stokes方程
上式中,t代表时间,矢量U代表方程求解的守恒变量,
步骤102、将流体的计算区域剖分成若干不重合的离散单元κ,在不同的离散单元内,流场中的守恒变量U可以用基函数的线性组合来表示,选择单项式等级基作为基函数,单元内的积分点根据网格单元的类型来确定。
步骤103、在方程两边同时乘以基函数φi并进行积分,得到弱形式下的控制方程(2)。
κ和
其中,
可以采用内罚方法进行计算。
步骤2、构建基于Newton/GMRES的隐式迭代框架。
其具体包括以下步骤:
步骤201、采用Euler后插方法对时间导数进行离散,从而得到控制方程的全离散形式。
上式中,Un代表守恒变量U在第n迭代步下的值,Δt代表迭代步之间的时间间隔。
步骤202、式(3)等效为关于守恒变量U的非线性系统R(U)=0,采用Newton迭代方法对非线性系统进行迭代求解。
通过Newton迭代,非线性系统转化为一系列线性系统的逼近,式(4)中,R(U)为公式(3)中的左侧部分,即
包含时间积分、体积分、面积分和提升算子积分项,
步骤3:采用混合求解方法计算线化矩阵,带入到控制方程,求解得到仿真计算结果。
其具体包括以下步骤:
步骤301、用基函数的线性组合代替守恒变量U,其具体形式为:
其中,φj代表基函数,j为基函数的下标,uj(t)代表对应基函数的系数,又被称作自由度,是一个关于时间的函数,后续的公式中用uj简化表示uj(t)。N代表总的基函数的个数,是一个与计算阶数直接相关的量。
步骤302、将线化矩阵中的守恒变量用(6)式来代替,分别可以得到时间积分、体积分、面积分和提升算子积分项的线化矩阵。
步骤303、计算时间积分项的线化矩阵。
步骤304、利用链式法则和全微分,对体积分项进行分解。
上式中,通量对守恒变量或者守恒变量的梯度的导数能够通过差分微小扰动的方式求出,守恒变量对自由度的偏导数能够转换成基函数,而守恒变量的梯度对自由度的偏导数可以转化为基函数的梯度。
步骤305、同样利用链式法则和全微分,对面积分项和提升算子项进行分解,并计算对应项的线化矩阵。
步骤305、将得到的线化矩阵带入步骤202,更写迭代,得到最终计算结果。
如图3-6所示,应用本发明计算圆柱问题和NACA0012翼型问题的具体实施例:
图4为计算二维圆柱绕流问题时,本发明与其他方法在计算量(迭代步数和CPU时间)的对比,纵坐标为计算残差。收敛到相同残差时,本发明在CPU计算时间上具有明显优势。
图6为计算NACA0012二维翼型问题时,本发明的收敛曲线与传统方法的对比。可以看到收敛到相同残差的情况下,本发明的计算时间上明显少于传统方法。
结论:在计算圆柱问题和NACA0012翼型问题时,相比其他方法,采用本发明方法所需的计算时间有效减少,计算效率得到了提升。
如表1所示,采用本发明的方法与手动解析方法和自动微分方法在不同精度下计算量的对比,从表1数据可以看出,本发明在求导计算量方面具备明显优势,且求导次数不会随着精度的提高而增加。
表1三维情况下本发明与传统方法线化矩阵计算量的对比
以上所述仅为本发明的较佳实施例而已,并不用以限制本发明,凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明的保护范围之内。
机译: 隐式油藏模拟矩阵方程的求解方法
机译: 隐式油藏模拟矩阵方程的求解方法
机译: 隐式储层模拟矩阵的求解方法