基于赫林格-賴斯納變分原理的一致高效無(wú)網(wǎng)格本質(zhì)邊界條件施加方法1)

吳俊超 吳新瑜 趙珧冰 王東東

*(華僑大學(xué)土木工程學(xué)院,福建省智慧基礎(chǔ)設(shè)施與監(jiān)測(cè)重點(diǎn)實(shí)驗(yàn)室,福建廈門 361021)

?(廈門大學(xué)土木工程系,福建省濱海土木工程數(shù)字仿真重點(diǎn)實(shí)驗(yàn)室,福建廈門 361005)

引言

無(wú)網(wǎng)格法[1-7]是指基于節(jié)點(diǎn)信息構(gòu)造離散形函數(shù)的一類數(shù)值方法的總稱.該類方法通常具有高階光滑、全域協(xié)調(diào)的形函數(shù),形函數(shù)構(gòu)造過(guò)程不依賴于網(wǎng)格單元的拓?fù)湫畔?適用于大變形分析[8]、薄板殼高階問(wèn)題[9]及裂紋擴(kuò)展模擬[10]等.然而,與傳統(tǒng)有限元法相比,高階連續(xù)光滑的特點(diǎn)導(dǎo)致無(wú)網(wǎng)格形函數(shù)在離散節(jié)點(diǎn)上通常不具有插值性,在求解過(guò)程中難以直接施加本質(zhì)邊界條件[11-12].其次,無(wú)網(wǎng)格形函數(shù)通常為有理式,并且形函數(shù)影響域高度重疊,導(dǎo)致形函數(shù)在背景積分單元上為分段的有理式.在伽遼金法的求解過(guò)程中,傳統(tǒng)高斯積分法無(wú)法精確數(shù)值積分由形函數(shù)組成的剛度矩陣和力向量,導(dǎo)致伽遼金弱形式不滿足積分約束條件或稱變分一致性條件[13-14],無(wú)法保證計(jì)算精度和最優(yōu)誤差收斂率.

為了使無(wú)網(wǎng)格法能夠直接施加本質(zhì)邊界條件,許多學(xué)者構(gòu)造了諸多具有插值性的無(wú)網(wǎng)格近似方法,例如奇異權(quán)函數(shù)法[15]、插值最小二乘法[16-17]、復(fù)變量移動(dòng)最小二乘法[18]、廣義移動(dòng)最小二乘法[19]、變換法[20]等.然而,這類方法不是建立在變分原理基礎(chǔ)上,并不能保證節(jié)點(diǎn)之間位移邊界條件施加精度和無(wú)網(wǎng)格法的變分一致性.此外,對(duì)于滿足積分約束條件的無(wú)網(wǎng)格數(shù)值積分方法,例如穩(wěn)定節(jié)點(diǎn)積分法[13]、一致性積分法[21-22]、變分一致積分法[23]、嵌套子域積分法[24]、再生光滑梯度積分法[25]等,在計(jì)算過(guò)程中采用形函數(shù)的光滑梯度替代傳統(tǒng)無(wú)網(wǎng)格形函數(shù)的導(dǎo)數(shù),在保證無(wú)網(wǎng)格法的計(jì)算精度和最優(yōu)誤差收斂率的同時(shí)提高了計(jì)算效率.但其本質(zhì)邊界條件仍需要具有變分一致性的方法進(jìn)行施加[26-27].

在伽遼金無(wú)網(wǎng)格法中,拉格朗日乘子法[1]是常采用的一種變分一致本質(zhì)邊界條件施加方法.該方法需要在整體剛度矩陣上增加額外的自由度離散拉格朗日乘子.當(dāng)采用滿足變分一致無(wú)網(wǎng)格數(shù)值積分方法時(shí),拉格朗日乘子的自由度需要與光滑梯度構(gòu)造過(guò)程中的數(shù)值積分點(diǎn)相一致.但是,過(guò)多的拉格朗日乘子自由度將增加整體剛度矩陣奇異性,不適用于高階基函數(shù)的無(wú)網(wǎng)格法.Lu等[28]根據(jù)拉格朗日乘子的物理意義,采用位移自由度離散拉格朗日乘子,無(wú)需增加額外自由度,并稱之為修正變分原理法.但該方法施加本質(zhì)邊界條件過(guò)程中的修正變分項(xiàng)降低了剛度矩陣的正定性,計(jì)算精度不高.Zhu和Atluri[29]將罰函數(shù)法引入無(wú)網(wǎng)格法中施加本質(zhì)邊界條件,該方法格式簡(jiǎn)潔易實(shí)現(xiàn),但其形式并不滿足變分的一致性,求解精度依賴人工參數(shù),穩(wěn)定性較低.Fernández-Méndez和Huerta[11]采用尼茲法(Nitsche method)施加無(wú)網(wǎng)格法本質(zhì)邊界條件,該方法可視為修正變分原理與罰函數(shù)法的結(jié)合,修正變分項(xiàng)保證了變分一致性,而罰函數(shù)項(xiàng)恢復(fù)了剛度矩陣的正定性,提高了求解的穩(wěn)定性,但該方法仍然需要選擇合適的人工參數(shù)來(lái)保證數(shù)值解的精度.

