第章有限元法課件

第3章有限元法基礎(chǔ)3.6有限元方程組的求解用有限元方法求解電磁場問題，在導(dǎo)出有限元方程組后，即歸結(jié)為求解一個(gè)方程組的問題，方程組階數(shù)等于未知節(jié)點(diǎn)參數(shù)的節(jié)點(diǎn)數(shù)。γ＝常數(shù)時(shí)：線性方程組γ≠常數(shù)時(shí)：非線性方程組3.6.1有限元線性方程組解法高斯-若爾當(dāng)解法注意有限元方程組系數(shù)矩陣的特點(diǎn)：(1)對稱性(2)稀疏性與正定性第3章有限元法基礎(chǔ)3.6有限元方程組的求解用有限1第3章有限元法基礎(chǔ)下面主要介紹牛頓-拉夫遜迭代法。3.6.2有限元非線性方程組解法

線性化方法：如逐次線性化方法，牛頓-拉夫遜迭代法及改進(jìn)型的牛頓-拉夫遜迭代法等。通過求解目標(biāo)函數(shù)的極小值點(diǎn)來獲得非線性方程組的解：如最速下降法，共軛梯度法等。

逐次線性化方法特點(diǎn)：比較簡便，有良好的收斂性，但收斂速度慢，適宜在求解方程的階數(shù)不是特別高時(shí)采用；

牛頓-拉夫遜迭代法特點(diǎn)：收斂速度快，但在形成方程組時(shí)需要很大的計(jì)算量，并要求有很好的初值，否則迭代過程可能不收斂。為了克服這個(gè)缺點(diǎn)，但又能發(fā)揮其特點(diǎn)，產(chǎn)生了多種改進(jìn)型的牛頓-拉夫遜迭代法。第3章有限元法基礎(chǔ)下面主要介紹牛頓-拉夫遜迭代法。3.621.牛頓-拉夫遜迭代法第3章有限元法基礎(chǔ)對矢量磁位A，使用直角坐標(biāo)系，一般認(rèn)為A是近似解，對于真解必有一個(gè)不為零的余矢量在此方程中，K和F都是A的函數(shù)，如果A(0)是的近似解，將F(A)用泰勒級數(shù)表示，取一次項(xiàng)，則有若用J(A(0))表示F(A)的雅可比矩陣在A(0)處的值，則有令此余矢量對于A(1)解時(shí)為0，有1.牛頓-拉夫遜迭代法第3章有限元法基礎(chǔ)對矢量磁位3第3章有限元法基礎(chǔ)式中求余矢量的雅可比矩陣，需要從一個(gè)單元的余矢量函數(shù)討論起。一個(gè)單元的余矢量為：第3章有限元法基礎(chǔ)式中求余矢量的雅可比矩陣，需要從一個(gè)單4第3章有限元法基礎(chǔ)單元雅可比矩陣為：以為例，看雅可比矩陣第ii項(xiàng)元素對單元e的貢獻(xiàn)。K是A的函數(shù)，，對A的偏導(dǎo)數(shù)均不為零，由(3.66)式可知因此有：第3章有限元法基礎(chǔ)單元雅可比矩陣為：以5第3章有限元法基礎(chǔ)下面求。由以下各式：第3章有限元法基礎(chǔ)下面求。由以下各式6第3章有限元法基礎(chǔ)可以求出：式中所以第3章有限元法基礎(chǔ)可以求出：式中所以7同理可求出：單元雅可比矩陣元素統(tǒng)一表示為：（γ、B

對各單元均不相同）求出每個(gè)單元的雅可比矩陣J(e)和余矢量F(A(k))后，將所有單元的分析結(jié)果進(jìn)行總體合成，得到：則牛頓-拉夫遜迭代公式可以歸結(jié)為下式：第3章有限元法基礎(chǔ)同理可求出：單元雅可比矩陣元素統(tǒng)一表示為：（γ、B對各單元8第3章有限元法基礎(chǔ)（3.185)在上式中，對于J(e)

中各元素，必須在求出后才能進(jìn)行計(jì)算，因此，由矢量磁位A求解時(shí)，由每次迭代解A(K)計(jì)算出單元磁通密度B，由B值查磁化曲線得對應(yīng)H值，再計(jì)算，磁化曲線的查取歸結(jié)為數(shù)值計(jì)算。對B=f(H)的數(shù)學(xué)處理，一般包括插值法和擬合法，比較簡單而有效的辦法是線性插值和拉格朗日插值等。當(dāng)磁場磁化曲線采用逐段線性插值函數(shù)逼近，H值在HK+1與HK之間時(shí)，有第3章有限元法基礎(chǔ)（3.185)在上式中，對于J9第3章有限元法基礎(chǔ)2.改進(jìn)型的牛頓-拉夫遜迭代法牛頓-拉夫遜迭代法的優(yōu)點(diǎn)：（1）收斂速度快，按平方律收斂，每經(jīng)一次迭代有效位數(shù)基本上增加一倍；（2）自校正功能，即A(K+1)僅依賴于F(A)與A(K)，前面迭代的舍入誤差不會一步步傳遞下去。其缺點(diǎn)：（1）每次迭代都要形成一次J(K)，而J(K)的計(jì)算時(shí)間往往比增加迭代的次數(shù)所需要的時(shí)間多；（2）對初值要求較高，選擇不當(dāng)，常會引起振蕩。為了克服所存在的缺點(diǎn)，現(xiàn)已形成多種改進(jìn)型牛頓-拉夫遜迭代法。第3章有限元法基礎(chǔ)2.改進(jìn)型的牛頓-拉夫遜迭代法牛頓10第3章有限元法基礎(chǔ)(1)修正的牛頓—拉夫遜迭代法為了克服每次迭代都需要計(jì)算J

