三維Wilson元及在近不可壓彈性問題中的應用

李真有,肖映雄

(湘潭大學土木工程與力學學院， 湖南湘潭411105)

?

三維Wilson元及在近不可壓彈性問題中的應用

李真有,肖映雄

(湘潭大學土木工程與力學學院， 湖南湘潭411105)

為克服三維近不可壓縮問題的體積閉鎖現象，建立了兩種基于六面體單元的Wilson非協調元計算格式，并將其應用于兩類含混合邊界條件的近不可壓縮彈性問題的求解。數值結果表明:Wilson非協調元能有效克服三維體積閉鎖現象，與相同規模下的協調元相比較，它具有更高的計算精度。在三維有限元分析中，剖分網格的質量將對計算精度影響很大，實際計算時若能采用各向同性網格，則對問題的分析將具有更好的收斂性。

近不可壓縮問題；體積閉鎖；Wilson非協調元；網格質量

利用通常的低階協調有限元方法處理三維近不可壓縮問題常出現體積閉鎖現象。克服這種閉鎖現象的方法很多，如混合有限元法[1-3]，非協調有限元法[4-6]，高階協調元法[7-8]及減縮積分法[9-10]等。減縮積分使單元剛度矩陣秩降低，常引起多余的零能模式。與基于H-R變分原理的混合元格式相比，基于能量泛函極小的離散變分問題更容易被解決，并且有限元方程的系數矩陣是正定的，將給求解帶來很大方便。對三維問題，若利用高階協調元方法，要求單元階次p≥8，這在實際計算中是不可取的。此時，常采用低階非協調元(如二次元、線性元)來克服三維近不可壓縮問題的閉鎖現象[5,11]，這種方法具有自由度少、精度高的特點，可解決三維問題計算規模過大的困難。

在三維有限元分析中，常采用六面體網格剖分，這種網格在計算精度、劃分數量、抗畸變程度等方面比四面體網格具有明顯的優勢。基于六面體單元的Wilson非協調元，它通過在單元內部設置附加自由度，以提高完全多項式的次數、改善計算精度。關于Wilson元的研究和應用，已有很多研究成果[12-15]。本文建立了兩種基于六面體單元的Wilson非協調元，并將其應用于懸臂梁和Cook膜兩類近不可壓縮彈性問題的求解。數值結果表明：Wilson非協調元能有效克服三維體積閉鎖現象，與相同規模下的協調元相比較，它具有更高的計算精度。還研究了有限元網格質量對Wilson元計算精度的影響，計算中若能采用各向同性網格，則對問題的分析將具有更好的收斂性。

1　近不可壓縮問題和體積閉鎖現象

考慮下述三維線彈性力學模型:

(1)

其中，Ω?R3，?Ω=Γ0∪Γ1且Γ0∩Γ1=Φ, u為位移向量，f為外力，g為邊界Γ1上的面力，n為邊界Γ1的單位外法線向量，拉梅常數λ和μ可分別用楊氏模量E和泊松比ν表示為:

(2)

這里，泊松比ν∈[0,0.5)，當ν→0.5時，稱相應的彈性力學問題為近不可壓縮問題。

記Hm(Ω)={v|?αv∈L2(Ω),|α|≤m}，則問題(1)對應的變分問題可描述為：

(3)

在有限元分析中，網格的劃分及單元的選取是影響其計算精度的重要因素。四面體網格和六面體網格是三維有限元分析中常采用的兩種網格類型，而六面體網格在計算精度、網格數量、抗畸變程度等方面比四面體網格具有明顯的優勢，已成為三維有限元分析中的首選網格。設Th是區域Ω上的六面體網格剖分，其中h為Th上所有剖分單元的最大直徑。引入如下p次拉格朗日有限元空間：

(4)

設變分問題(4)的p次有限元離散化線性系統的矩陣形式為:

(5)

對近不可壓縮問題，通常的低階協調元(如三線性元、三二次元)解不再收斂到原問題的解或達不到最優收斂階，即出現所謂的體積閉鎖現象。下面，舉兩個例子來進行驗證。

圖1　懸臂梁幾何結構示意圖

nx×ny×nzv=0.3v=0.4v=0.49v=0.499uAzuBzuAzuBzuAzuBzuAzuBz16×16×16-0.2186-0.6157-0.2113-0.6023-0.1584-0.4900-0.0739-0.272424×24×24-0.2370-0.6663-0.2308-0.6545-0.1893-0.5680-0.1070-0.371240×4×4-0.2455-0.6895-0.2376-0.6710-0.1906-0.5545-0.1257-0.396580×8×8-0.2527-0.7087-0.2485-0.6994-0.2244-0.6418-0.1681-0.5008160×16×16-0.2547-0.7140-0.2518-0.7079-0.2410-0.6832-0.2061-0.5964