本文針對(duì)滿足積分約束條件的伽遼金無(wú)網(wǎng)格法,提出了一種具有變分一致且高效的本質(zhì)邊界條件施加方法.該方法以赫林格-賴斯納(Hellinger-Reissner,HR)變分原理[30]為基礎(chǔ),采用混合離散的方式近似弱形式中的位移和應(yīng)力.其中,位移采用傳統(tǒng)無(wú)網(wǎng)格形函數(shù)進(jìn)行離散,而應(yīng)力則在每個(gè)背景積分單元中近似為局部多項(xiàng)式.此時(shí),近似應(yīng)力中多項(xiàng)式的系數(shù)可在背景積分單元中采用位移自由度表示,無(wú)需額外增加整體自由度.不僅如此,應(yīng)力的系數(shù)表達(dá)式中包含再生光滑梯度表達(dá)式,在再生光滑梯度的理論框架下,光滑梯度構(gòu)造過(guò)程中可采用優(yōu)化的數(shù)值積分點(diǎn)以減少全域無(wú)網(wǎng)格形函數(shù)的計(jì)算量,計(jì)算效率高.整體的離散平衡方程具有與尼茲法相類似的形式,其中尼茲法中的修正變分項(xiàng)分別引入了再生光滑梯度和無(wú)網(wǎng)格形函數(shù)離散應(yīng)力和位移,整體求解過(guò)程無(wú)需計(jì)算傳統(tǒng)無(wú)網(wǎng)格形函數(shù)的導(dǎo)數(shù).而在赫林格-賴斯納變分原理的框架下,弱形式中自帶穩(wěn)定項(xiàng),無(wú)需增加額外的穩(wěn)定項(xiàng),且穩(wěn)定項(xiàng)中不存在人工參數(shù),消除了對(duì)人工參數(shù)的依賴性.文中通過(guò)二維彈性力學(xué)算例驗(yàn)證了所提方法的精度、效率和誤差收斂性.

1 赫林格-賴斯納變分原理

不失一般性,在求解域 Ω 內(nèi)考慮如下彈性力學(xué)控制方程

式中,W為彈性體的余能密度函數(shù),其與應(yīng)力之間的關(guān)系為

其中Cijkl為四階彈性張量.

對(duì)上式進(jìn)行變分可得與之相對(duì)應(yīng)的弱形式

2 位移離散與再生核近似

在赫林格-賴斯納變分原理的弱形式中,位移與應(yīng)力可分別采用不同的近似方案進(jìn)行離散.這里,位移采用基于再生核近似[29]的無(wú)網(wǎng)格形函數(shù)進(jìn)行離散.首先,在求解域 Ω 上布置一系列無(wú)網(wǎng)格節(jié)點(diǎn).此時(shí),位移分量ui的近似表達(dá)式可表示為

式中,diI為無(wú)網(wǎng)格節(jié)點(diǎn)xI上的位移分量節(jié)點(diǎn)系數(shù),ΨI為與之相對(duì)應(yīng)的無(wú)網(wǎng)格形函數(shù),根據(jù)再生核近似理論[2],無(wú)網(wǎng)格形函數(shù)具有如下表達(dá)式

其中,p為基函數(shù)向量,以二維p階基函數(shù)為例,此時(shí)基函數(shù)向量為

