計(jì)算流體力學(xué)課件：求代數(shù)方程組及網(wǎng)格生成

代數(shù)方程的求解及網(wǎng)格生成知識點(diǎn)：

1代數(shù)方程組的求解網(wǎng)格生成2知識回顧：有限體積法基本流程無粘項(xiàng)常用方法（流過AB邊的通量）：a.利用周圍點(diǎn)的值，計(jì)算出(I+1/2,J)點(diǎn)處的物理量；

直接利用“差分格式”b.利用該處的物理量，計(jì)算出流過AB邊的流通量

迎風(fēng)型方法需利用“通量分裂技術(shù)”FVS類：FDS類：利用Riemann解Reimann解：Godunov,Roe,HLL,HLLC利用坐標(biāo)變換，轉(zhuǎn)化為一維Riemann問題i,ji-1,ji+1,ji,j+1i,j+1知識回顧：隱格式代數(shù)方程組的求解4微分方程（組）代數(shù)方程組數(shù)值解（離散解）差分有限體積大部分計(jì)算量9.1.1Gauss消去法消元為了計(jì)算穩(wěn)定，通常使用主元消去法列主元消去法；全主元消去法計(jì)算量：乘法：

加法：9.1代數(shù)方程組求解的直接法優(yōu)點(diǎn)：簡單精確，缺點(diǎn)：計(jì)算量大上三角矩陣59.1.2LU分解法……00Step1Step2Stepk對角線上不能有0，計(jì)算之前先交換矩陣A的元素，將主值交換到對角線上6回代過程計(jì)算量：分解O(n3/3),回代O(n2)優(yōu)點(diǎn)：1）重復(fù)求解，……時(shí)，僅需一次LU分解，計(jì)算量小；2）LU分解不破壞帶狀稀疏矩陣的性質(zhì)，可大幅減小計(jì)算量。L帶寬的帶狀矩陣：LU分解:O(nL)回代：O(nL)79.1.3帶狀矩陣求解的追趕法追趕法：等價(jià)于帶狀矩陣的LU分解例：三對角矩陣一般項(xiàng)：邊界項(xiàng)：追趕法令：Step1:Step2:

Step3:Step4:

計(jì)算量：9n次（乘法）A為固定值時(shí)，3n次（乘法）簡單易用，計(jì)算量小89.2代數(shù)方程組求解的迭代法直接法迭代法優(yōu)點(diǎn)算法簡便，準(zhǔn)確（未知數(shù)少時(shí)）計(jì)算量小，誤差容易控制缺點(diǎn)計(jì)算量大O(n3)舍入誤差積累，不易控制快速收斂的算法設(shè)計(jì)較為復(fù)雜9.2.1Jocabi及Gauss-Seidel迭代解出對角元素Jocabi迭代Gauss-Seidel迭代“對角占優(yōu)”99.2.2松弛迭代——超松弛（SOR）、亞松弛Step1:采用Jocabian或Gauss-Seidel迭代產(chǎn)生新的值Step2:進(jìn)行松弛含義：改變步長超松弛精確解“步子邁大一些”，加快收斂亞松弛“步子邁小一些”，穩(wěn)定性好收斂性：

對角占優(yōu)矩陣，Jocabian及Gauss-Seidel迭代可收斂10舉例：Laplace方程的求解五點(diǎn)格式Jacobi迭代Gauss-Seidel迭代缺點(diǎn):每迭代一步，信息只傳遞到周圍網(wǎng)格點(diǎn)，n很大時(shí)收斂較慢n+1nnnnn+1nn+1n+111對稱Gauss-Seidel迭代（SGS）n+1nn+1n+1n+1nnnn+1n+1Step1Step2特點(diǎn)：兩次掃描，反復(fù)迭代129.2.3交替方向迭代（ADI）方法例：Step1:認(rèn)為已知（使用上一步的值），求解三對角方程,得到中間步的值Step2:代入中間步的值，求解三對角方程，得到n+1步的值