(K)，采取全部求解過程都用J

(0)，即每次迭代不需要形成新的系數(shù)矩陣，保持其斜率不變。這將使迭代次數(shù)增加而影響收斂速度，但如初值選取較好，總體上將節(jié)省運(yùn)算時(shí)間。(2)采用欠松弛因子的方法欠松弛迭代的方法為式中ω為收斂因子，0<ω<1；ΔA(K)為每次迭代近似解的誤差。此外，還有牛頓-拉夫遜迭代法與最速下降法的組合求解法，加速的全牛頓法，采用阻尼因子等，具體方法請參閱文獻(xiàn)［2］。第3章有限元法基礎(chǔ)(1)修正的牛頓—拉夫遜迭代法11否讀入所需系數(shù)及數(shù)組建立右端矢量R建立雅可比矩陣及右端剩余矢量ΔR計(jì)算迭代誤差計(jì)算γ的初值計(jì)算各單元b，c，Δ是開始牛頓－拉夫遜迭代迭代誤差是否小于控制值？結(jié)束非線性方程組計(jì)算框圖：根據(jù)邊界條件修正方程組計(jì)算各單元A，B解方程組得到各節(jié)點(diǎn)A計(jì)算各單元γ及否讀入所需系數(shù)及數(shù)組建立右端矢量R建立雅可比矩陣及右端剩余123.8有限元素的自動(dòng)剖分采取自動(dòng)剖分的必要性：第3章有限元法基礎(chǔ)

在單元剖分中，為了盡可能的壓縮存儲，減小計(jì)算量，提高精度，必須注意以下問題：

(1)三角形各邊不要相差太懸殊，避免出現(xiàn)尖銳的三角形；(2)一個(gè)節(jié)點(diǎn)周圍，不宜集中過多的三角形單元，為壓縮存儲創(chuàng)造條件；(3)三角形單元內(nèi)物理參數(shù)(γ或μ)變化連續(xù)，即媒質(zhì)交界面應(yīng)與單元的邊界重合；(4)精度要求不同的區(qū)域，元素的密度應(yīng)不同；(5)元素節(jié)點(diǎn)編排應(yīng)規(guī)格化。本節(jié)討論平面內(nèi)的幾種自動(dòng)剖分方法。3.8有限元素的自動(dòng)剖分采取自動(dòng)剖分的必要性：第3章13第3章有限元法基礎(chǔ)3.8.1直線內(nèi)插法適宜于對以直線段為邊界的的場域進(jìn)行自動(dòng)剖分。只需給出x,y方向上兩端點(diǎn)的坐標(biāo)，就可算出所有節(jié)點(diǎn)的坐標(biāo)。1.確定x方向、y方向節(jié)點(diǎn)數(shù)及總節(jié)點(diǎn)數(shù)x方向節(jié)點(diǎn)數(shù)y方向節(jié)點(diǎn)數(shù)節(jié)點(diǎn)總數(shù)2.確定各節(jié)點(diǎn)的坐標(biāo)x和y方向節(jié)點(diǎn)坐標(biāo)最小值為x1、y1，最大值為xm、ym，則節(jié)點(diǎn)坐標(biāo)增量分別為，第3章有限元法基礎(chǔ)3.8.1直線內(nèi)插法適宜于對以14第3章有限元法基礎(chǔ)可推得第Ni列Nj行的第Nk個(gè)節(jié)點(diǎn)的編號為式中第Nk個(gè)節(jié)點(diǎn)的坐標(biāo)值為例3.1設(shè)x1=0,y1=0,xm=7cm,ym=5cm，取n=7,m=5,可以求出x方向節(jié)點(diǎn)數(shù)Nx=7+1=8，y方向節(jié)點(diǎn)數(shù)Ny=5+1=6，總節(jié)點(diǎn)數(shù)N0=8×6=48。求第2列第3行節(jié)點(diǎn)數(shù)編號及坐標(biāo)。節(jié)點(diǎn)編號(Ni=2,Nj=3)，按(3.192)式，有Nk=(2-1)6+3=9其坐標(biāo)按(3.193)式為（3.193)第3章有限元法基礎(chǔ)可推得第Ni列Nj行的第Nk153.確定三角形單元編號及各三角形單元3頂點(diǎn)(節(jié)點(diǎn))按規(guī)定順序的的總體編號第3章有限元法基礎(chǔ)對于圖(a)中類型的三角形單元，其單元編號E與行、列的關(guān)系為