φ為核函數(shù),其影響域也決定了無(wú)網(wǎng)格形函數(shù)的影響域.在二維情況下,核函數(shù)的影響域通常為圓形或矩形,本文采用基于三次樣條函數(shù)、具有矩形影響域的核函數(shù)

式中,rx=|xI-x|/sx,ry=|yI-y|/sy,sx和sy分別為x和y方向的影響域尺寸.φ為三次樣條函數(shù)[1]

無(wú)網(wǎng)格形函數(shù)表達(dá)式(8)中的c為修正函數(shù)向量,該表達(dá)式可通過(guò)滿足再生條件[2]確定

將式無(wú)網(wǎng)格形函數(shù)表達(dá)式(8) 代入再生條件(13)中,可得

式中M稱為矩量矩陣

將式(14)、式(15)代入式(8)中,可得到再生核無(wú)網(wǎng)格形函數(shù)的表達(dá)式

圖1為二維二次基函數(shù)無(wú)網(wǎng)格形函數(shù),從圖中可以看出,無(wú)網(wǎng)格形函數(shù)在全域上連續(xù)光滑.但在無(wú)網(wǎng)格節(jié)點(diǎn)處,無(wú)論節(jié)點(diǎn)位于域內(nèi)或邊界處,形函數(shù)都不具有插值性,這是伽遼金無(wú)網(wǎng)格法難以施加本質(zhì)邊界的主要原因.同時(shí),相較于有限元形函數(shù),無(wú)網(wǎng)格形函數(shù)并不具有顯示表達(dá)式,形函數(shù)計(jì)算復(fù)雜耗時(shí),無(wú)網(wǎng)格形函數(shù)的導(dǎo)數(shù)更是如此.

圖1 無(wú)網(wǎng)格形函數(shù)Fig.1 Meshfree shape functions

3 應(yīng)力離散與再生光滑梯度近似

表述方便起見(jiàn),引入沃伊特標(biāo)記(Voigt notation),將近似應(yīng)力 σh采用矩陣向量形式改寫為

4 赫林格-賴斯納變分原理下的本質(zhì)邊界條件施加方法

基于赫林格-賴斯納變分原理的弱形式(4)中已考慮了本質(zhì)邊界條件.前兩節(jié)已分別采用再生核無(wú)網(wǎng)格形函數(shù)和背景積分單元上的分段多項(xiàng)式近似式(4)中的位移和應(yīng)力,進(jìn)一步將近似的位移分量表達(dá)式(7)和應(yīng)力分量表達(dá)式(17)代入弱形式(6)中,有

其中,δd為位移變分節(jié)點(diǎn)系數(shù)向量,f為外力向量.將式(19)代入式(33)中,式(33)等式左邊可化簡(jiǎn)為

式中最后一個(gè)等式引入了關(guān)系式(20)～式(26)進(jìn)行化簡(jiǎn),詳細(xì)的推導(dǎo)過(guò)程可參考文末附錄A.將簡(jiǎn)化后的式(34)代入式(33),即可得到離散的平衡方程

從式(35)和式(36)可以看出,在赫林格-賴斯納變分原理的框架下,本質(zhì)邊界條件的施加過(guò)程與尼茲法[11]具有相類似的格式.清晰起見(jiàn),這里給出基于尼茲法的無(wú)網(wǎng)格法離散平衡方程

其中,K和f與式(36) 中的剛度矩陣K和力向量f相同.Kv和fv是基于變分原理施加本質(zhì)邊界條件產(chǎn)生的剛度和力向量,兩者保證了邊界條件施加方法的變分一致性,其具體表達(dá)式為

再者,式(40)中的Kp和fp對(duì)應(yīng)于尼茲法中采用罰函數(shù)法施加穩(wěn)定項(xiàng)的剛度矩陣和力向量

式中 1為單位對(duì)角矩陣

該部分通過(guò)設(shè)置合適的罰因子 α[31]來(lái)保證剛度矩陣的正定性,提高求解的穩(wěn)定性.

