平面問題高階單元

平面問題高階單元第1頁，課件共54頁，創(chuàng)作于2023年2月可能達(dá)到較高的精度。圖7-1懸臂梁分別采用高、低階單元計(jì)算就是一個典型的例子。h4hPAB懸臂深梁解析解：

A=1.0

B=1.0常應(yīng)變單元：

A=0.866

B=0.619高階單元：

A=0.99

B=0.99??????圖7-1第2頁，課件共54頁，創(chuàng)作于2023年2月選擇位移模式時，第2章提到要考慮解的收斂性，即要考慮到位移模式的完備性和協(xié)調(diào)性。實(shí)際操作中，一般應(yīng)考慮位移模式的對稱性。這是因?yàn)?，有限元位移模式的選擇實(shí)際是以帕斯卡（Pascal）三解形基礎(chǔ)上的（如圖7-2所示），由低價至高階，順序選取，組成多項(xiàng)式。多項(xiàng)式中的項(xiàng)數(shù)等于單元節(jié)點(diǎn)自由度數(shù)。如三節(jié)點(diǎn)三角形單元，位移模式取完全一次式，共3項(xiàng)。六節(jié)點(diǎn)三角形單元，位移模式取完全二次式共6項(xiàng)。如果某一階次不能全取，則應(yīng)按對稱性原則適當(dāng)選取。

1xyx2xyy2x3x2yxy2y3

x4x3yx2y2xy3y4

圖7-2多項(xiàng)式選擇的怕斯卡三角形第3頁，課件共54頁，創(chuàng)作于2023年2月

1xyx2xyy2x3x2yxy2y3

x4x3yx2y2xy3y4

圖7-2多項(xiàng)式選擇的怕斯卡三角形例如在下節(jié)將要討論的四結(jié)點(diǎn)矩形單元中，位移模式不能取1，x，y，x2四項(xiàng)，也不能取1，x，y，y2四項(xiàng)，而應(yīng)取1，x，y，xy四項(xiàng)。7.2四節(jié)點(diǎn)矩形單元圖7-3示出的矩形單元，邊長分別為2a和2b。取4個角點(diǎn)為節(jié)點(diǎn)，編號為i，j，l，m。將x軸和y軸置于單元的對稱軸上。單元的位移函數(shù)可取為：1、位移函數(shù)第4頁，課件共54頁，創(chuàng)作于2023年2月在上式表示的位移模式中，a1,a2,a3,a5,a6,a7，a8反映了單元的剛體位移和常應(yīng)變。在單元的邊界(x=±a或y=±a)上（或），位移是按線性分布的。因此，相鄰單元在公共邊上的位移是連續(xù)的。這樣，位移模式滿足了解答收斂性的充分條件。ijlmxyaabb圖7-3在式（7-1）中代入節(jié)點(diǎn)位移和節(jié)點(diǎn)坐標(biāo)后，可解出（7-1）第5頁，課件共54頁，創(chuàng)作于2023年2月式中形函數(shù)為：（7-3）（7-2）各待定系數(shù)(a1…a8)。將這些系數(shù)再代入式（7-1），可得：第6頁，課件共54頁，創(chuàng)作于2023年2月則式（7-3）可簡寫為（7-4）將位移函數(shù)寫成矩陣形式，即有與式（2-20）相同的形式（7-5）式中（7-6）令在節(jié)點(diǎn)上的值為：第7頁，課件共54頁，創(chuàng)作于2023年2月（7-7）其中，I為二階單位矩陣。2、應(yīng)變矩陣根據(jù)幾何方程，可得與式（2-25）同樣的形式（7-8）把應(yīng)變矩陣[B]寫成子矩陣形式（7-9）其中（7-10）第8頁，課件共54頁，創(chuàng)作于2023年2月由此可見，[B]是

、