E＝2(Ni-1)(Ny-1)+Nj

式中Ni＝1，2，3…，(Nx-1);

Nj＝1，2，3…，(Ny-1)。圖(a)中第E個(gè)三角形單元的3節(jié)點(diǎn)，i,j,m的總體編號與行、列的關(guān)系分別為同理，對圖(b)類型的三角形單元有3.確定三角形單元編號及各三角形單元3頂點(diǎn)(節(jié)點(diǎn))按規(guī)定順序16第3章有限元法基礎(chǔ)對于圖(b)中類型的三角形單元，其單元編號E與行、列的關(guān)系為

E＝2(Ni-1)(Ny-1)+Nj

式中Ni＝1，2，3…，(Nx-1);Nj＝1，2，3…，(Ny-1)。3節(jié)點(diǎn)i,j,m的總體編號與行、列的關(guān)系分別為例3.2(1)求上例中1列3行（b）類型單元的編號E和3個(gè)節(jié)點(diǎn)i,j,m的總體節(jié)點(diǎn)編號。(2)求出1列3行（a）類型單元的編號E的3個(gè)節(jié)點(diǎn)i,j,m的總體節(jié)點(diǎn)編號。第3章有限元法基礎(chǔ)對于圖(b)中類型的三角形單元，其單17當(dāng)場域形狀復(fù)雜，要求場域內(nèi)單元密度分布不均勻，又希望密度是平滑過渡時(shí)，這種簡單的直線內(nèi)插法就不適用了，廣泛應(yīng)用的是等勢剖分的方法。第3章有限元法基礎(chǔ)3.8.2等勢剖分1962年，W.P.Crowley1.微分方程的推導(dǎo)等勢剖分是把網(wǎng)格線看成平面上的兩個(gè)勢函數(shù)的等勢線，定義兩個(gè)函數(shù)

(x,y)及ψ(x,y)，它們在整個(gè)計(jì)算區(qū)域內(nèi)部滿足拉普拉斯方程：(3.198)當(dāng)場域形狀復(fù)雜，要求場域內(nèi)單元密度分布不均勻，又希望密度18第3章有限元法基礎(chǔ)當(dāng)邊界和ψs給定后，可以解出場域內(nèi)和ψ的分布，求出等勢線。兩族曲線的交點(diǎn)就是剖分所要確定的節(jié)點(diǎn)，但是得到的節(jié)點(diǎn)坐標(biāo)是