對(duì)照本文所提方法和尼茲法的無(wú)網(wǎng)格離散方程,即式(35)和式(40),可見(jiàn)式(35)中的對(duì)應(yīng)于尼茲法中利用變分原理施加本質(zhì)邊界條件,但其應(yīng)變離散過(guò)程采用式(25)中的光滑梯度替代傳統(tǒng)無(wú)網(wǎng)格形函數(shù)梯度,而剛度矩陣的正定性則是通過(guò)施加得以保證.與尼茲法相比,本文所提方法的穩(wěn)定項(xiàng)表達(dá)式相對(duì)較為復(fù)雜,但在求解全過(guò)程中,無(wú)需計(jì)算復(fù)雜耗時(shí)的無(wú)網(wǎng)格形函數(shù)梯度,整體過(guò)程計(jì)算效率優(yōu)于尼茲法,且穩(wěn)定項(xiàng)也不包括任何人工參數(shù).

在引入數(shù)值積分求解式(35)的過(guò)程中,為保證離散弱形式的變分一致性,數(shù)值積分在求解的各個(gè)過(guò)程也需保持一致.如圖2 所示,所提方法需要兩套數(shù)值積分點(diǎn),和中的數(shù)值積分方案需要與,,f,和保持一致,為了保證理論誤差收斂率,該數(shù)值積分方案的積分精度在文獻(xiàn)[25]有詳細(xì)討論.另一方面,變分一致性要求式(20)中G與式(35)中的K采用相同數(shù)值積分方案.在所提方法中,由于再生光滑梯度已滿足積分約束條件,再生光滑梯度也為低階多項(xiàng)式,此時(shí)數(shù)值積分方案采用傳統(tǒng)有限元中對(duì)應(yīng)階次的高斯積分即可保證數(shù)值精度.

圖2 背景積分單元示意圖和優(yōu)化的數(shù)值積分方案Fig.2 Illustration of background integration cells and optimized quadrature rules

此外,構(gòu)造再生光滑梯度需要計(jì)算積分點(diǎn)處無(wú)網(wǎng)格形函數(shù).為了減小形函數(shù)的計(jì)算量、提高計(jì)算效率,這里采用無(wú)網(wǎng)格再生光滑梯度積分法[25],利用積分點(diǎn)在背景積分單元間的共享特性,從全域上優(yōu)化整體求解過(guò)程中積分點(diǎn)的數(shù)量,提高計(jì)算效率.本文采用二維二次、三次基函數(shù)的優(yōu)化數(shù)值積分方案,具體的數(shù)值積分點(diǎn)的位置與權(quán)重詳見(jiàn)附錄B.

5 數(shù)值算例

本節(jié)首先通過(guò)分片試驗(yàn),驗(yàn)證了采用傳統(tǒng)高斯積分法和再生光滑梯度積分法的不同本質(zhì)邊界條件施加方法下是否滿足積分約束條件.然后,通過(guò)典型的懸臂梁?jiǎn)栴}和帶孔無(wú)限大平板問(wèn)題來(lái)驗(yàn)證采用二次和三次基函數(shù)時(shí),所提算法的計(jì)算精度和效率.在精度分析中,采用如下的位移誤差(L2-Err)和能量誤差(He-Err)

表述清晰起見(jiàn),下面的討論中用“GI”和“RKGSI”表示高斯積分法和再生光滑梯度積分法,用“penalty”、“LM”和“Nitsche”表示罰函數(shù)法、拉格朗日乘子法和尼茲法這三種施加本質(zhì)邊界條件的方法.例如,本文所提的基于赫林格-賴斯納變分原理的本質(zhì)邊界條件施加方法簡(jiǎn)寫為“RKGSI-HR”.為了保證拉格朗日乘子法求解的穩(wěn)定性,拉格朗日乘子均采用線性有限元形函數(shù)進(jìn)行離散.當(dāng)采用二次基函數(shù)時(shí),GI 的域內(nèi)積分采用13 點(diǎn)高斯積分,邊界積分采用3 點(diǎn)高斯積分;當(dāng)采用三次基函數(shù)時(shí),GI 采用16 點(diǎn)高斯積分,邊界積分采用5 點(diǎn)高斯積分.

5.1 分片試驗(yàn)

首先,采用線性、二次和三次彈性力學(xué)分片試驗(yàn)驗(yàn)證不同無(wú)網(wǎng)格本質(zhì)邊界條件施加方案的變分一致性.分片試驗(yàn)的求解域?yàn)檫呴L(zhǎng)等于1 的正方形,求解域的四邊施加本質(zhì)邊界條件,邊界條件的預(yù)設(shè)值根據(jù)如下所示的精確解施加