三對角方程采用追趕法求解，效率較高在每一條線上采用直接法，信息快速傳遞，有利于收斂Step3:重復(fù)以上兩個(gè)步驟，直至收斂因追趕法實(shí)際上是LU分解法，因此又稱LU-ADI方法n+1n+1n+1nn+1n+1n+1nnn139.2.4近似解法——LU-SGS方法

SGS方法信息傳遞速度仍較慢，需要加速

近似LU分解Step1:Step2:優(yōu)點(diǎn)：不含任何迭代過程，兩步掃描即可完成，效率很高；缺點(diǎn)：近似LU分解，結(jié)果不夠準(zhǔn)確OK，

不迭代是LU分解的SGS方法，因此成為LU-SGS近似解法149.2.5加速收斂的多重網(wǎng)格法Gauss-Seidel迭代含義：線性系統(tǒng)，誤差滿足同樣的方程定義誤差：1）收斂速度的Fourier分析增長（收斂）因子含義：極端高波數(shù)情況迭代一次，誤差減小一半極端低波數(shù)情況收斂速度趨近于015策略：多重網(wǎng)格——粗網(wǎng)格加速低波數(shù)擾動(dòng)收斂，細(xì)網(wǎng)格加速高波數(shù)收斂細(xì)網(wǎng)格粗網(wǎng)格使用多重網(wǎng)格法求解方程：迭代方程：以Jacobian迭代為例修正方程Step1:在細(xì)網(wǎng)格上迭代一定步數(shù)（無需收斂），得到中間步的值Step2:將修正項(xiàng)插值到粗網(wǎng)格上，并迭代求解Step3:將求解后的修正項(xiàng)插值到細(xì)網(wǎng)格，并計(jì)算出細(xì)網(wǎng)格上新的值Step4:重復(fù)Step1-3直到收斂修正項(xiàng)16常用方法：V型及W型迭代細(xì)網(wǎng)格粗網(wǎng)格更粗網(wǎng)格細(xì)網(wǎng)格粗網(wǎng)格更粗網(wǎng)格V型迭代W型迭代179.2.6共軛梯度法1.求解對稱正定矩陣的共軛梯度法化代數(shù)方程組問題為極值問題（設(shè)A為對稱正定矩陣）的最小值問題可用最速下降法之類方法求解例：的分布最小值點(diǎn)（0,1）18求極值問題的最速下降法思路：沿當(dāng)?shù)靥荻确较蚯斑M(jìn)1）根據(jù)當(dāng)前位置，計(jì)算當(dāng)?shù)靥荻确较颍?/p>

2）沿該方向前進(jìn)，使得達(dá)到極小值方向步長“殘余向量”特點(diǎn)：沿當(dāng)?shù)靥荻确较蚯斑M(jìn)，直到不能前進(jìn)為止，然后以按照新的梯度方向前進(jìn)；相鄰路徑方向正交；缺點(diǎn)：局部最速下降路線并非全局最速下降路線，因而收斂速度并非最優(yōu)。“最速下降法”示意圖19最速下降法的改進(jìn)：共軛梯度法最速下降法：每步在一維空間求最優(yōu)解：改進(jìn)：在二維空間尋求最優(yōu)解（不再沿當(dāng)?shù)靥荻确较颍┡f路徑方向及當(dāng)前梯度方向所張成的二維空間尋找該平面內(nèi)的極小值解出新線路：修正最終得：特點(diǎn)：相鄰兩次方向關(guān)于矩陣A正交（稱為共軛）。202.求解一般非奇異方程組的共軛梯度法A為一般滿秩陣（非對稱正定陣）的情況正則化方法：對稱正定陣Step1:設(shè)置初值