和ψ值，我們需要的是等勢線交點(diǎn)在(x,y)平面上的坐標(biāo)！若將x,y看做(,ψ)平面上的函數(shù)x(,ψ)及y(,ψ)，通過函數(shù)的變換，方程式(3.198)可化為所要求解的形式，簡單步驟如下：設(shè)兩端對x求導(dǎo)第3章有限元法基礎(chǔ)當(dāng)邊界和ψs給定后，可以解19第3章有限元法基礎(chǔ)可解得令可得(3.204)第3章有限元法基礎(chǔ)可解得令可得(3.204)20第3章有限元法基礎(chǔ)又因?yàn)閮啥藢求導(dǎo)將式(3.204)兩式展開，與上式比較可得(3.208)(3.209)(3.210)(3.211)第3章有限元法基礎(chǔ)又因?yàn)閮啥藢求導(dǎo)將式(3.20421第3章有限元法基礎(chǔ)將式(3.208)、(3.209)代入，則有其中由式可得(3.213)(3.212)(3.214)第3章有限元法基礎(chǔ)將式(3.208)、(3.209)代22將式(3.213)中展開為第3章有限元法基礎(chǔ)將式(3.203)、(3.208)、(3.209)、(3.210)、(3.211)、(3.213)、(3.214)及表達(dá)式代入式(3.212)，化簡可得(3.215)同理可以導(dǎo)出(3.216)(3.217)將式(3.213)中展開為第3章23第3章有限元法基礎(chǔ)其中齊次方程式(3.216)、式(3.217)有唯一解的條件為J≠0，因此有若在(,ψ)平面得到式(3.219)的兩組數(shù)值解x(,ψ)和y(,ψ)，則對于定域內(nèi)每一個(gè)確定的、ψi，就有一對xi和yi，這就是要求的某兩條網(wǎng)格線相交處的節(jié)點(diǎn)坐標(biāo)。(3.218)(3.219)第3章有限元法基礎(chǔ)其中齊次方程式(3.216)、式(3.24在(,ψ)平面上網(wǎng)格與節(jié)點(diǎn)數(shù)與(x,y)平面上是相等的，但在(,ψ)平面上是等距的網(wǎng)格，在(x,y)平面上卻是用相應(yīng)邊界確定的，各單元密度過渡是平滑的。但并不是等距網(wǎng)格。第3章有限元法基礎(chǔ)2.等勢剖分計(jì)算公式的推導(dǎo)圖3.19網(wǎng)格與邊界的劃分xy平面區(qū)域邏輯網(wǎng)格（1）邊界值的確定在(,ψ)平面上網(wǎng)格與節(jié)點(diǎn)數(shù)與(x,y)平面上25第3章有限元法基礎(chǔ)（2）三角形網(wǎng)格的邏輯坐標(biāo)在空間坐標(biāo)中用三角形網(wǎng)格剖分時(shí)，在(,ψ)坐標(biāo)中，等線和等ψ線有多種取法，常用的有如圖3.20所示3種。圖3.20不同的邏輯網(wǎng)格方式(a)(b)(c)（3）對等邊三角形邏輯網(wǎng)格求(3.219)式的差分方程第3章有限元法基礎(chǔ)（2）三角形網(wǎng)格的邏輯坐標(biāo)在空間坐標(biāo)中26第3章有限元法基礎(chǔ)令等邊三角形,，推導(dǎo)的關(guān)鍵是對三角形網(wǎng)格求出，，，，，。利用積分中值定理有設(shè)函數(shù)格林公式同理，有第3章有限元法基礎(chǔ)令等邊三角形,27第3章有限元法基礎(chǔ)對于二次積分，同樣應(yīng)用積分中值定理和格林公式，可得現(xiàn)把替換成x(,ψ)或y(,ψ)，可以計(jì)算出每個(gè)節(jié)點(diǎn)(l,m)的，，與。第3章有限元法基礎(chǔ)對于二次積分，同樣應(yīng)用積分中值定理和格28第3章有限元法基礎(chǔ)圖3.21等邊三角形節(jié)點(diǎn)編號如圖示，節(jié)點(diǎn)(l,m)由6個(gè)等邊三角形包圍，求的回路積分是沿6個(gè)三角形的外邊界積分，即第3章有限元法基礎(chǔ)圖3.21等邊三角形節(jié)點(diǎn)編號如圖示29第3章有限元法基礎(chǔ)同理可得，與上面求法完全相同，即將xi換成yi以上式子對y全部成立。由此可進(jìn)一步求出α，β，γ。為求出，，，，，，要利用圖3.21中虛線網(wǎng)格，即對質(zhì)心構(gòu)成包圍(l,m)節(jié)點(diǎn)的六邊形進(jìn)行積分。第3章有限元法基礎(chǔ)同理可得，與上面求30第3章有限元法基礎(chǔ)首先求出各三角形內(nèi)平均的和。利用上式，求出圖3.21式中對應(yīng)三角形的。第3章有限元法基礎(chǔ)首先求出各三角形內(nèi)平均的31第3章有限元法基礎(chǔ)再沿圖中虛線做環(huán)路積分，求出對應(yīng)圖示各點(diǎn)(m,l)：1’[(m+1/3),(l+1/3)]I[m,(l+1/2)]2’[(m-1/3),(l+1/3)]II[(m-1/2),(l+1/2)]3’[(m-2/3),(l+1/3)]III[(m-1/2),l]4’[(m-1/3),(l-1/3)]IV[m,(l-1/2)]5’[(m+1/3),(l-2/3)]V[(m+1/2),(l-1/2)]6’[(m+2/3),(l-1/3)]VI[(m+1/2),l]同理可得第3章有限元法基礎(chǔ)再沿圖中虛線做環(huán)路積分，求出對應(yīng)圖示各32第3章有限元法基礎(chǔ)同理可得到，，類似的關(guān)系。將它們代入式(3.219)：整理后可得第3章有限元法基礎(chǔ)同理可得到，33上式整理后可得第3章有限元法基礎(chǔ)歸納成統(tǒng)一形式為：式中wi是α，β，和γ的函數(shù)。上面得到的值是對等邊三角形的情況。對于不同的邏輯網(wǎng)格，它們?nèi)≈祵⒉煌?3.236)上式整理后可得第3章有限元法基礎(chǔ)歸納成統(tǒng)一形式為：式中343.差分方程組的求解第3章有限元法基礎(chǔ)方程式(3.236)是非線性代數(shù)方程，n個(gè)內(nèi)部節(jié)點(diǎn)有2n個(gè)坐標(biāo)，有2n個(gè)方程，一般采用高斯-賽德爾選代。計(jì)算的步驟為：(1)邊界頂點(diǎn)坐標(biāo)已知，其它節(jié)點(diǎn)賦初值。(2)從左下角開始計(jì)算J，α，β