其中,n=1,2,3 代表線性、二次和三次分片試驗(yàn).分片試驗(yàn)采用如圖3 所示的非均布節(jié)點(diǎn)離散求解域,針對(duì)二次基函數(shù)的無(wú)網(wǎng)格近似,采用線性和二次分片試驗(yàn)進(jìn)行測(cè)試,核函數(shù)相對(duì)影響域?yàn)?.5;而三次基函數(shù)的無(wú)網(wǎng)格近似采用二次和三次分片試驗(yàn),核函數(shù)的相對(duì)影響域?yàn)?.5.

圖3 分片試驗(yàn)無(wú)網(wǎng)格離散模型Fig.3 Meshfree discretization for the patch test

表1、表2 分別為具有二次、三次基函數(shù)的無(wú)網(wǎng)格法分片試驗(yàn)結(jié)果,從表中可以看出,傳統(tǒng)高斯積分法由于不滿足積分約束條件,即使采用高階高斯積分也不能通過(guò)分片試驗(yàn).采用滿足積分約束條件的再生光滑梯度積分法時(shí),罰函數(shù)法不具有變分一致性,無(wú)法通過(guò)分片試驗(yàn).而拉格朗日乘子法由于其拉格朗日乘子采用線性形函數(shù)進(jìn)行離散,無(wú)法與再生光滑梯度相匹配,也不能通過(guò)分片試驗(yàn).當(dāng)采用尼茲法和本文所提基于赫林格-賴斯納變分原理的本質(zhì)邊界施加方法時(shí),RKGSI 可通過(guò)分片試驗(yàn).

表1 二次基函數(shù)無(wú)網(wǎng)格法分片試驗(yàn)結(jié)果Table 1 The results of patch test with quadratic basis functions

表2 三次基函數(shù)無(wú)網(wǎng)格法分片試驗(yàn)結(jié)果Table 2 The results of patch test with cubic basis functions

5.2 懸臂梁?jiǎn)栴}

考慮圖4 所示懸臂梁?jiǎn)栴},懸臂梁求解域的長(zhǎng)和寬分別為L(zhǎng)=48,D=12,懸臂梁楊氏模量E=3×106、泊松比 ν=0.3.懸臂梁的右端沿y軸正方向施加荷載P,而左端為固定支座.根據(jù)平面應(yīng)力假設(shè)和圣維南原理,該問(wèn)題的解析解[32]為

圖4 懸臂梁?jiǎn)栴}模型Fig.4 Description of the cantilever beam problem

與之相應(yīng)的應(yīng)力分量為

為保證一致性,該問(wèn)題的本質(zhì)邊界條件和自然邊界條件分別根據(jù)解析解式(47)和式(48)進(jìn)行施加.

懸臂梁求解域采用圖5 所示的四個(gè)疏密程度不同的節(jié)點(diǎn)離散研究計(jì)算誤差收斂特性,算例采用具有二次基函數(shù)的無(wú)網(wǎng)格近似,相應(yīng)的核函數(shù)相對(duì)影響域?yàn)?.5.圖6為該問(wèn)題的位移誤差和能量誤差收斂率結(jié)果,結(jié)果表明所提RKGSI-HR 方法的精度與RKGSI-Nitsche 相當(dāng),均達(dá)到理論收斂率,但RKGSIHR 方法無(wú)需人工參數(shù).傳統(tǒng)高斯積分法、RKGSILM和RKGSI-penalty 均不具有變分一致性,不能達(dá)到理論誤差收斂率.

圖5 懸臂梁?jiǎn)栴}節(jié)點(diǎn)離散Fig.5 Meshfree discretizations of the cantilever beam problem

圖6 懸臂梁?jiǎn)栴}誤差對(duì)比Fig.6 Error comparison for the cantilever beam problem