的函數(shù)，即是x、y的函數(shù)。因此單元中的應(yīng)變不再是常數(shù)。3、應(yīng)力矩陣根據(jù)應(yīng)力-應(yīng)變關(guān)系，可以計(jì)算單元中的應(yīng)力，得到式（2-28）同樣形式（7-11）應(yīng)力矩陣[S]具有與式（2-29）同樣形式（7-12）將[S]寫成子矩陣形式（7-13）第9頁，課件共54頁，創(chuàng)作于2023年2月其中（7-14）上式對應(yīng)平面應(yīng)力情形。對于平面應(yīng)變情形，只需將其中的E，

作相應(yīng)的改變即可。4、單元剛度矩陣單元剛度矩陣可采用式（2-33a）進(jìn)行計(jì)算（2-33a）第10頁，課件共54頁，創(chuàng)作于2023年2月在四節(jié)點(diǎn)矩形單元中，[k]是一個8×8的矩陣。將[k]寫成分塊形式：（7-16）其中的子矩陣[krs]2×2可由下式計(jì)算第11頁，課件共54頁，創(chuàng)作于2023年2月（7-17）上式對應(yīng)平面應(yīng)力情形。對于平面應(yīng)變情形，只須將上式中的E、

作相應(yīng)的改變。5、等價節(jié)點(diǎn)力單元體積力和表面力引起的節(jié)點(diǎn)力仍可用式（2-45）和（2-46）進(jìn)行計(jì)算。第12頁，課件共54頁，創(chuàng)作于2023年2月對本問題給定的位移函數(shù)，若體積力是重力的情形（設(shè)重度為

），單元等價節(jié)點(diǎn)載荷列陣為：（2-45）（2-46）（7-18）有了對單元的上述結(jié)果，便可應(yīng)用第5章的方法組集結(jié)構(gòu)剛度矩陣和節(jié)點(diǎn)荷載向量；求解節(jié)點(diǎn)位移；計(jì)算內(nèi)力和應(yīng)力。

四節(jié)點(diǎn)矩形單元采用較高階的位移模式，具有比三節(jié)點(diǎn)三角形單元較高的計(jì)算精度。但矩形單元也有缺點(diǎn)，第13頁，課件共54頁，創(chuàng)作于2023年2月在三角形單元i,j,m的各邊中點(diǎn)增設(shè)一個節(jié)點(diǎn)，使每個單元具有6個節(jié)點(diǎn)，得到圖7-4所示的六節(jié)點(diǎn)三角形單元。這種單元具有12個自由度，可以采用完全二次多項(xiàng)式的位移模式：一是不能適應(yīng)斜線及曲線邊界，二是不便于采用大小不同的單元。7.3六節(jié)點(diǎn)三角形單元1、位移模式

???ijmi

j

m

xy圖7-4第14頁，課件共54頁，創(chuàng)作于2023年2月（7-20）所取位移模式反映了單元的剛體位移和常應(yīng)變；單元內(nèi)部是連續(xù)的；在單元邊界上位移分量按拋物線變化，而每條公共邊界上有3個公共結(jié)點(diǎn)，可以保證相鄰兩單元位移的連續(xù)性。因此，上述位移模式滿足收斂的必要和充分條件。上述位移模式確定之后，可以用分析三節(jié)點(diǎn)三角形單元和四節(jié)點(diǎn)矩形單元相同的方法進(jìn)行分析。得到形函數(shù)、應(yīng)變矩陣、應(yīng)力矩陣、單元剛度矩陣、等價節(jié)點(diǎn)力向量。但其過程十分繁復(fù)，采用面積坐標(biāo)可以大大簡化計(jì)算。第15頁，課件共54頁，創(chuàng)作于2023年2月2、面積坐標(biāo)