，和γ。(3)計(jì)算(3.236)式，求x,y，經(jīng)過多次迭代，最后得到各節(jié)點(diǎn)坐標(biāo)值。3.差分方程組的求解第3章有限元法基礎(chǔ)方程式(3.235

網(wǎng)格自動(dòng)剖分在有限元方法中占有很重要的地位，它直接決定了生成矩陣的特性。目前二維有限元網(wǎng)格剖分算法已經(jīng)趨于成熟和完善，幾乎可以處理任意復(fù)雜的場域，已不需要一般性工作。但網(wǎng)格的自適應(yīng)細(xì)分仍然是一個(gè)熱門課題，目的是使網(wǎng)格分別自動(dòng)適應(yīng)于場域結(jié)構(gòu)或場量分布，使場域中每個(gè)單元都能給出幾乎相同的計(jì)算精度，這樣就要求程序本身能自行判斷何處的單元需要細(xì)分，細(xì)分到何種程度。第3章有限元法基礎(chǔ)3.8.3自適應(yīng)剖分技術(shù)網(wǎng)格自動(dòng)剖分在有限元方法中占有很重要的地位，36第3章有限元法基礎(chǔ)自適應(yīng)技術(shù)是靠網(wǎng)格細(xì)分與場量計(jì)算的循環(huán)過程來實(shí)現(xiàn)，具體步驟為：(1)生成開端網(wǎng)格，即采用網(wǎng)格自動(dòng)細(xì)分軟件將場域剖成很粗的網(wǎng)格；(2)求解場量，即形成有限元方程與求解方程；(3)分析場結(jié)果計(jì)算誤差(對每個(gè)單元進(jìn)行誤差分析)；(4)根據(jù)誤差分析，確定需要細(xì)分的網(wǎng)格單元；(5)細(xì)分局部網(wǎng)格；(6)返回到(2)，求解場量；(7)計(jì)算最后結(jié)果。第3章有限元法基礎(chǔ)自適應(yīng)技術(shù)是靠網(wǎng)格細(xì)分與場量計(jì)算的37第3章有限元法基礎(chǔ)否場域信息輸入是精度滿足要求？開端網(wǎng)格生成場量計(jì)算誤差分析網(wǎng)格自動(dòng)細(xì)分后處理圖3.22自適應(yīng)剖分框圖第3章有限元法基礎(chǔ)否場域信息輸入是精度滿足要求？開端網(wǎng)格38網(wǎng)格初始剖分采用自適應(yīng)技術(shù)后網(wǎng)格一次細(xì)分網(wǎng)格二次細(xì)分第3章有限元法基礎(chǔ)Matlab提供PDETOOL軟件包就具有網(wǎng)格自動(dòng)剖分和自適應(yīng)細(xì)分功能，并且它采用CAD輸入方式，既簡單又直觀，所以可選用它作為有限元計(jì)算的前處理程序。網(wǎng)格初始剖分采用自適應(yīng)技術(shù)后網(wǎng)格一次細(xì)分網(wǎng)格二次細(xì)分第39第章有限元法課件40第章有限元法課件41第章有限元法課件42第章有限元法課件43第3章有限元法基礎(chǔ)3.6有限元方程組的求解用有限元方法求解電磁場問題，在導(dǎo)出有限元方程組后，即歸結(jié)為求解一個(gè)方程組的問題，方程組階數(shù)等于未知節(jié)點(diǎn)參數(shù)的節(jié)點(diǎn)數(shù)。γ＝常數(shù)時(shí)：線性方程組γ≠常數(shù)時(shí)：非線性方程組3.6.1有限元線性方程組解法高斯-若爾當(dāng)解法注意有限元方程組系數(shù)矩陣的特點(diǎn)：(1)對稱性(2)稀疏性與正定性第3章有限元法基礎(chǔ)3.6有限元方程組的求解用有限44第3章有限元法基礎(chǔ)下面主要介紹牛頓-拉夫遜迭代法。3.6.2有限元非線性方程組解法

線性化方法：如逐次線性化方法，牛頓-拉夫遜迭代法及改進(jìn)型的牛頓-拉夫遜迭代法等。通過求解目標(biāo)函數(shù)的極小值點(diǎn)來獲得非線性方程組的解：如最速下降法，共軛梯度法等。

逐次線性化方法特點(diǎn)：比較簡便，有良好的收斂性，但收斂速度慢，適宜在求解方程的階數(shù)不是特別高時(shí)采用；

牛頓-拉夫遜迭代法特點(diǎn)：收斂速度快，但在形成方程組時(shí)需要很大的計(jì)算量，并要求有很好的初值，否則迭代過程可能不收斂。為了克服這個(gè)缺點(diǎn)，但又能發(fā)揮其特點(diǎn)，產(chǎn)生了多種改進(jìn)型的牛頓-拉夫遜迭代法。第3章有限元法基礎(chǔ)下面主要介紹牛頓-拉夫遜迭代法。3.6451.牛頓-拉夫遜迭代法第3章有限元法基礎(chǔ)對矢量磁位A，使用直角坐標(biāo)系，一般認(rèn)為A是近似解，對于真解必有一個(gè)不為零的余矢量在此方程中，K和F都是A的函數(shù)，如果A(0)是的近似解，將F(A)用泰勒級數(shù)表示，取一次項(xiàng)，則有若用J(A(0))表示F(A)的雅可比矩陣在A(0)處的值，則有令此余矢量對于A(1)解時(shí)為0，有1.牛頓-拉夫遜迭代法第3章有限元法基礎(chǔ)對矢量磁位46第3章有限元法基礎(chǔ)式中求余矢量的雅可比矩陣，需要從一個(gè)單元的余矢量函數(shù)討論起。一個(gè)單元的余矢量為：第3章有限元法基礎(chǔ)式中求余矢量的雅可比矩陣，需要從一個(gè)單47第3章有限元法基礎(chǔ)單元雅可比矩陣為：以為例，看雅可比矩陣第ii項(xiàng)元素對單元e的貢獻(xiàn)。K是A的函數(shù)，，對A的偏導(dǎo)數(shù)均不為零，由(3.66)式可知因此有：第3章有限元法基礎(chǔ)單元雅可比矩陣為：以48第3章有限元法基礎(chǔ)下面求。由以下各式：第3章有限元法基礎(chǔ)下面求。由以下各式49第3章有限元法基礎(chǔ)可以求出：式中所以第3章有限元法基礎(chǔ)可以求出：式中所以50同理可求出：單元雅可比矩陣元素統(tǒng)一表示為：（γ、B

對各單元均不相同）求出每個(gè)單元的雅可比矩陣J(e)和余矢量F(A(k))后，將所有單元的分析結(jié)果進(jìn)行總體合成，得到：則牛頓-拉夫遜迭代公式可以歸結(jié)為下式：第3章有限元法基礎(chǔ)同理可求出：單元雅可比矩陣元素統(tǒng)一表示為：（γ、B對各單元51第3章有限元法基礎(chǔ)（3.185)在上式中，對于J(e)