圖7為該問(wèn)題的效率對(duì)比,其中圖7(a)為整體求解過(guò)程中的效率分析.從圖7(a) 中可以看出,RKGSI 整體的計(jì)算效率要優(yōu)于GI,GI 所用CPU 耗時(shí)大概為RKGSI 的2.5 倍.這是由于RKGSI 整體無(wú)需計(jì)算傳統(tǒng)無(wú)網(wǎng)格形函數(shù)導(dǎo)數(shù),同時(shí)采用優(yōu)化的數(shù)值積分點(diǎn)減小形函數(shù)的計(jì)算量,兩方面共同作用的結(jié)果.此外,從該圖可以看出拉格朗日乘子法施加邊界條件時(shí),整體的計(jì)算效率要略高于其他方法,這是由于拉格朗日乘子法需要增加自由度離散拉格朗日乘子,導(dǎo)致整體線性代數(shù)的計(jì)算量增加.為了更好說(shuō)明本文所提RKGSI-HR 的計(jì)算效率,圖7(b)對(duì)比了不同本質(zhì)邊界方法施加過(guò)程的效率對(duì)比.該圖對(duì)比了施加過(guò)程計(jì)算形函數(shù)及梯度和組裝相應(yīng)的剛度矩陣和力向量所用時(shí)間,從圖中可以看出,罰函數(shù)法和拉格朗日乘子法由于計(jì)算過(guò)程僅需要無(wú)網(wǎng)格形函數(shù)本身,兩種方法的計(jì)算形函數(shù)耗時(shí)相同且優(yōu)于其余兩種方法.相比于傳統(tǒng)尼茲法,基于赫林格-賴斯納變分原理施加本質(zhì)邊界條件過(guò)程也無(wú)需計(jì)算形函數(shù)梯度,但需要構(gòu)造相應(yīng)的光滑梯度,此部分所用施加大約為罰函數(shù)法和拉格朗日乘子法的1.6 倍.而在組裝力向量方面,赫林格-賴斯納變分原理法和尼茲法具有相當(dāng)?shù)男?拉格朗日乘子法由于采用有限元形函數(shù)離散拉格朗日乘子,其計(jì)算效率最高.需要指出的是,罰函數(shù)法和拉格朗日乘子法并不能保證理論誤差收斂率,且從整體看拉格朗日乘子法的計(jì)算效率并不優(yōu)于其他方法.綜合本算例的精度分析和效率分析可得,RKGSI-HR 不僅保證計(jì)算精度和最優(yōu)理論收斂性,并且相較于傳統(tǒng)尼茲法具有更高的計(jì)算效率.

圖7 懸臂梁?jiǎn)栴}效率對(duì)比Fig.7 Efficiency comparison for the cantilever beam problem

5.3 帶孔無(wú)限大平板問(wèn)題

考慮圖8 所示帶孔無(wú)限大平板問(wèn)題,其中板的中心有一半徑為a=1 的小孔,平板的無(wú)窮遠(yuǎn)處沿x軸方向施加均布荷載T=1000.平板材料參數(shù)為:楊氏模量E=3×106、泊松比 ν=0.3.根據(jù)Michell解[32]可得該問(wèn)題的解析解為

其中 κ和μ為

與之相對(duì)應(yīng)的應(yīng)力表達(dá)式為

根據(jù)帶孔無(wú)限大平板的對(duì)稱性,如圖8 所示,取邊長(zhǎng)為b=5 的四分之一方形域進(jìn)行研究,求解域的左端與下端(圖中紅線邊界)約束法向位移,求解域上端、右端和圓孔邊界(圖中藍(lán)線邊界)根據(jù)解析應(yīng)力表達(dá)式(51)施加自然邊界條件.該問(wèn)題的誤差收斂率測(cè)試采用圖9 所示的112 個(gè)、403 個(gè)、1525 個(gè)、5929 個(gè)節(jié)點(diǎn)進(jìn)行離散,無(wú)網(wǎng)格近似采用三次基函數(shù),核函數(shù)的相對(duì)影響域?yàn)?.5.

圖8 帶孔無(wú)限大平板問(wèn)題模型Fig.8 Description of the plate with a hole problem

圖9 帶孔無(wú)限大平板問(wèn)題無(wú)網(wǎng)格離散Fig.9 Meshfree discretizations of the plate with a hole problem

圖10為位移和能量誤差收斂的結(jié)果,結(jié)果再次驗(yàn)證了所提赫林格-賴斯納變分原理法能夠保證理論誤差收斂率.同樣地,具有變分一致性的RKGSINitsche 也達(dá)到了理論誤差收斂率,而傳統(tǒng)高斯積分法即使在節(jié)點(diǎn)較為稀疏時(shí),也難以保證收斂效果.該問(wèn)題中,雖然RKGSI-LM 不能滿足變分一致性,但拉格朗日乘子法保持著較高的計(jì)算精度,其誤差并未影響該方法達(dá)到理論誤差收斂率.圖11為該問(wèn)題的效率對(duì)比,與同樣具有變分一致性的RKGSINitsche 方法相比,RKGSI-HR 無(wú)需計(jì)算傳統(tǒng)形函數(shù)導(dǎo)數(shù),具有更高的計(jì)算效率.最后,圖12為該問(wèn)題的應(yīng)力云圖,圖中清晰表明具有RKGSI 的計(jì)算精度要高于GI,GI 的應(yīng)力云圖均出現(xiàn)不同程度的振蕩.