對于一個三角形ijm（圖7-5），三角形內(nèi)任一點(diǎn)P（x,y）的位置，可以用如下的三個比值來確定：ijmxy圖7-5·P（7-21）AiAjAm（1）定義第16頁，課件共54頁，創(chuàng)作于2023年2月其中A為三角形ijm的面積，Ai,Aj,Am分別為三角形的Pjm,Pmi,Pijd的面積。這三個比值Li,Lj,Lm稱為P點(diǎn)的面積坐標(biāo)。由于則（7-22）由此可見，P點(diǎn)的三個面積坐標(biāo)不是獨(dú)立的。同時，面積坐標(biāo)只是用以確定三角形內(nèi)部某點(diǎn)的位置，因而是一種局部坐標(biāo)。下面進(jìn)一步給出面積坐標(biāo)的幾個性質(zhì)。（2）面積坐標(biāo)與直角坐標(biāo)的關(guān)系在圖7-5中，三角形Pjm的面積為第17頁，課件共54頁，創(chuàng)作于2023年2月（7-23）由式（7-23），式（7-21）化為（7-24）將式（7-24）、（7-23a）和式（2-18）、（2-17）對比，可知，面積坐標(biāo)就是三節(jié)點(diǎn)三角形單元的形函數(shù)（7-23a）Ni、Nj、Nm。第18頁，課件共54頁，創(chuàng)作于2023年2月將式（7-24）的3個式子分別乘以xi,xj,xm，然后相加，并利用關(guān)系式（7-23a），有同理（7-25）（3）面積坐標(biāo)的導(dǎo)數(shù)公式根據(jù)面積坐標(biāo)與直角坐標(biāo)的關(guān)系，由復(fù)合函數(shù)的求導(dǎo)公式，有第19頁，課件共54頁，創(chuàng)作于2023年2月（7-26）（4）面積坐標(biāo)的積分公式下面給出面積坐標(biāo)的冪函數(shù)積分公式。它們在計(jì)算單元剛度矩陣和等效結(jié)點(diǎn)載荷時有用。第20頁，課件共54頁，創(chuàng)作于2023年2月在三角形單元上進(jìn)行積分時，有（7-27）在三角形某一邊（設(shè)ij邊，邊長為l）上進(jìn)行積分時，有（7-28）3、用面積坐標(biāo)表示六節(jié)點(diǎn)三角形單元計(jì)算公式

對應(yīng)如圖7-4所示的六節(jié)點(diǎn)三角形單元，形函數(shù)可用面積坐標(biāo)表示為

（1）形函數(shù)和位移表達(dá)式第21頁，課件共54頁，創(chuàng)作于2023年2月???ijmi

j

m

xy圖7-6現(xiàn)利用形函數(shù)的性質(zhì)檢驗(yàn)式（7-29）的正確性。先考慮三角形的角點(diǎn)，例如圖7-6中的i點(diǎn)，有由式（7-21）（P16），有代入式（7-29），有（7-29）第22頁，課件共54頁，創(chuàng)作于2023年2月再考慮三角形的邊中點(diǎn)，例如i

點(diǎn)，面積劃分如圖7-7所示。顯然有：???ijmi

j

m

xy圖7-7由式（7-21）(P16)，有代入式（7-29）(P16)，進(jìn)一步說明式（7-29）所表示的形函數(shù)的正確性。說明形函數(shù)Ni在i點(diǎn)等于1，在其它節(jié)點(diǎn)等于0，因此是正確的。第23頁，課件共54頁，創(chuàng)作于2023年2月形函數(shù)確定后，單元中任意一點(diǎn)的位移可以表示為：（7-30）其中（7-31）（7-32）其中I為二階單位陣，形函數(shù)由式（7-29）確定。（2）應(yīng)變矩陣單元中的應(yīng)變?nèi)钥杀硎緸椋海?-33）第24頁，課件共54頁，創(chuàng)作于2023年2月式中應(yīng)變矩陣[B]為：（7-34）其中（7-35）第25頁，課件共54頁，創(chuàng)作于2023年2月單元中的應(yīng)力仍可表示為：（3）應(yīng)力矩陣（7-36）式中[D]是彈性矩陣，由式（2-9）確定；應(yīng)變矩陣由式（7-34）、（7-35）確定。根據(jù)矩陣乘法，可以給出用面積坐標(biāo)表示的應(yīng)力矩陣[S]（4）單元剛度矩陣單元剛度矩陣仍可表示為：（7-37）根據(jù)[B]、[D]的表達(dá)式以及面積坐標(biāo)的積分公式（7-27），可以求出[k]中元素的顯式表示。由于較為繁復(fù)，這里就不列出詳細(xì)結(jié)果。