中各元素，必須在求出后才能進(jìn)行計(jì)算，因此，由矢量磁位A求解時(shí)，由每次迭代解A(K)計(jì)算出單元磁通密度B，由B值查磁化曲線得對應(yīng)H值，再計(jì)算，磁化曲線的查取歸結(jié)為數(shù)值計(jì)算。對B=f(H)的數(shù)學(xué)處理，一般包括插值法和擬合法，比較簡單而有效的辦法是線性插值和拉格朗日插值等。當(dāng)磁場磁化曲線采用逐段線性插值函數(shù)逼近，H值在HK+1與HK之間時(shí)，有第3章有限元法基礎(chǔ)（3.185)在上式中，對于J52第3章有限元法基礎(chǔ)2.改進(jìn)型的牛頓-拉夫遜迭代法牛頓-拉夫遜迭代法的優(yōu)點(diǎn)：（1）收斂速度快，按平方律收斂，每經(jīng)一次迭代有效位數(shù)基本上增加一倍；（2）自校正功能，即A(K+1)僅依賴于F(A)與A(K)，前面迭代的舍入誤差不會一步步傳遞下去。其缺點(diǎn)：（1）每次迭代都要形成一次J(K)，而J(K)的計(jì)算時(shí)間往往比增加迭代的次數(shù)所需要的時(shí)間多；（2）對初值要求較高，選擇不當(dāng)，常會引起振蕩。為了克服所存在的缺點(diǎn)，現(xiàn)已形成多種改進(jìn)型牛頓-拉夫遜迭代法。第3章有限元法基礎(chǔ)2.改進(jìn)型的牛頓-拉夫遜迭代法牛頓53第3章有限元法基礎(chǔ)(1)修正的牛頓—拉夫遜迭代法為了克服每次迭代都需要計(jì)算J

(K)，采取全部求解過程都用J

(0)，即每次迭代不需要形成新的系數(shù)矩陣，保持其斜率不變。這將使迭代次數(shù)增加而影響收斂速度，但如初值選取較好，總體上將節(jié)省運(yùn)算時(shí)間。(2)采用欠松弛因子的方法欠松弛迭代的方法為式中ω為收斂因子，0<ω<1；ΔA(K)為每次迭代近似解的誤差。此外，還有牛頓-拉夫遜迭代法與最速下降法的組合求解法，加速的全牛頓法，采用阻尼因子等，具體方法請參閱文獻(xiàn)［2］。第3章有限元法基礎(chǔ)(1)修正的牛頓—拉夫遜迭代法54否讀入所需系數(shù)及數(shù)組建立右端矢量R建立雅可比矩陣及右端剩余矢量ΔR計(jì)算迭代誤差計(jì)算γ的初值計(jì)算各單元b，c，Δ是開始牛頓－拉夫遜迭代迭代誤差是否小于控制值？結(jié)束非線性方程組計(jì)算框圖：根據(jù)邊界條件修正方程組計(jì)算各單元A，B解方程組得到各節(jié)點(diǎn)A計(jì)算各單元γ及否讀入所需系數(shù)及數(shù)組建立右端矢量R建立雅可比矩陣及右端剩余553.8有限元素的自動(dòng)剖分采取自動(dòng)剖分的必要性：第3章有限元法基礎(chǔ)

在單元剖分中，為了盡可能的壓縮存儲，減小計(jì)算量，提高精度，必須注意以下問題：