圖10 帶孔無(wú)限大平板問(wèn)題誤差對(duì)比Fig.10 Error comparison for the plate with a hole problem

圖11 帶孔無(wú)限大平板問(wèn)題效率對(duì)比Fig.11 Efficiency comparison for the plate with a hole problem

圖12 帶孔無(wú)限大平板問(wèn)題 σxx應(yīng)力云圖Fig.12 Contour plot of stress σxxfor the plate with a hole problem

6 結(jié)論

本文以赫林格-賴斯納變分原理為基礎(chǔ),提出了一種滿足積分約束條件的變分一致高效本質(zhì)邊界條件施加方法.該方法采用混合離散近似赫林格-賴斯納變分原理弱形式中位移和應(yīng)力,其中位移采用傳統(tǒng)無(wú)網(wǎng)格形函數(shù)進(jìn)行離散,而應(yīng)力則在每個(gè)背景積分單元上近似為對(duì)應(yīng)階次的多項(xiàng)式.在赫林格-賴斯納變分原理的框架下,該方法的離散平衡方程具有與傳統(tǒng)尼茲法相類似的格式,可視為與再生光滑梯度積分法相配套的新型尼茲法.與傳統(tǒng)尼茲法相比,所提方法的修正變分項(xiàng)采用無(wú)網(wǎng)格形函數(shù)和再生光滑梯度進(jìn)行混合離散,在保證了變分一致性的同時(shí)避免了復(fù)雜耗時(shí)的形函數(shù)導(dǎo)數(shù)計(jì)算,明顯提高了計(jì)算效率;而對(duì)應(yīng)于尼茲法中的穩(wěn)定項(xiàng)則直接源于赫林格-賴斯納變分原理弱形式,無(wú)需額外增加穩(wěn)定項(xiàng),更重要的是穩(wěn)定項(xiàng)中不包含任何人工參數(shù),有效消除了尼茲法中的人工參數(shù)依賴性問(wèn)題.文中通過(guò)典型算例系統(tǒng)地驗(yàn)證了所提基于赫林格-賴斯納變分原理施加本質(zhì)邊界條件方法的變分一致性、計(jì)算精度和計(jì)算效率.

附錄A

本附錄中詳細(xì)推導(dǎo)了第4 小節(jié)式(34)中最后一個(gè)等式的推導(dǎo)過(guò)程,該等式為

引入關(guān)系式(20)～式(26),式(A1)中等式的左端各項(xiàng)簡(jiǎn)化過(guò)程如下所示

綜合式(A2)～式(A6),可得等式(A1)成立.需要注意的是,從式(A4)的化簡(jiǎn)過(guò)程可以看出,的表達(dá)式雖然直觀上看不對(duì)稱,但從原表達(dá)式中可知其實(shí)是對(duì)稱矩陣.

附錄 B

本附錄中列出了本文所提方法在數(shù)值實(shí)現(xiàn)過(guò)程中所采用的積分點(diǎn)位置和權(quán)重,其中附表B1和附表B2 分別為二次基函數(shù)情況和三次基函數(shù)情況的優(yōu)化數(shù)值積分方案.表中“”為KIJ,G所采用的數(shù)值積分點(diǎn),所采用的數(shù)值積分點(diǎn).并且 ξ,η和γ為三角形參數(shù)空間坐標(biāo),w和wB分別為三角形域內(nèi)積分權(quán)重和邊界積分權(quán)重.

表B1 二次基函數(shù)無(wú)網(wǎng)格法優(yōu)化的數(shù)值積分方案Table B1 The optimized quadrature rules for quadratic basis function

表B2 三次基函數(shù)無(wú)網(wǎng)格法優(yōu)化的數(shù)值積分方案Table B2 The optimized quadrature rules for cubic basis function