Step2:迭代推進(jìn)直到殘余向量足夠小為止219.3代數(shù)網(wǎng)格生成法基本思路：通過代數(shù)方程計(jì)算出網(wǎng)格點(diǎn)的位置優(yōu)點(diǎn)：靈活、計(jì)算量小缺點(diǎn)：光滑性差，過于依賴人工如圖：葉柵通道已知計(jì)算域上邊界（紅線）及下邊界（藍(lán)線）的方程為：和則網(wǎng)格為：其中可控制法向的疏密分布均勻分布；在下壁面處密集分布上下壁面兩側(cè)加密229.4橢圓形方程網(wǎng)格生成法ABCDEFA’B’C’D’E’F’對于如圖單聯(lián)通的計(jì)算域，可通過坐標(biāo)變換變換到圖示矩形計(jì)算域物理空間邊界

計(jì)算空間邊界物理空間內(nèi)點(diǎn)計(jì)算空間內(nèi)點(diǎn)物理空間計(jì)算空間通常：給定邊界點(diǎn)的對應(yīng)關(guān)系（代數(shù)方法）

通過求解方程獲得內(nèi)點(diǎn)的對應(yīng)關(guān)系方程的邊值問題橢圓型方程

邊值問題拋物型方程雙曲型方程初邊值問題橢圓型方程23通常：或含義：給物理空間的每個(gè)點(diǎn)找到計(jì)算空間的對應(yīng)位置。注：

由于拓?fù)涞膶?yīng)性，物理空間必須是單聯(lián)通域如果是多聯(lián)通的，可通過切割，形成單聯(lián)通域24求解1.形式變換，改寫成以為自變量便于進(jìn)行差分求解

25離散化：中心差分離散方程迭代求解：Jacobi,Seidel,SOR,LU-ADI,LU-SGS,多重網(wǎng)格，……26Step1:確定邊界網(wǎng)格通常采用代數(shù)方法生成——一維網(wǎng)格，容易生成注意：1）邊界對應(yīng)關(guān)系容易出錯(cuò)2）考慮網(wǎng)格的疏密分布（翼型尾緣區(qū)，激波區(qū)，近壁區(qū)……)Step2:利用上頁的離散方程，解出全部網(wǎng)格坐標(biāo)

不足：內(nèi)部區(qū)域的網(wǎng)格分布不易控制

無法做到指定區(qū)域網(wǎng)格加密

無法保證網(wǎng)格正交邊界網(wǎng)格可控，內(nèi)部網(wǎng)格只能“聽天由命”方案1：源項(xiàng)P,Q為0，求解Laplace方程27方案2：設(shè)定源項(xiàng)P,Q求解Poison方程源項(xiàng)P,Q對網(wǎng)格的影響

數(shù)值實(shí)驗(yàn)發(fā)現(xiàn)，在某點(diǎn)處加入點(diǎn)源P:P<0使方向網(wǎng)格線匯聚P>0使方向網(wǎng)格線發(fā)散

點(diǎn)源PP>0P<0在某點(diǎn)處加入點(diǎn)源Q,可對方向網(wǎng)格線產(chǎn)生同樣效果1）網(wǎng)格線的匯聚啟發(fā)：在某條網(wǎng)格線上加入負(fù)的源項(xiàng)，可令網(wǎng)格匯聚使網(wǎng)格匯聚于282）邊界網(wǎng)格的正交（并指定邊界網(wǎng)格間距）令：源項(xiàng)在邊界處，內(nèi)部衰減利用邊界處網(wǎng)格的正交性及網(wǎng)格間距要求確定系數(shù)P和Q指定值基本思路：以壁面線處為例

網(wǎng)格線正交指定法向網(wǎng)格間距計(jì)算第2層網(wǎng)格線上的坐標(biāo)通過差分計(jì)算邊界處的如需要利用的信息，可用上一迭代步的值計(jì)算出邊界處的P,Q根據(jù)指數(shù)衰減原則給出全場的P,Q具體公式見：傅德薰《計(jì)算空氣動(dòng)力學(xué)》284-28629習(xí)題9.1網(wǎng)格生成

通過解橢圓型方程生成NACA0