(1)三角形各邊不要相差太懸殊，避免出現(xiàn)尖銳的三角形；(2)一個(gè)節(jié)點(diǎn)周圍，不宜集中過多的三角形單元，為壓縮存儲創(chuàng)造條件；(3)三角形單元內(nèi)物理參數(shù)(γ或μ)變化連續(xù)，即媒質(zhì)交界面應(yīng)與單元的邊界重合；(4)精度要求不同的區(qū)域，元素的密度應(yīng)不同；(5)元素節(jié)點(diǎn)編排應(yīng)規(guī)格化。本節(jié)討論平面內(nèi)的幾種自動(dòng)剖分方法。3.8有限元素的自動(dòng)剖分采取自動(dòng)剖分的必要性：第3章56第3章有限元法基礎(chǔ)3.8.1直線內(nèi)插法適宜于對以直線段為邊界的的場域進(jìn)行自動(dòng)剖分。只需給出x,y方向上兩端點(diǎn)的坐標(biāo)，就可算出所有節(jié)點(diǎn)的坐標(biāo)。1.確定x方向、y方向節(jié)點(diǎn)數(shù)及總節(jié)點(diǎn)數(shù)x方向節(jié)點(diǎn)數(shù)y方向節(jié)點(diǎn)數(shù)節(jié)點(diǎn)總數(shù)2.確定各節(jié)點(diǎn)的坐標(biāo)x和y方向節(jié)點(diǎn)坐標(biāo)最小值為x1、y1，最大值為xm、ym，則節(jié)點(diǎn)坐標(biāo)增量分別為，第3章有限元法基礎(chǔ)3.8.1直線內(nèi)插法適宜于對以57第3章有限元法基礎(chǔ)可推得第Ni列Nj行的第Nk個(gè)節(jié)點(diǎn)的編號為式中第Nk個(gè)節(jié)點(diǎn)的坐標(biāo)值為例3.1設(shè)x1=0,y1=0,xm=7cm,ym=5cm，取n=7,m=5,可以求出x方向節(jié)點(diǎn)數(shù)Nx=7+1=8，y方向節(jié)點(diǎn)數(shù)Ny=5+1=6，總節(jié)點(diǎn)數(shù)N0=8×6=48。求第2列第3行節(jié)點(diǎn)數(shù)編號及坐標(biāo)。節(jié)點(diǎn)編號(Ni=2,Nj=3)，按(3.192)式，有Nk=(2-1)6+3=9其坐標(biāo)按(3.193)式為（3.193)第3章有限元法基礎(chǔ)可推得第Ni列Nj行的第Nk583.確定三角形單元編號及各三角形單元3頂點(diǎn)(節(jié)點(diǎn))按規(guī)定順序的的總體編號第3章有限元法基礎(chǔ)對于圖(a)中類型的三角形單元，其單元編號E與行、列的關(guān)系為

E＝2(Ni-1)(Ny-1)+Nj

式中Ni＝1，2，3…，(Nx-1);

Nj＝1，2，3…，(Ny-1)。圖(a)中第E個(gè)三角形單元的3節(jié)點(diǎn)，i,j,m的總體編號與行、列的關(guān)系分別為同理，對圖(b)類型的三角形單元有3.確定三角形單元編號及各三角形單元3頂點(diǎn)(節(jié)點(diǎn))按規(guī)定順序59第3章有限元法基礎(chǔ)對于圖(b)中類型的三角形單元，其單元編號E與行、列的關(guān)系為

E＝2(Ni-1)(Ny-1)+Nj

式中Ni＝1，2，3…，(Nx-1);Nj＝1，2，3…，(Ny-1)。3節(jié)點(diǎn)i,j,m的總體編號與行、列的關(guān)系分別為例3.2(1)求上例中1列3行（b）類型單元的編號E和3個(gè)節(jié)點(diǎn)i,j,m的總體節(jié)點(diǎn)編號。(2)求出1列3行（a）類型單元的編號E的3個(gè)節(jié)點(diǎn)i,j,m的總體節(jié)點(diǎn)編號。第3章有限元法基礎(chǔ)對于圖(b)中類型的三角形單元，其單60當(dāng)場域形狀復(fù)雜，要求場域內(nèi)單元密度分布不均勻，又希望密度是平滑過渡時(shí)，這種簡單的直線內(nèi)插法就不適用了，廣泛應(yīng)用的是等勢剖分的方法。第3章有限元法基礎(chǔ)3.8.2等勢剖分1962年，W.P.Crowley1.微分方程的推導(dǎo)等勢剖分是把網(wǎng)格線看成平面上的兩個(gè)勢函數(shù)的等勢線，定義兩個(gè)函數(shù)

(x,y)及ψ(x,y)，它們在整個(gè)計(jì)算區(qū)域內(nèi)部滿足拉普拉斯方程：(3.198)當(dāng)場域形狀復(fù)雜，要求場域內(nèi)單元密度分布不均勻，又希望密度61第3章有限元法基礎(chǔ)當(dāng)邊界和ψs給定后，可以解出場域內(nèi)和ψ的分布，求出等勢線。兩族曲線的交點(diǎn)就是剖分所要確定的節(jié)點(diǎn)，但是得到的節(jié)點(diǎn)坐標(biāo)是