第26頁，課件共54頁，創(chuàng)作于2023年2月（5）等價節(jié)點(diǎn)力向量由于位移模式是非線性的，因此體積力和表面力引起的節(jié)點(diǎn)力向量不能采用靜力等效原理進(jìn)行分配，而應(yīng)采用相應(yīng)公式進(jìn)行計(jì)算。單元體積力引起的等價節(jié)點(diǎn)力計(jì)算公式仍為：（7-38）

將由式（7-29）、（7-32）表示的[N]代入，并應(yīng)用積分式（7-27），可以計(jì)算FV

e。例如對于重力引起的FV

e，有

它表示各邊中點(diǎn)承擔(dān)單元重力的1/3。第27頁，課件共54頁，創(chuàng)作于2023年2月單元表面力引起的結(jié)點(diǎn)力計(jì)算公式仍為：

（7-39）設(shè)在ij邊上受有x方向的均勻分布力ps，對應(yīng)的等價節(jié)點(diǎn)力向量為（圖7-8）pslh/6pslh/64pslh/6???ijmi

j

m

xy圖7-9lps???ijmi

j

m

xy圖7-8lpspslh/62pslh/6第28頁，課件共54頁，創(chuàng)作于2023年2月如在ij邊上受到x方向的三角形分布面力，其集度在i點(diǎn)為ps，在j點(diǎn)為0。對應(yīng)的節(jié)點(diǎn)力向量為（圖3-9）它表示邊中點(diǎn)承擔(dān)載荷的2/3，載荷集度大的角節(jié)點(diǎn)承擔(dān)1/3。六結(jié)點(diǎn)三角形單元中的應(yīng)變、應(yīng)力不為常量，因此可以應(yīng)用于應(yīng)力梯度較大的地方，精度較高。顯然，其計(jì)算也較復(fù)雜。

7.4四節(jié)點(diǎn)四邊形等參數(shù)單元1、等參數(shù)單元的概念

第29頁，課件共54頁，創(chuàng)作于2023年2月現(xiàn)在，我們從任意四邊形單元著手，介紹等參數(shù)單元的概念。1234xy圖7-10任意四邊形單元前面講到的四節(jié)點(diǎn)矩形單元雖然比較簡單，但難以應(yīng)用于斜線邊界。圖7-10所示四節(jié)點(diǎn)任意四邊形單元容易適應(yīng)這種邊界，但要在整體坐標(biāo)系內(nèi)，寫出它的統(tǒng)一的形函數(shù)又是相當(dāng)復(fù)雜和困難的。第30頁，課件共54頁，創(chuàng)作于2023年2月但是若能找到它與一個規(guī)則正方形的關(guān)系，就能寫出它的統(tǒng)一的位移模式，這可以通過坐標(biāo)變換來解決。在圖7-10所示四邊形單元上，用等分四邊的兩族直線分割該四邊形，以兩族曲線的中心（=0、=0）為原點(diǎn)，沿、增大的方向作軸和軸，并令四邊的=±1、=±1，就得出一組新坐標(biāo)系（圖7-11）。1234xy圖11實(shí)際單元=-1=1=1=-1

這里，、是一種局部（單元）坐標(biāo)，它只應(yīng)用于單元范圍內(nèi)。而x，y是整體（結(jié)構(gòu)）坐標(biāo)，它適用于所有的單元。圖中的任意四邊形單元是研究對象，稱為實(shí)際單元。第31頁，課件共54頁，創(chuàng)作于2023年2月參照式（7-2）和（7-3）P6，此基本單元位移函數(shù)可寫為：（7-40）