表2　懸臂梁軸線上A點和B點20節點二次元解Tab.2　The triquadratic element solutions of two typical points A and B on the cantilever axis

算例2(Cook膜問題)考慮如圖2所示的三維Cook膜問題，其左側面(即x=0)固定，右側面上面力g=(0,0,t)T，其中，t=10N/mm2，彈性模量E=103N/mm2，不計自重。采用六面體網格進行剖分，分別采用8節點線三性元和20節點三二次元進行計算，相應的數值結果如表3和表4所示。

圖2　三維Cook膜問題Fig.2　Three dimensional Cook’s membrane problems表3　結構網格下特征點C在z方向上的8節點線性元解Tab.3　The trilinear element solutions of on the structured meshes

表4　結構網格下特征點C在z方向上的20節點二次元解Tab.4　The triquadratic element solutions of on the structured meshes

由上述數值結果可知：對于可壓問題，隨著網格規模的不斷增加，特征點處的8節點線性元解均能收斂于理論解，而對近不可壓縮問題，當ν→0.5時，8節點線性元解已不再收斂，20節點二次元解盡管收斂，但達不到最佳收斂結果，出現所謂的體積閉鎖現象。另外，計算中所使用網格單元的質量也會影響計算結果的精度，如采用各向同性網格，則具有更好的計算精度。

克服這種閉鎖現象的方法很多，如混合有限元法，非協調有限元法，高階協調元法及減縮積分法等。對三維問題，常采用低階非協調元(如二次元、線性元)來克服閉鎖現象。這種方法具有自由度少、精度高的優點，可有效解決三維實際問題計算量過大的困難。

2　Wilson非協調元

六面體單元對應的位移模式不是完全多項式，而決定有限元解精度的是完全多項式的次數，非完全的高次項對改善精度不起作用，有時還可能會起相反的作用。Wilson非協調元通過在單元內部設置附加自由度，達到提高完全多項式的次數、改善計算精度的目的。這種單元的位移模式中包含了剛體位移和常量應變，但在相鄰單元邊界上不再保持連續性。筆者對如圖3所示的局部坐標系o-ξηζ下的8節點和20節點六面體單元，分別建立Wilson非協調元計算格式。

首先，設8節點六面體單元的非協調元位移模式為：

(6)

類似地，可得到20節點六面體單元形如式(6)所示的位移模式，其中：

利用最小位能原理，可建立如下形式的單元剛度方程：

(7)

消去內部附加自由度αe，可得凝聚后的單元剛度方程為：

Keae=Pe,

(8)

將上述8節點和20節點Wilson非協調元應用于兩類近不可壓縮問題：懸臂梁和Cook膜，以驗證其有效性。

首先，對懸臂梁問題，當ν→0.5時，在不同網格規模和不同網格特性(各向同性網格和各向異性網格)下列出了懸臂梁軸線上A點和B點的兩種Wilson非協調元解，總結如表5和表6所示。

表5懸臂梁軸線上A點和B點8節點Witson元解

Tab.5The 8-node Wilson element solutions of two typical pointsAandB

ν0.40.450.490.499nx×ny×nzuAzuBzuAzuBzuAzuBzuAzuBz16×16×16-0.2440-0.6935-0.2387-0.6828-0.2333-0.6718-0.2308-0.666724×24×24-0.2478-0.7008-0.2441-0.6932-0.2402-0.6853-0.2389-0.682640×4×4-0.2505-0.7061-0.2479-0.7008-0.2448-0.6946-0.2409-0.686580×8×8-0.2522-0.7093-0.2503-0.7054-0.2483-0.7014-0.2474-0.6994160×16×16-0.2529-0.7107-0.2512-0.7073-0.2496-0.7039-0.2491-0.7030

表6　懸臂梁軸線上A點和B點20節點Witson元解Tab.6　The 20-node Wilson element solutions of two typical points A and B

表7結構化網格下特征C在z方向上的8節點WILSON解

ν0.400.490.4950.4990.499516×16×163.954773.936903.932953.911783.9039124×24×243.977763.699743.967953.962543.957828×2×83.886813.829643.811363.725023.6493416×2×163.957153.941213.937593.922513.9083632×4×323.989303.985683.984773.982863.98116

表8　結構化網格下特征點C在z方向上的20節點Wilson解Tab.8　The 20-node Wilson element solution of on the structured meshes

(a) 336個單元

(b) 11 856個單元

圖4三維Cook膜問題兩種非結構化網格剖分

Fig.4Two unstructured meshes used for Cook’s membrane

表9非結構網格下特征C在z方向上的8節點Wilson解

表10　非結構網格下特征C在z方向上的20節點Wilson解Tab.10　The 20-node Wilson element solution of on the structured meshes

3　結　語