和ψ值，我們需要的是等勢線交點(diǎn)在(x,y)平面上的坐標(biāo)！若將x,y看做(,ψ)平面上的函數(shù)x(,ψ)及y(,ψ)，通過函數(shù)的變換，方程式(3.198)可化為所要求解的形式，簡單步驟如下：設(shè)兩端對x求導(dǎo)第3章有限元法基礎(chǔ)當(dāng)邊界和ψs給定后，可以解62第3章有限元法基礎(chǔ)可解得令可得(3.204)第3章有限元法基礎(chǔ)可解得令可得(3.204)63第3章有限元法基礎(chǔ)又因?yàn)閮啥藢求導(dǎo)將式(3.204)兩式展開，與上式比較可得(3.208)(3.209)(3.210)(3.211)第3章有限元法基礎(chǔ)又因?yàn)閮啥藢求導(dǎo)將式(3.20464第3章有限元法基礎(chǔ)將式(3.208)、(3.209)代入，則有其中由式可得(3.213)(3.212)(3.214)第3章有限元法基礎(chǔ)將式(3.208)、(3.209)代65將式(3.213)中展開為第3章有限元法基礎(chǔ)將式(3.203)、(3.208)、(3.209)、(3.210)、(3.211)、(3.213)、(3.214)及表達(dá)式代入式(3.212)，化簡可得(3.215)同理可以導(dǎo)出(3.216)(3.217)將式(3.213)中展開為第3章66第3章有限元法基礎(chǔ)其中齊次方程式(3.216)、式(3.217)有唯一解的條件為J≠0，因此有若在(,ψ)平面得到式(3.219)的兩組數(shù)值解x(,ψ)和y(,ψ)，則對于定域內(nèi)每一個(gè)確定的、ψi，就有一對xi和yi，這就是要求的某兩條網(wǎng)格線相交處的節(jié)點(diǎn)坐標(biāo)。(3.218)(3.219)第3章有限元法基礎(chǔ)其中齊次方程式(3.216)、式(3.67在(,ψ)平面上網(wǎng)格與節(jié)點(diǎn)數(shù)與(x,y)平面上是相等的，但在(,ψ)平面上是等距的網(wǎng)格，在(x,y)平面上卻是用相應(yīng)邊界確定的，各單元密度過渡是平滑的。但并不是等距網(wǎng)格。第3章有限元法基礎(chǔ)2.等勢剖分計(jì)算公式的推導(dǎo)圖3.19網(wǎng)格與邊界的劃分xy平面區(qū)域邏輯網(wǎng)格（1）邊界值的確定在(,ψ)平面上網(wǎng)格與節(jié)點(diǎn)數(shù)與(x,y)平面上68第3章有限元法基礎(chǔ)（2）三角形網(wǎng)格的邏輯坐標(biāo)在空間坐標(biāo)中用三角形網(wǎng)格剖分時(shí)，在(,ψ)坐標(biāo)中，等線和等ψ線有多種取法，常用的有如圖3.20所示3種。圖3.20不同的邏輯網(wǎng)格方式(a)(b)(c)（3）對等邊三角形邏輯網(wǎng)格求(3.219)式的差分方程第3章有限元法基礎(chǔ)（2）三角形網(wǎng)格的邏輯坐標(biāo)在空間坐標(biāo)中69第3章有限元法基礎(chǔ)令等邊三角形,，推導(dǎo)的關(guān)鍵是對三角形網(wǎng)格求出，，，，，。利用積分中值定理有設(shè)函數(shù)格林公式同理，有第3章有限元法基礎(chǔ)令等邊三角形,70第3章有限元法基礎(chǔ)對于二次積分，同樣應(yīng)用積分中值定理和格林公式，可得現(xiàn)把替換成x(,ψ)或y(,ψ)，可以計(jì)算出每個(gè)節(jié)點(diǎn)(l,m)的，，與。第3章有限元法基礎(chǔ)對于二次積分，同樣應(yīng)用積分中值定理和格71第3章有限元法基礎(chǔ)圖3.21等邊三角形節(jié)點(diǎn)編號如圖示，節(jié)點(diǎn)(l,m)由6個(gè)等邊三角形包圍，求的回路積分是沿6個(gè)三角形的外邊界積分，即第3章有限元法基礎(chǔ)圖3.21等邊三角形節(jié)點(diǎn)編號如圖示72第3章有限元法基礎(chǔ)同理可得，與上面求法完全相同，即將xi換成yi以上式子對y全部成立。由此可進(jìn)一步求出α，β，γ。為求出，，，，，，要利用圖3.21中虛線網(wǎng)格，即對質(zhì)心構(gòu)成包圍(l,m)節(jié)點(diǎn)的六邊形進(jìn)行積分。第3章有限元法基礎(chǔ)同理可得，與上面求73第3章有限元法基礎(chǔ)首先求出各三角形內(nèi)平均的和。利用上式，求出圖3.21式中對應(yīng)三角形的。第3章有限元法基礎(chǔ)首先求出各三角形內(nèi)平均的74第3章有限元法基礎(chǔ)再沿圖中虛線做環(huán)路積分，求出對應(yīng)圖示各點(diǎn)(m,l)：1’[(m+1/3),(l+1/3)]I[m,(l+1/2)]2