1234=-1=1=-1=1圖7-12基本單元為了得出實(shí)際單元的位移模式和局部坐標(biāo)與整體坐標(biāo)之間的變換關(guān)系，引入一個四節(jié)點(diǎn)的正方形單元，稱基本單元（圖7-12）。第32頁，課件共54頁，創(chuàng)作于2023年2月其中，形函數(shù)應(yīng)為：引入新變量

i、i（i=1,2,3,4)基本單元的形函數(shù)被寫成：第33頁，課件共54頁，創(chuàng)作于2023年2月（7-41）

現(xiàn)在，把基本單元的位移模式（7-40）和形函數(shù)式（7-41）移用于圖（7-11）所示的實(shí)際單元，則實(shí)際單元的位移模式取為：

（7-40）在結(jié)點(diǎn)處：在其它結(jié)點(diǎn)處：第34頁，課件共54頁，創(chuàng)作于2023年2月且，式中的形函數(shù)Ni仍由式（7-41）確定。而把式（7-41）中的、理解為圖7-11所示實(shí)際單元的局部坐標(biāo)，i、i便是實(shí)際單元中節(jié)點(diǎn)i的局部坐標(biāo)。（7-42）利用形函數(shù)的上述性質(zhì)，可以將任意四邊形的整體坐標(biāo)寫成：任意四邊形單元中結(jié)點(diǎn)的整體坐標(biāo)，如果它已知，那么(7-42)表示了局部坐標(biāo)與整體坐標(biāo)的變換第35頁，課件共54頁，創(chuàng)作于2023年2月另一方面，式（7-42）表明了實(shí)際單元中局部坐標(biāo)（、）與整體坐標(biāo)（x、y）的一一對應(yīng)關(guān)系，是一個坐標(biāo)變換式。實(shí)際單元是任意四邊形四節(jié)點(diǎn)單元，基本單元是正方形單元，可以認(rèn)為：實(shí)際單元是對基本單元通過變換得來的。由于實(shí)際單元的位移模式中采用了基本單元等同的形函數(shù)，這個實(shí)際單元就稱為等參數(shù)單元。類似于本章3.2節(jié)進(jìn)行的四結(jié)點(diǎn)矩形單元的特性分析，可以建立等參單元的應(yīng)變矩陣、應(yīng)力矩陣、剛度矩陣、節(jié)點(diǎn)力向量等的計(jì)算公式。與前面不同之處在于，在等參數(shù)單元法中，要將對整體坐標(biāo)x、y的導(dǎo)數(shù)計(jì)算和積分計(jì)算轉(zhuǎn)換為對局部坐標(biāo)、的導(dǎo)數(shù)計(jì)算和積分計(jì)算。第36頁，課件共54頁，創(chuàng)作于2023年2月例：實(shí)際單元的結(jié)點(diǎn)整體坐標(biāo)如圖(a)中括號內(nèi)數(shù)字所示，基本單元的結(jié)點(diǎn)局部坐標(biāo)如圖(b)中括號內(nèi)數(shù)字所示。圖(a)實(shí)際單元1(0,0)(1)試驗(yàn)證基本單元上的結(jié)點(diǎn)局部坐標(biāo)與實(shí)際單元上對應(yīng)點(diǎn)的整體坐標(biāo)的對應(yīng)關(guān)系。(2)求基本單元的局部坐標(biāo)原點(diǎn)(),在實(shí)際單元上的整體坐標(biāo)(x,y)是多少？4(0,1)3(1,2)2(2,0)xy0(3/4,3/4)01(－1,－1)2(1,－1)3(1,1)4(－1,1)圖(b)基本單元第37頁，課件共54頁，創(chuàng)作于2023年2月解（1）以3結(jié)點(diǎn)為例，根據(jù)形函數(shù)的性質(zhì)：在3結(jié)點(diǎn)處：應(yīng)用式（7-42）說明了由基本單元上結(jié)點(diǎn)的局部坐標(biāo)可映射出實(shí)際單元上對應(yīng)的結(jié)點(diǎn)整體坐標(biāo)（2）由于，由（7-40）P33得：第38頁，課件共54頁，創(chuàng)作于2023年2月代入（7-42）由上例可知：利用（7-42）在基本單元上任意一點(diǎn)，都可以在實(shí)際單元上找到一個對應(yīng)點(diǎn)的坐標(biāo)(x,y),這樣就把實(shí)際單元與基本單元緊密地聯(lián)系起來。反之，則比較困難，這是因?yàn)樾魏瘮?shù)是一個二次函數(shù)。為了避開這個困難，一般都假定基本單元上已知點(diǎn)去求實(shí)際單元上的對應(yīng)點(diǎn)。第39頁，課件共54頁，創(chuàng)作于2023年2月