對三維近不可壓問題，若利用通常的低階協調元方法求解會出現所謂的體積閉鎖現象。Wilson非協調元法是克服三維體積閉鎖現象的一種有效方法，它具有自由度少、精度高的優點，可解決三維問題計算規模大的難題。Wilson非協調元對網格質量依賴性強，如何設計適用于任意六面體網格單元的低階非協調元方法，將是我們進一步研究的問題之一。另外，針對這種非協調有限元分析中形成的大型的、稀疏的和高度病態的正定方程組，如何設計高效求解算法是一個非常困難的問題，將決定其有限元分析的整體效率，這也是我們今后將進一步研究的問題。

[1]CHEUNG Y K, CHEN W J.Hybrid element method for incompressible and nearly incompressible materials[J]. International Journal of Solids and Structures, 1989, 25(5): 483-495.

[2]MORLEY M.A mixed family of elements for linear elasticity[J]. Numerische Mathematik, 1989(55):633-666.

[3]CHAMA A, REDDY B D.New stable mixed finite element approximations for problems in linear elasticity[J]. Computer Methods in Applied Mechanics and Engineering, 2013, 256: 211-223.

[4]WANG L H, QI H.A locking-free scheme of nonconforming rectangular finite element for the planar elasticity[J]. Journal of Computational Mathematics, 2004, 22: 641-650.

[5]FALK R S.Nonconforming finite element methods for the equations of linear elasticity[J]. Mathematics of Computation, 1991, 57: 529- 556.

[6]BRENNER S C.A nonconforming mixed multigrid method for the pure traction problem in planar linear elasticity[J]. Mathematics of Computation, 1994, 63: 435-460.

[7]SCOTT L R, VOGELIUS M.Conforming finite element methods for incompressible and nearly incompressible continua[C]//Large Scale Computations in Fluid Mechanics. Providence RI. American Mathematical Socity, USA: 1985, 22: 221-244.

[8]STENBERG R, SURI M.Mixed h-p finite element methods for problems in elasticity and Stokes flow[J]. Numerische Mathematik, 1996, 72: 367-390.

[9]HUGHES T J R.Equivalence of finite elements for nearly incompressible elasticity[J]. Journal of Applied Mechanics, 1977, 44: 181-183.

[10]MALKUS D S, HUGHES T J R.Mixed finite element methods-reduced and selective integration techniques: a unification of concepts[J]. Computer Methods in Applied Mechanics and Engineering, 1978, 15 (1): 63-81.

[11]QI H, WANG L H, ZHENG W Y.On locking-free finite element schemes for three dimensional elasticity[J]. Journal of Computational Mathematics, 2005, 23:101-112.

[12]鹿曉陽.Wilson 非協調元構造機理研究[J]. 工程力學，2000, 17(5): 58-62.

[13]張春生，龍馭球，須寅.三維內參型附加非協調位移基本項[J]. 工程力學，2000, 18(5): 58-62.

[14]陳震，林府標.三維Wilson元與特征值下界[J]. 南昌大學學報(理科版)，2010, 34(1): 19-23.

[15]王芬玲，石東偉，石東洋.擬線性Sobolev方程Wilson元解的超收斂分析及外推[J]. 工程數學學報，2012, 29(5): 720- 724.

(責任編輯唐漢民梁碧芬)

Wilson element and its application in nearly incompressible elasticity problems in three dimensions

LI Zhen-you, XIAO Ying-xiong

(Civil Engineering and Mechanics College, Xiangtan University, Xiangtan 411105, China)

In order to overcome the volume locking phenomenon of nearly incompressible problems in three dimensions, two types of Wilson nonconforming finite elements have been presented based on the hexahedral elements, and the resulting methods are then applied to the solution of two nearly incompressible problems with mixed boundary conditions.The numerical results have been shown that Wilson elements can effectively overcome the locking phenomenon, and it has higher accuracy compared with the conforming elements under the same mesh size.In three-dimensional finite element analysis, the quality of the mesh will have a great effect on the accuracy, and if the isotropic grids can be used in the practical calculations, the method will have better convergence.

nearly incompressible problem；volume locking；Wilson nonconforming element；mesh quality

2016-04-22；

2016-05-01

國家自然科學基金資助項目(10972191);湖南省自然科學基金資助項目(14JJ2063);湖南省教育廳資助科研項目(15A183)

肖映雄(1970—)，男，湖南城步人，湘潭大學教授；E-mail: xyx610xyx@xtu.edu.cn。

10.13624/j.cnki.issn.1001-7445.2016.1271

O343.2

A

1001-7445(2016)04-1271-08

引文格式：李真有,肖映雄.三維Wilson元及在近不可壓彈性問題中的應用[J].廣西大學學報(自然科學版)，2016，41(4)：1271-1278.