CN110619160B - 一种基于伴随残差排序的隐式解法 - Google Patents

一种基于伴随残差排序的隐式解法 Download PDF

Info

Publication number
CN110619160B
CN110619160B CN201910822819.9A CN201910822819A CN110619160B CN 110619160 B CN110619160 B CN 110619160B CN 201910822819 A CN201910822819 A CN 201910822819A CN 110619160 B CN110619160 B CN 110619160B
Authority
CN
China
Prior art keywords
residual
calculation
flow field
influence
influence quantity
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
Application number
CN201910822819.9A
Other languages
English (en)
Other versions
CN110619160A (zh
Inventor
胡国风
吴泓宇
刘可
杨坤
王喆
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sichuan Tengdun Technology Co Ltd
Original Assignee
Sichuan Tengdun Technology Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Sichuan Tengdun Technology Co Ltd filed Critical Sichuan Tengdun Technology Co Ltd
Priority to CN201910822819.9A priority Critical patent/CN110619160B/zh
Publication of CN110619160A publication Critical patent/CN110619160A/zh
Application granted granted Critical
Publication of CN110619160B publication Critical patent/CN110619160B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Abstract

本发明公开了一种基于伴随残差排序的隐式解法,包括以下步骤:(1)针对待测对象建立流场控制方程;(2)设置边界条件和初始流场条件计算各单元格的残差;(3)根据流场控制方程的伴随方程计算残差影响量估计值;(4)根据残差影响量估计值对计算结果的影响进行排序计算。本发明能保证计算全部是有效计算,并且随着信号的传播过程,可保证计算区域始终覆盖信号传播的路径,最大限度地减少冗余计算量;计算域始终由扰动信号的传播驱动,自动区分有效区域和无效区域;由于扰动不向上游传播,可以有效缩减不影响计算结果的下游区域的无效计算;还可利用伴随矩阵获得最终计算结果的精度,并对结果进行修正,以较少的计算量获得较高的计算精度。

Description

一种基于伴随残差排序的隐式解法
技术领域
本发明涉及计算流体动力学技术领域,尤其涉及一种基于伴随残差排序的隐式解法。
背景技术
在计算流体动力学技术(CFD)领域,大部分计算模型都可以用如下矩阵方程求解:
Ax+b=0 (1)
其中矩阵A为代表控制方程的矩阵,b表示边界条件,x为待求向量。
通常解流场控制方程(1)的方法是采用时间相关法,对待求向量x设一个初值x0,经历一个随时间变化的历程
Figure GDA0002282195720000011
最终得到式(1)的解,其中R为残差。
稳态计算的目标是使残差趋向于0,获得流场目标解,R的大小为流场收敛的主要标志。CFD计算的一般过程是设初始流场,在物面(或其他限制性边界)以外的流场区域残差为0或接近于0,限制性边界附近由边界条件限制而产生非零的残差,在计算迭代过程中向外扩散、在边界上反射、可能还有振荡过程,然后逐渐减小幅值并趋向于0。残差的变化实际上代表了由限制性边界条件产生的扰动信号在流动方程的控制下的传播与相互作用过程。
对于超声速/高超声速流场,存在扰动的有限影响区(通常是以激波为分界线),扰动的影响局限在影响区内,影响区外的流场永远不变;同时还存在信号传播的方向性,下游的变化不会越过马赫锥影响到上游区域。基于这些特点,提出了一些方法以减少扰动影响区外不必要的计算,如基于预先分析对计算域进行限定,或将控制方程简化成抛物型(PNS),进行空间推进计算等。
对于亚声速流场,扰动的影响会传播到整个流场,但仍然存在传播方向和衰减的规律,相应的有交替方向扫描(ADI)等针对信号传播方向性的加速方法。
但这些方法均存在局限性:计算域限定不能用一个方案对所有计算状态都达到最高的效率;PNS则对方程进行了改变,影响精度,并只能适用于超声速区域;ADI也存在一定的盲目性,没有精确跟踪信号的传播路径,存在较大程度的计算量冗余。在计算过程中,也常会出现流场中有局部区域振荡收敛性差,而使全流场进行长时间迭代计算的情况,也有用全流场收敛准则判断收敛,但局部仍存在较大偏差的情况。这些都会大大影响计算效率和计算精度。
基于残差排序的隐式解法,对全流场所有网格的残差进行排序管理,优先计算残差排序靠前的单元。但是在流场中不同区域的残差对最终计算结果的影响是不一样的,特别是在超声速流场中,由于扰动的影响不能向上游传播,下游的残差对计算结果完全没有影响,没有计算的必要,用残差排序并不能完全避免这一部分的重复计算。
发明内容
为了解决上述问题,减少重复的计算量,提高流场控制方程求解效率,即以最优的计算步骤和顺序计算方程,最大限度加速方程的收敛过程,本发明提出一种基于伴随残差排序的隐式解法,具体的,包括以下步骤:
步骤1:针对待测对象建立流场控制方程:Ax+b=0;其中A为流场控制方程的矩阵;b为向量,表示根据待测对象建立的边界条件;x为向量,表示待求的流场参数;
步骤2:设置边界条件和初始流场条件,一般地,由于初始流场条件下各单元格均不满足所给的矩阵方程,具有残差,故计算各单元格的残差R;
步骤3:根据流场控制方程的伴随方程
Figure GDA0002282195720000021
得到每个单元残差变化量δR与计算结果变化量关系的估计式:
Figure GDA0002282195720000022
其中F为待求的气动力,Λ为伴随矩阵,最终计算收敛的结果与当前状态相比,残差变化量δR=-R,因此最终结果与当前结果的变化量估计值为:
δF=-ΛδR=ΛR
从而可以得到每个单元残差对最终计算结果的影响量估计值,即残差影响量估计值ΛR;
步骤4:根据残差影响量估计值ΛR对计算结果的影响进行排序计算。
进一步的,所述步骤4还包括以下子步骤:
步骤4-1:按残差影响量估计值的大小对各单元格排序建立一个求解队列,残差影响量估计值最大的单元格位于队首;
步骤4-2:从队首取出残差影响量估计值最大的单元格进行控制方程迭代求解,直到该残差影响量估计值小于当前队首单元格的残差影响量估计值;
步骤4-3:计算更新当前单元格临近单元格的残差和残差影响量估计值,将当前单元格及受当前单元格变化而影响残差发生变化的单元格,按残差影响量估计值大小放入队列中相应位置;
步骤4-4:重复步骤4-2和步骤4-3,直到残差影响量估计值最大值减少一个数量级,然后重复步骤3和步骤4-1;
步骤4-5:重复步骤4-4,直到总的残差影响量估计值小于预先给定的收敛判据值,从而得到最终的流场计算收敛结果;
步骤4-6:求最后剩余残差对计算结果的影总响量,并加到最后计算结果中,以提高最后计算结果的精度。
进一步的,若计算模型物面在某单元的影响区外,其流场变化扰动对物面附近的流场没有影响,则该单元的伴随矩阵值就等于0。
本发明的有益效果在于:保证计算全部是有效计算,并且随着信号的传播过程,可保证计算区域始终覆盖信号传播的路径,最大限度地减少冗余计算量;计算域始终由扰动信号的传播驱动,自动区分有效区域和无效区域;特别对于超声速流场的计算,由于扰动不向上游传播,可以有效缩减不影响计算结果的下游区域的无效计算;可以适应所有计算状态,无需针对不同计算状态进行特殊处理;该方法只控制计算域,对算法层没有约束和特殊要求,可以适用于任何算法,不会牺牲精度或限制算法的适用范围。最后,可以利用伴随矩阵获得最终计算结果的精度,并对结果进行修正,以较少的计算量获得较高的计算精度。
具体实施方式
为了对本发明的技术特征、目的和效果有更加清楚的理解,现说明本发明的具体实施方式。应当理解,此处所描述的具体实施例仅用以解释本发明,并不用于限定本发明,即所描述的实施例仅仅是本发明一部分实施例,而不是全部的实施例。基于本发明的实施例,本领域技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例,都属于本发明保护的范围。
实施例1
伴随矩阵方法主要用于求外形变化对气动力影响的方法,针对方程(1),设伴随矩阵Λ,满足:
Figure GDA0002282195720000031
利用伴随矩阵Λ可以获得流场或边界条件的变化对气动力求解结果的影响,常用于获得改变边界条件对气动力特性的影响,以进行外形的参数优化。
因此,本实施例提供了一种基于伴随残差排序的隐式解法,包括以下步骤:
步骤1:针对待测对象建立流场控制方程:Ax+b=0;其中A为流场控制方程的矩阵;b为向量,表示根据待测对象建立的边界条件;x为向量,表示待求的流场参数;
步骤2:设置边界条件和初始流场条件,一般地,由于初始流场条件下各单元格均不满足所给的矩阵方程,具有残差,故计算各单元格的残差R;
步骤3:根据流场控制方程的伴随方程
Figure GDA0002282195720000041
得到每个单元残差变化量δR与计算结果变化量关系的估计式:
Figure GDA0002282195720000042
其中F为待求的气动力,Λ为伴随矩阵,最终计算收敛的结果与当前状态相比,残差变化量δR=-R,因此最终结果与当前结果的变化量估计值为:
δF=-ΛδR=ΛR
从而可以得到每个单元残差对最终计算结果的影响量估计值,即残差影响量估计值ΛR;
步骤4:根据残差影响量估计值ΛR对计算结果的影响进行排序计算。
其中,步骤4还包括以下子步骤:
步骤4-1:按残差影响量估计值的大小对各单元格排序建立一个求解队列,残差影响量估计值最大的单元格位于队首;
步骤4-2:从队首取出残差影响量估计值最大的单元格进行控制方程迭代求解,直到该残差影响量估计值小于当前队首单元格的残差影响量估计值;
步骤4-3:计算更新当前单元格临近单元格的残差和残差影响量估计值,将当前单元格及受当前单元格变化而影响残差发生变化的单元格,按残差影响量估计值大小放入队列中相应位置;
步骤4-4:重复步骤4-2和步骤4-3,直到残差影响量估计值最大值减少一定程度,优选地为减少一个数量级,然后重复步骤3和步骤4-1;
步骤4-5:重复步骤4-4,直到总的残差影响量估计值小于预先给定的收敛判据值,从而得到最终的流场计算收敛结果;
步骤4-6:求最后剩余残差对计算结果的影总响量,并加到最后计算结果中,以提高最后计算结果的精度。
此外,若计算模型物面在某单元的影响区外,其流场变化扰动对物面附近的流场没有影响,则该单元的伴随矩阵值就等于0。
以上所述仅是本发明的优选实施方式,应当理解本发明并非局限于本文所披露的形式,不应看作是对其他实施例的排除,而可用于各种其他组合、修改和环境,并能够在本文所述构想范围内,通过上述教导或相关领域的技术或知识进行改动。而本领域人员所进行的改动和变化不脱离本发明的精神和范围,则都应在本发明所附权利要求的保护范围内。

Claims (2)

1.一种基于伴随残差排序的隐式解法,其特征在于,包括以下步骤:
步骤1:针对待测对象建立流场控制方程:Ax+b=0;其中A为流场控制方程的矩阵;b为向量,表示根据待测对象建立的边界条件;x为向量,表示待求的流场参数;
步骤2:设置边界条件和初始流场条件,由于初始流场条件下各单元格均不满足所给的矩阵方程,具有残差,故计算各单元格的残差R;
步骤3:根据流场控制方程的伴随方程
Figure FDA0003873331800000011
得到每个单元残差变化量δR与计算结果变化量关系的估计式:
Figure FDA0003873331800000012
其中F为待求的气动力,Λ为伴随矩阵,最终计算收敛的结果与当前状态相比,残差变化量δR=-R,因此最终结果与当前结果的变化量估计值为:
δF=-ΛδR=ΛR
从而可以得到每个单元残差对最终计算结果的影响量估计值,即残差影响量估计值ΛR;
步骤4:根据残差影响量估计值ΛR对计算结果的影响进行排序计算;
所述步骤4还包括以下子步骤:
步骤4-1:按残差影响量估计值的大小对各单元格排序建立一个求解队列,残差影响量估计值最大的单元格位于队首;
步骤4-2:从队首取出残差影响量估计值最大的单元格进行控制方程迭代求解,直到该残差影响量估计值小于当前队首单元格的残差影响量估计值;
步骤4-3:计算更新当前单元格临近单元格的残差和残差影响量估计值,将当前单元格及受当前单元格变化而影响残差发生变化的单元格,按残差影响量估计值大小放入队列中相应位置;
步骤4-4:重复步骤4-2和步骤4-3,直到残差影响量估计值最大值减少一个数量级,然后重复步骤3和步骤4-1;
步骤4-5:重复步骤4-4,直到总的残差影响量估计值小于预先给定的收敛判据值,从而得到最终的流场计算收敛结果;
步骤4-6:求最后剩余残差对计算结果的影总响量,并加到最后计算结果中,以提高最后计算结果的精度。
2.根据权利要求1所述的一种基于伴随残差排序的隐式解法,其特征在于,若计算模型物面在某单元的影响区外,其流场变化扰动对物面附近的流场没有影响,则该单元的伴随矩阵值就等于0。
CN201910822819.9A 2019-09-02 2019-09-02 一种基于伴随残差排序的隐式解法 Active CN110619160B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910822819.9A CN110619160B (zh) 2019-09-02 2019-09-02 一种基于伴随残差排序的隐式解法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910822819.9A CN110619160B (zh) 2019-09-02 2019-09-02 一种基于伴随残差排序的隐式解法

Publications (2)

Publication Number Publication Date
CN110619160A CN110619160A (zh) 2019-12-27
CN110619160B true CN110619160B (zh) 2022-12-02

Family

ID=68922976

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910822819.9A Active CN110619160B (zh) 2019-09-02 2019-09-02 一种基于伴随残差排序的隐式解法

Country Status (1)

Country Link
CN (1) CN110619160B (zh)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US3299255A (en) * 1962-10-18 1967-01-17 Sperry Rand Corp Fluid data comparator
CN107256290A (zh) * 2017-05-19 2017-10-17 四川腾盾科技有限公司 一种基于残差摄动法的边界条件变化对流场影响差量的计算方法
CN107357976A (zh) * 2017-06-27 2017-11-17 四川腾盾科技有限公司 一种飞行器的动导数的计算方法
CN108563843A (zh) * 2018-03-26 2018-09-21 北京航空航天大学 定常可压缩流动的扰动区域更新方法
CN108932554A (zh) * 2017-05-26 2018-12-04 西安交通大学 一种风电场流场量测点的配置优化方法及装置

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7117108B2 (en) * 2003-05-28 2006-10-03 Paul Ernest Rapp System and method for categorical analysis of time dependent dynamic processes
US8380473B2 (en) * 2009-06-13 2013-02-19 Eric T. Falangas Method of modeling dynamic characteristics of a flight vehicle

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US3299255A (en) * 1962-10-18 1967-01-17 Sperry Rand Corp Fluid data comparator
CN107256290A (zh) * 2017-05-19 2017-10-17 四川腾盾科技有限公司 一种基于残差摄动法的边界条件变化对流场影响差量的计算方法
CN108932554A (zh) * 2017-05-26 2018-12-04 西安交通大学 一种风电场流场量测点的配置优化方法及装置
CN107357976A (zh) * 2017-06-27 2017-11-17 四川腾盾科技有限公司 一种飞行器的动导数的计算方法
CN108563843A (zh) * 2018-03-26 2018-09-21 北京航空航天大学 定常可压缩流动的扰动区域更新方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
An implicit,exact dual adjoint solution method for turbulent flows on unstructured grids;Eric J Nielsen et al.;《Computers & Fluid》;20041130;1131-1155 *
基于非结构混合网格的离散伴随优化方法研究;李彬;《中国博士学位论文全文数据库 工程科技II辑》;20161215;C031-4 *

Also Published As

Publication number Publication date
CN110619160A (zh) 2019-12-27

Similar Documents

Publication Publication Date Title
CN108304600B (zh) 一种高超声速飞行器转捩位置预测方法
CA2887795C (en) Closed loop control of aircraft control surfaces
Baiihaus et al. Implicit approximate-factorization schemes for steady transonic flow problems
CN107977491B (zh) 一种非稳态情况下飞行器空气舵缝隙的气动热评估方法
Prananta et al. Two-dimensional transonic aeroelastic analysis using thin-layer Navier–Stokes method
CN114611416B (zh) 导弹非线性非定常气动特性ls-svm建模方法
CN114444216A (zh) 基于数值模拟的高空条件下飞行器姿态控制方法及系统
CN110619160B (zh) 一种基于伴随残差排序的隐式解法
CN110659447B (zh) 一种基于影响因子残差排序的隐式解法
CN110647716B (zh) 一种基于残差排序的隐式解法
CN110674607B (zh) 一种基于残差量级排序的隐式解法
CN112214835A (zh) 一种旋翼悬停状态气动噪声工程估算方法
Srivastava et al. Flutter analysis of a transonic fan
CN114489098B (zh) 一种飞行器的姿态控制方法及飞行器
CN113569360B (zh) 一种风力机叶片抗颤振翼型簇设计方法
CN106777918A (zh) 基于导波和模糊算法的功能梯度结构材料特性的反演方法
CN111177645A (zh) 一种基于大规模点云数据的大型高速回转装备误差混合评定方法
CN109752086A (zh) 基于bellhop的快速声场计算方法
CN109359372B (zh) 考虑结构-载荷-边界耦合影响的结构拓扑优化设计方法
CN114676378B (zh) 基于rad求解器的激波计算方法
CN112949222B (zh) 一种准乘波体构型横航向稳定性优化方法、系统及升力体
CN112162563B (zh) 基于自适应弱敏无迹Kalman滤波的直升机状态估计方法
Ambrosi et al. Full potential and Euler solutions for transonic unsteady flow
Guo et al. A New Adaptive Tracking Algorithm for Near-Space Hypersonic Target
Birhane et al. A Computational Framework for the Aerodynamic Shape Optimization of Long-Span Bridge Decks

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