2、應(yīng)變矩陣單元的幾何方程與式（7-8）、（7-9）相同，即：（7-8）（7-9）（7-43）式中（7-44）?第40頁，課件共54頁，創(chuàng)作于2023年2月這里采用記號由于形函數(shù)式（3-41）是用局部坐標(biāo)、給出的，將、看作x、y的函數(shù)，根據(jù)復(fù)合函數(shù)的求導(dǎo)規(guī)則，有：上式可記為：（7-45）第41頁，課件共54頁，創(chuàng)作于2023年2月上式右邊第一個矩陣稱為雅可比（Jacobi）矩陣：其逆矩陣為：式中|J|為雅可比行列式（7-48）由式（7-42）p35，有(7-46)(7-47)第42頁，課件共54頁，創(chuàng)作于2023年2月（7-49）由式（7-41）p34，有（7-50）由式（7-45）p41，有（7-51）第43頁，課件共54頁，創(chuàng)作于2023年2月式中分別由式（7-47）和（7-50）確定。從而由式（3-43）、（3-44）確定出應(yīng)變矩陣[B]。和3、應(yīng)力矩陣應(yīng)力矩陣仍由下式得到（7-52）

4、單元剛度矩陣單元剛度矩陣是一個8×8的矩陣，仍為第44頁，課件共54頁，創(chuàng)作于2023年2月（7-53）由于[B]是用局部坐標(biāo)系、給出的，坐標(biāo)變換時有面積微元公式因此，[k]可由下式計(jì)算5、等價節(jié)點(diǎn)力向量（1）體積力設(shè)單元的體積力是，則等價節(jié)點(diǎn)力公式為第45頁，課件共54頁，創(chuàng)作于2023年2月（7-54）(2)表面力設(shè)單元的某邊（如對應(yīng)的=±1）上作用有表面力ps=[psxpsy]T，則該邊上節(jié)點(diǎn)的節(jié)點(diǎn)力為（7-55）第46頁，課件共54頁，創(chuàng)作于2023年2月

上式中，用到坐標(biāo)變換時線積分微元公式(當(dāng)=常量時)p37對于=±1的邊界表面力，等價節(jié)點(diǎn)力的計(jì)算過程完全一樣。6、關(guān)于高斯（Gauss）積分在計(jì)算單元剛度矩陣式（7-53）和等價節(jié)點(diǎn)載荷向量式（7-54）、（7-55）時，由于被積函數(shù)比較復(fù)雜，通?？刹捎脭?shù)值積分。即在單元上選擇某些點(diǎn)，稱為第47頁，課件共54頁，創(chuàng)作于2023年2月積分點(diǎn)，求出被積函數(shù)在這些積分點(diǎn)上的數(shù)值，再用一些權(quán)函數(shù)乘這些函數(shù)值后求和，就可得到近似積分值。高斯積分法是數(shù)值積分法中具有較高精度的方法?，F(xiàn)簡要介紹如下：（1）一維高斯積分公式

（7-56）式中Hi對應(yīng)于積分點(diǎn)

i的權(quán)函數(shù)，對于n=2~4個積分點(diǎn)，其坐標(biāo)

i和權(quán)函數(shù)Hi列于表7-1中。第48頁，課件共54頁，創(chuàng)作于2023年2月n±

iHi20.57735026921.000000000030.77459666920.00000000000.55555555560.888888888940.86113631160.3399810436