CFD2015-第9讲-求代数方程组及网格生成

计算流体力学讲义第九讲代数方程的求解及网格生成李新亮lixl；力学所主楼214；82543801,知识点：,1,CopyrightbyLiXinliang,代数方程组的求解网格生成,课件下载：,CopyrightbyLiXinliang,2,知识回顾：有限体积法基本流程,无粘项常用方法（流过AB边的通量）：a.利用周围点的值，计算出(I+1/2,J)点处的物理量；直接利用“差分格式”b.利用该处的物理量，计算出流过AB边的流通量迎风型方法需利用“通量分裂技术”,FVS类：,FDS类：,利用Riemann解,Reimann解：Godunov,Roe,HLL,HLLC,利用坐标变换，转化为一维Riemann问题,i,j,i-1,j,i+1,j,i,j+1,i,j+1,知识回顾：隐格式,代数方程组的求解,CopyrightbyLiXinliang,4,微分方程（组）,代数方程组,数值解（离散解）,差分有限体积,大部分计算量,9.1.1Gauss消去法,消元,为了计算稳定，通常使用主元消去法列主元消去法；全主元消去法计算量：乘法：加法：,9.1代数方程组求解的直接法,优点：简单精确，缺点：计算量大,上三角矩阵,CopyrightbyLiXinliang,5,9.1.2LU分解法,0,0,Step1,Step2,Stepk,对角线上不能有0，计算之前先交换矩阵A的元素，将主值交换到对角线上,CopyrightbyLiXinliang,6,回代过程,计算量：分解O(n3/3),回代O(n2),优点：1）重复求解，时，仅需一次LU分解，计算量小；2）LU分解不破坏带状稀疏矩阵的性质，可大幅减小计算量。,L带宽的带状矩阵：,LU分解:O(nL)回代：O(nL),CopyrightbyLiXinliang,7,9.1.3带状矩阵求解的追赶法,追赶法：等价于带状矩阵的LU分解,例：三对角矩阵,一般项：,边界项：,追赶法,令：,Step1:Step2:Step3:Step4:,计算量：9n次（乘法）A为固定值时，3n次（乘法）,简单易用，计算量小,CopyrightbyLiXinliang,8,9.2代数方程组求解的迭代法,9.2.1Jocabi及Gauss-Seidel迭代,解出对角元素,Jocabi迭代,Gauss-Seidel迭代,“对角占优”,CopyrightbyLiXinliang,9,9.2.2松弛迭代超松弛（SOR）、亚松弛,Step1:采用Jocabian或Gauss-Seidel迭代产生新的值,Step2:进行松弛,含义：改变步长,超松弛,精确解,“步子迈大一些”，加快收敛,亚松弛,“步子迈小一些”，稳定性好,收敛性：对角占优矩阵，Jocabian及Gauss-Seidel迭代可收敛,CopyrightbyLiXinliang,10,举例：Laplace方程的求解,五点格式,Jacobi迭代,Gauss-Seidel迭代,缺点:每迭代一步，信息只传递到周围网格点，n很大时收敛较慢,n+1,n,n,n,n,n+1,n,n+1,n+1,CopyrightbyLiXinliang,11,对称Gauss-Seidel迭代（SGS）,n+1,n,n+1,n+1,n+1,n,n,n,n+1,n+1,Step1,Step2,特点：两次扫描，反复迭代,CopyrightbyLiXinliang,12,9.2.3交替方向迭代（ADI）方法,例：,Step1:认为已知（使用上一步的值），求解三对角方程,得到中间步的值,Step2:代入中间步的值，求解三对角方程，得到n+1步的值,三对角方程采用追赶法求解，效率较高在每一条线上采用直接法，信息快速传递，有利于收敛,Step3:重复以上两个步骤，直至收敛,因追赶法实际上是LU分解法，因此又称LU-ADI方法,n+1,n+1,n+1,n,n+1,n+1,n+1,n,n,n,CopyrightbyLiXinliang,13,9.2.4近似解法LU-SGS方法,SGS方法信息传递速度仍较慢，需要加速,近似LU分解,Step1:,Step2:,优点：不含任何迭代过程，两步扫描即可完成，效率很高；缺点：近似LU分解，结果不够准确,OK，不迭代,是LU分解的SGS方法，因此成为LU-SGS近似解法,CopyrightbyLiXinliang,14,9.2.5加速收敛的多重网格法,Gauss-Seidel迭代,含义：线性系统，误差满足同样的方程,定义误差：,1）收敛速度的Fourier分析,增长（收敛）因子,含义：,极端高波数情况,迭代一次，误差减小一半,极端低波数情况,收敛速度趋近于0,CopyrightbyLiXinliang,15,策略：多重网格粗网格加速低波数扰动收敛，细网格加速高波数收敛,细网格,粗网格,使用多重网格法求解方程：,迭代方程：,以Jacobian迭代为例,修正方程,Step1:在细网格上迭代一定步数（无需收敛），得到中间步的值Step2:将修正项插值到粗网格上，并迭代求解Step3:将求解后的修正项插值到细网格，并计算出细网格上新的值Step4:重复Step1-3直到收敛,修正项,CopyrightbyLiXinliang,16,常用方法：V型及W型迭代,细网格,粗网格,更粗网格,细网格,粗网格,更粗网格,V型迭代W型迭代,CopyrightbyLiXinliang,17,9.2.6共轭梯度法,1.求解对称正定矩阵的共轭梯度法,化代数方程组问题为极值问题（设A为对称正定矩阵）,的最小值问题,可用最速下降法之类方法求解,例：,的分布,最小值点（0,1）,CopyrightbyLiXinliang,18,求极值问题的最速下降法,思路：沿当地梯度方向前进,1）根据当前位置，计算当地梯度方向：,2）沿该方向前进，使得达到极小值,方向,步长,“残余向量”,特点：沿当地梯度方向前进，直到不能前进为止，然后以按照新的梯度方向前进；相邻路径方向正交；缺点：局部最速下降路线并非全局最速下降路线，因而收敛速度并非最优。,“最速下降法”示意图,CopyrightbyLiXinliang,19,最速下降法的改进：共轭梯度法,最速下降法：每步在一维空间求最优解：,改进：在二维空间寻求最优解（不再沿当地梯度方向）,旧路径方向及当前梯度方向所张成的二维空间,寻找该平面内的极小值,解出,新线路：,修正,最终得：,特点：相邻两次方向关于矩阵A正交（称为共轭）。,CopyrightbyLiXinliang,20,2.求解一般非奇异方程组的共轭梯度法,A为一般满秩阵（非对称正定阵）的情况,正则化方法：,对称正定阵,Step1:设置初值,Step2:迭代推进,直到残余向量足够小为止,CopyrightbyLiXinliang,21,9.3代数网格生成法,基本思路：通过代数方程计算出网格点的位置优点：灵活、计算量小缺点：光滑性差，过于依赖人工,如图：叶栅通道,已知计算域上边界（红线）及下边界（蓝线）的方程为：和,则网格为：,其中可控制法向的疏密分布,均匀分布；,在下壁面处密集分布,上下壁面两侧加密,CopyrightbyLiXinliang,22,9.4椭圆形方程网格生成法,A,B,C,D,E,F,A,B,C,D,E,F,对于如图单联通的计算域，可通过坐标变换,变换到图示矩形计算域,物理空间边界计算空间边界物理空间内点计算空间内点,物理空间,计算空间,通常：给定边界点的对应关系（代数方法）通过求解方程获得内点的对应关系,方程的边值问题,椭圆型方程边值问题抛物型方程双曲型方程,初边值问题,椭圆型方程,CopyrightbyLiXinliang,23,通常：,或,含义：给物理空间的每个点找到计算空间的对应位置。,注：由于拓扑的对应性，物理空间必须是单联通域,如果是多联通的，可通过切割，形成单联通域,CopyrightbyLiXinliang,24,求解,1.形式变换，改写成以为自变量,便于进行差分求解,CopyrightbyLiXinliang,25,离散化：中心差分,离散方程,迭代求解：Jacobi,Seidel,SOR,LU-ADI,LU-SGS,多重网格，,CopyrightbyLiXinliang,26,Step1:确定边界网格,通常采用代数方法生成一维网格，容易生成注意：1）边界对应关系容易出错2）考虑网格的疏密分布（翼型尾缘区，激波区，近壁区),Step2:利用上页的离散方程，解出全部网格坐标,不足：内部区域的网格分布不易控制无法做到指定区域网格加密无法保证网格正交,边界网格可控，内部网格只能“听天由命”,方案1：源项P,Q为0，求解Laplace方程,CopyrightbyLiXinliang,27,方案2：设定源项P,Q求解Poison方程,源项P,Q对网格的影响数值实验发现，在某点处加入点源P:P0使方向网格线发散,点源P,P0,P0,在某点处加入点源Q,可对方向网格线产生同样效果,1）网格线的汇聚,启发：在某条网格线上加入负的源项，可令网格汇聚,使网格汇聚于,CopyrightbyLiXinliang,28,2）边界网格的正交（并指定边界网格间距）,令：,源项在边界处，内部衰减,利用边界处网格的正交性及网格间距要求确定系数P和Q,指定值,基本思路：以壁面线处为例,网格线正交,指定法向网格间距,计算第2层网格线上的坐标,通过差分计算边界处的,如需要利用的信息，可用上一迭代步的值,计算出边界处的P,Q,根据指数衰减原则给出全场的P,Q,具体公式见：傅德薰计算空气动力学284-286,CopyrightbyLiXinliang,29,习题9.1网格生成,通过解椭圆型方程生成NA