基于等幾何方法的功能梯度Mindlin板結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計(jì)

劉宏亮，姚力銘，王培金

(1.沈陽航空航天大學(xué)遼寧省飛行器復(fù)合材料結(jié)構(gòu)分析與仿真重點(diǎn)實(shí)驗(yàn)室，沈陽 110136；2.大連理工大學(xué)工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室，大連 116023)

0 引言

隨著高新技術(shù)的發(fā)展，飛行器通常需要在惡劣極端的環(huán)境下長時間工作[1]，對結(jié)構(gòu)和材料性能的要求越來越高。由于比剛度和比強(qiáng)度等方面的優(yōu)點(diǎn)，復(fù)合材料板在航空航天工業(yè)的應(yīng)用愈發(fā)普遍。與傳統(tǒng)的復(fù)合材料板不同，功能梯度復(fù)合板具有材料相之間平滑過渡且連續(xù)分布的特點(diǎn)，可以在很大程度上避免材料界面的出現(xiàn)對結(jié)構(gòu)力學(xué)性能的影響，例如降低應(yīng)力集中和減少殘余應(yīng)力。這實(shí)際上也是功能梯度結(jié)構(gòu)相比傳統(tǒng)復(fù)合材料結(jié)構(gòu)的一個優(yōu)勢[2-3]。

在航空航天產(chǎn)品設(shè)計(jì)過程中，利用優(yōu)化方法獲得輕質(zhì)高承載結(jié)構(gòu)是實(shí)現(xiàn)裝備輕量化和高性能化的一個重要且有效途徑。相比尺寸和形狀優(yōu)化，結(jié)構(gòu)拓?fù)鋬?yōu)化的設(shè)計(jì)潛力更大，對結(jié)構(gòu)的綜合性能提升也更加明顯，但同時其挑戰(zhàn)性也更大[4-5]。目前，各種拓?fù)鋬?yōu)化方法不斷發(fā)展，例如變密度方法[6]、漸進(jìn)結(jié)構(gòu)優(yōu)化方法[7]、水平集方法[8]、移動變形組件/孔洞方法[9-10]等。而且，相較于單相材料拓?fù)鋬?yōu)化，考慮多相材料的結(jié)構(gòu)拓?fù)鋬?yōu)化研究近年來受到了更多的關(guān)注[11-13]。然而，多材料結(jié)構(gòu)拓?fù)鋬?yōu)化需要一個合適的優(yōu)化模型，可以有效描述迭代過程中設(shè)計(jì)域內(nèi)不同材料的分布，同時應(yīng)該方便于設(shè)計(jì)優(yōu)化過程的靈敏度分析計(jì)算[14-16]。

功能梯度結(jié)構(gòu)的拓?fù)鋬?yōu)化可以是單材料優(yōu)化，也可以是多材料優(yōu)化。當(dāng)功能梯度結(jié)構(gòu)的材料漸變方式確定時，只優(yōu)化結(jié)構(gòu)的拓?fù)錁?gòu)型屬于單材料優(yōu)化，例如，材料沿著厚度或者面內(nèi)方向漸變且其分布確定的功能梯度板結(jié)構(gòu)拓?fù)鋬?yōu)化[17]。當(dāng)梯度材料漸變與結(jié)構(gòu)拓?fù)渫絻?yōu)化以實(shí)現(xiàn)功能梯度結(jié)構(gòu)性能最優(yōu)時，其本質(zhì)上是一種多材料結(jié)構(gòu)拓?fù)鋬?yōu)化[18-19]。與傳統(tǒng)多材料結(jié)構(gòu)拓?fù)鋬?yōu)化不同的是，這樣的設(shè)計(jì)在材料相之間具有平滑且連續(xù)分布的特點(diǎn)。值得指出，現(xiàn)有的功能梯度結(jié)構(gòu)拓?fù)鋬?yōu)化研究主要考慮平面問題，尚未考慮工程中普遍存在的板結(jié)構(gòu)設(shè)計(jì)。事實(shí)上，板結(jié)構(gòu)的多材料拓?fù)鋬?yōu)化設(shè)計(jì)研究也是近年來才受到關(guān)注[20-21]。因此，本文研究的重點(diǎn)是一種沿著面內(nèi)方向材料漸變的功能梯度板結(jié)構(gòu)拓?fù)渑c材料分布協(xié)同設(shè)計(jì)優(yōu)化問題。

對于板結(jié)構(gòu)的分析和設(shè)計(jì)，目前主要采用基于有限元方法的Kirchhoff板單元和Mindlin板單元。其中，Kirchhoff板單元適合薄板模型，而Mindlin板單元適合中厚板模型。薄板模型的適用范圍有限，因此Mindlin板單元在工程中的應(yīng)用范圍更廣。然而，Mindlin板單元進(jìn)行有限元分析時可能會存在數(shù)值閉鎖問題，嚴(yán)重影響數(shù)值計(jì)算的精度。為了實(shí)現(xiàn)更具工程應(yīng)用性的板結(jié)構(gòu)拓?fù)鋬?yōu)化，克服數(shù)值求解難題，本文基于等幾何分析[22]發(fā)展一種功能梯度Mindlin板結(jié)構(gòu)的拓?fù)鋬?yōu)化設(shè)計(jì)方法。采用具有高階連續(xù)性的非均勻有理B樣條(non-uniform rational B-spline，NURBS)基函數(shù)作為功能梯度板結(jié)構(gòu)建模、分析和設(shè)計(jì)的基礎(chǔ)，獲得了梯度材料與結(jié)構(gòu)拓?fù)鋮f(xié)同優(yōu)化設(shè)計(jì)方案。

1 功能梯度Mindlin板等幾何拓?fù)鋬?yōu)化設(shè)計(jì)

1.1 Mindlin板結(jié)構(gòu)等幾何分析

本文考慮的Mindlin板基于一階剪切變形理論，假設(shè)板中面的法向纖維在變形過程中保持平直，長度不變，但不一定保持垂直于中面。因此，橫向剪切變形不可忽略，Mindlin板的應(yīng)變能由下式給出

(1)

式中，σf和εf是彎曲應(yīng)力和應(yīng)變，σc和εc是橫向剪切應(yīng)力和應(yīng)變，α為剪切修正系數(shù)。線彈性應(yīng)力應(yīng)變關(guān)系可定義為

σf=Dfεf

σc=Dcεc

(2)

其中，Df和Dc分別表示彎曲項(xiàng)與剪切項(xiàng)的彈性常數(shù)矩陣，即

(3)

式中，E是彈性模量，ν是泊松比，G則是剪切模量。將式(2)和式(3)帶入Mindlin板應(yīng)變能式中可得

(4)

基于等幾何分析，利用NURBS曲面表示結(jié)構(gòu)的幾何模型，即

(5)

(6)

式中，Ni，p和Mj，q分別表示定義在節(jié)點(diǎn)矢量{ξ0，ξ1，…，ξm1}和{η0，η1，…，ηm1}上的B樣條基函數(shù)，wi，j是控制點(diǎn)對應(yīng)的權(quán)重系數(shù)。

在Mindlin板結(jié)構(gòu)等幾何分析中，NURBS基函數(shù)不僅用于結(jié)構(gòu)幾何模型的定義，同時也作為等參元的形函數(shù)。離散后的Mindlin板線彈性控制方程可表示為

KU=F

(7)

式中，U為位移矢量，K表示結(jié)構(gòu)的整體剛度陣，它由單元剛度陣Ke組裝而成，表示如下

(8)

式中，Bf和Bc分別為彎曲項(xiàng)和剪切項(xiàng)的應(yīng)變位移矩陣，h表示板的厚度，ui，j為結(jié)構(gòu)的位移場。

本文利用等幾何分析的k細(xì)分方案可以實(shí)現(xiàn)NURBS基函數(shù)跨單元高階的Cp-1次連續(xù)，在很大程度上可以緩解閉鎖作用的影響。對于較薄的板模型，將采用NURBS等幾何分析單元結(jié)合縮減積分的策略解決數(shù)值閉鎖問題。

為了構(gòu)建一個基于連續(xù)密度的梯度漸變結(jié)構(gòu)優(yōu)化問題，式(8)中的彈性系數(shù)矩陣近似為不同材料相的彈性常數(shù)，這將在下面討論的插值方案中實(shí)現(xiàn)。

1.2 功能梯度拓?fù)鋬?yōu)化的插值方案

在功能梯度結(jié)構(gòu)中，材料之間不存在明確的界面，材料性能平滑變化。借鑒基于SIMP的變密度方法插值方案，假設(shè)板由兩種材料構(gòu)成，彈性模量分別為E1和E2，通過設(shè)計(jì)域內(nèi)某一點(diǎn)上兩種材料所占的百分比，即體積分?jǐn)?shù)來表示這種連續(xù)變化材料的性質(zhì)。

功能梯度材料性質(zhì)的計(jì)算通常采用數(shù)值均勻化法和混合法。本文基于混合法，假設(shè)兩種材料均勻混合，沒有微觀結(jié)構(gòu)，如圖1所示。目前，被廣泛采用的幾種混合模型包括Voigt模型、Halpin-Tsai復(fù)合模型、Voigt-Reuss界限和Hashin-Shtrikman界限[23]。本文研究重點(diǎn)是采用等幾何方法實(shí)現(xiàn)優(yōu)化設(shè)計(jì)，因此選擇使用各向同性復(fù)合材料的Hashin-Shtrikman上界和下界的平均值作為材料模型。

圖1 結(jié)構(gòu)中某一點(diǎn)材料混合物

當(dāng)兩種材料的泊松比等于1/3時，Hashin-Shtrikam的上限和下限可簡化為楊氏模量條件

(9)

當(dāng)y→1時，Hashin-Shtrikam模型上下限均為E1，當(dāng)y→0時，Hashin-Shtrikam模型上下限均為E2。假設(shè)E1>E2，當(dāng)y=1時，此時結(jié)構(gòu)中某一處全為強(qiáng)材料，當(dāng)y=0時即為弱材料。當(dāng)結(jié)構(gòu)中一點(diǎn)是兩種材料均勻混合時，采用Hashin-Shtrikam模型上下限均值近似表示為

(10)

借鑒SIMP方法，各向同性均質(zhì)材料的性質(zhì)被描述為E(x)=xpE0，設(shè)計(jì)變量為0時表示此處為優(yōu)化所產(chǎn)生的孔洞。為了實(shí)現(xiàn)材料屬性和結(jié)構(gòu)拓?fù)涞膮f(xié)同設(shè)計(jì)，功能梯度材料的插值模型被定義為

(11)

式中，p為懲罰系數(shù)。值得指出，對于式(11)中的設(shè)計(jì)變量y并不需要懲罰，因?yàn)樵O(shè)計(jì)目標(biāo)是獲得光滑漸變的功能梯度結(jié)構(gòu)。為了避免數(shù)值奇異，插值格式中采用一個足夠小的量綱為1的數(shù)Emin，即

(12)

1.3 功能梯度結(jié)構(gòu)拓?fù)渑c材料協(xié)同優(yōu)化列式

考慮體積約束下的功能梯度Mindlin板結(jié)構(gòu)最小柔順度優(yōu)化問題，當(dāng)功能梯度材料是由兩種材料(例如金屬與陶瓷)構(gòu)成時，其對應(yīng)的拓?fù)鋬?yōu)化列式可以表示為

minimize:J(xi,yi)=FTU(xi,yi)

(13)

式中，J表示目標(biāo)函數(shù)即結(jié)構(gòu)的柔順度，xi和yi表示控制點(diǎn)處的設(shè)計(jì)變量，分別對應(yīng)總材料的體分比和強(qiáng)材料的體分比，Ve表示NURBS單元體積，V0表示結(jié)構(gòu)總體積，fvx和fvy分別為材料用量的限制分?jǐn)?shù)，即體積分?jǐn)?shù)。fvx限制結(jié)構(gòu)的總材料用量，fvy則控制有材料處的強(qiáng)材料占比。ym和yM為給定設(shè)計(jì)變量yi的上下限，由式(9)可知，若ym和yM取0和1則對應(yīng)為Hashin-Shtrikam模型的上下限，為了保證實(shí)現(xiàn)強(qiáng)弱材料的光滑漸變，本文中ym和yM分別設(shè)置為0.1和0.9。

對于本文基于等幾何方法的功能梯度Mindlin板結(jié)構(gòu)拓?fù)鋬?yōu)化，用于幾何建模和結(jié)構(gòu)分析的NURBS基函數(shù)也用于表示相對密度場，即

(14)

其中，Θ表示設(shè)計(jì)變量x和y組成的集合。由于NURBS基函數(shù)的非負(fù)性和規(guī)范性，當(dāng)設(shè)計(jì)變量在給定的設(shè)計(jì)區(qū)間變化時，所對應(yīng)的相對密度分布也會滿足規(guī)定的區(qū)間分布，從而建立起合理的控制點(diǎn)密度設(shè)計(jì)變量和材料相對密度分布之間的映射關(guān)系。此外，NURBS基函數(shù)的高階連續(xù)性可以避免棋盤格模式等數(shù)值問題[24]。

2 設(shè)計(jì)優(yōu)化靈敏度分析

本文采用基于梯度的移動漸進(jìn)線法[25](MMA)來更新設(shè)計(jì)變量，進(jìn)而求解建立的優(yōu)化列式(13)。這里需要計(jì)算目標(biāo)函數(shù)和約束函數(shù)的靈敏度。功能梯度板的幾何建模、結(jié)構(gòu)分析和材料分布場采用統(tǒng)一的參數(shù)化，NURBS基函數(shù)的特性保證了解析靈敏度的有效計(jì)算。目標(biāo)函數(shù)對于設(shè)計(jì)變量xe的靈敏度可表示為

(15)

式中單元剛度陣對設(shè)計(jì)變量xe的導(dǎo)數(shù)表達(dá)式為

(16)

根據(jù)鏈?zhǔn)角髮?dǎo)法則可得到

(17)

式(17)中Mindlin板單元剛度陣由兩項(xiàng)組成，分別為彎曲項(xiàng)與剪切項(xiàng)，則?E/?xe的計(jì)算表達(dá)式為

(18)

由式(15)和式(18)可得

(19)

目標(biāo)函數(shù)對于設(shè)計(jì)變量ye的靈敏度表達(dá)式為

(20)

式中，?Ke/?ye可以參考式(16)和式(17)，而彈性模量對設(shè)計(jì)變量ye的一階導(dǎo)數(shù)為

(21)

式中，E′(y)表示Hashin-Shtrikam界限對設(shè)計(jì)變量ye求導(dǎo)，可以表示為

(22)

根據(jù)式(21)和式(22)可得

(23)

另外，體積約束的靈敏度可以直接求導(dǎo)得出。

3 數(shù)值算例與討論

通過兩個數(shù)值算例證明本文方法可以有效獲得體積約束下的功能梯度Mindlin板結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計(jì)。優(yōu)化算例中采用二次NURBS基函數(shù)，為避免數(shù)值閉鎖現(xiàn)象對Mindlin板結(jié)構(gòu)等幾何分析精度的影響，采用縮減的兩點(diǎn)積分策略。剪切修正系數(shù)設(shè)置為α=5/6，材料用量限制比值為0.5。當(dāng)連續(xù)兩次迭代的目標(biāo)值相對變化小于0.01或達(dá)到最大迭代次數(shù)50時，優(yōu)化過程終止。算例中材料、幾何以及載荷參數(shù)選擇為量綱為1的參數(shù)。

3.1 四邊固定的方板

第一個算例考慮一個四邊固定的矩形方板(a×a=1×1)，在板中心處施加一個集中載荷(F=1)，如圖2所示。板的厚度為0.1，參數(shù)設(shè)置參考了文獻(xiàn)[20]，用紅色表示強(qiáng)材料E1=4，藍(lán)色表示弱材料E2=2，兩種材料泊松比為0.3?？傮w積分?jǐn)?shù)設(shè)為0.5，強(qiáng)材料體積分?jǐn)?shù)占比V1從0.4下降到0.1，弱材料體積分?jǐn)?shù)占比V2從0.1上升到0.4，J表示目標(biāo)函數(shù)即柔順度。將設(shè)計(jì)域細(xì)分為100×100的單元時，優(yōu)化結(jié)果如圖3所示。

(a)V1=0.4，V2=0.1，J=30.52

在相同厚度和邊界條件下，可以觀察到材料填充面積與體積分?jǐn)?shù)具有明顯的相關(guān)性，隨著硬質(zhì)材料(紅色)的體積分?jǐn)?shù)逐漸下降，結(jié)構(gòu)的柔順度逐漸增加?？梢酝茢喑觯趯Y(jié)構(gòu)剛度貢獻(xiàn)較小的地方，軟質(zhì)材料傾向于替代硬質(zhì)材料，而強(qiáng)材料總是分配在中心區(qū)域(集中載荷作用區(qū))和應(yīng)力集中區(qū)。4種不同強(qiáng)材料體積分?jǐn)?shù)所對應(yīng)的優(yōu)化迭代曲線如圖4所示。

圖4 4種不同強(qiáng)材料體積分?jǐn)?shù)下的迭代曲線

圖3中的4種優(yōu)化結(jié)果，拓?fù)湓O(shè)計(jì)構(gòu)型幾乎一致，不同點(diǎn)在于材料分布，即強(qiáng)弱材料之間的過渡區(qū)域。圖3右側(cè)顏色欄中的0.1和0.9對應(yīng)式(13)設(shè)定的Hashin-Shtrikam模型上下限。本文協(xié)同優(yōu)化體現(xiàn)在兩種設(shè)計(jì)變量分別控制結(jié)構(gòu)拓?fù)渑c材料配比，如圖5所示。右側(cè)“0”表示結(jié)構(gòu)中的孔洞，即圖5黑白設(shè)計(jì)中的白材料。由于圖3中可視化點(diǎn)的數(shù)量有限，功能梯度材料的漸變細(xì)節(jié)表示得并不是很清楚。由圖4可知，其收斂速度較快，在20多步就已經(jīng)收斂。在不改變設(shè)計(jì)優(yōu)化計(jì)算效率的情況下，為了更清楚地顯示梯度漸變細(xì)節(jié)，在圖3優(yōu)化結(jié)果的基礎(chǔ)上，將用于描述最終優(yōu)化結(jié)果的可視化點(diǎn)擴(kuò)大100倍，如圖6所示，可見優(yōu)化結(jié)果具有清晰的拓?fù)湓O(shè)計(jì)和平滑的材料過渡區(qū)域。

圖5 體積約束下的功能梯度板協(xié)同優(yōu)化設(shè)計(jì)

(a)強(qiáng)材料體積分?jǐn)?shù)為0.4

3.2 四邊簡支的方板

第2個算例考慮一個四邊簡支的矩形板，厚度為0.01，設(shè)計(jì)域參考圖2，板的中心仍然施加一個集中載荷。強(qiáng)材料(紅色相)的彈性模量E1=2，弱材料的彈性模量E2=1，兩種材料泊松比為0.3，總體積分?jǐn)?shù)為0.5，其中強(qiáng)材料占比為0.3。該算例的設(shè)計(jì)參考了文獻(xiàn)[20]在多材料薄板拓?fù)鋬?yōu)化設(shè)計(jì)中的算例。將設(shè)計(jì)域離散為50×50的單元，最終的優(yōu)化結(jié)果如圖7所示，目標(biāo)函數(shù)的迭代曲線如圖8所示。

(a)結(jié)構(gòu)拓?fù)?/p>

圖8 四邊簡支功能梯度板的迭代曲線

功能梯度結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計(jì)也可以視作一種多材料結(jié)構(gòu)拓?fù)鋬?yōu)化設(shè)計(jì)，但是傳統(tǒng)多材料結(jié)構(gòu)拓?fù)湓O(shè)計(jì)的不同材料相之間存在明顯的界面。圖7(a)中的優(yōu)化設(shè)計(jì)構(gòu)型與參考文獻(xiàn)[20]中單材料結(jié)構(gòu)設(shè)計(jì)類似，但通過功能梯度材料結(jié)構(gòu)協(xié)同設(shè)計(jì)優(yōu)化，最終結(jié)構(gòu)的目標(biāo)函數(shù)值為137 319.91，在材料參數(shù)設(shè)置一致的情況下，結(jié)構(gòu)性能優(yōu)于參考文獻(xiàn)中的雙材料板拓?fù)鋬?yōu)化設(shè)計(jì)，目標(biāo)函數(shù)下降了11.28%，證明了本文設(shè)計(jì)方案的有效性和優(yōu)越性。另外，可視化點(diǎn)的數(shù)量不依賴建模、分析和設(shè)計(jì)變量的數(shù)量，可以在保證功能梯度材料相之間連續(xù)且平滑漸變的同時兼顧設(shè)計(jì)優(yōu)化的效率。

4 結(jié)束語

本文研究了一種在等幾何分析框架下的功能梯度Mindlin板結(jié)構(gòu)拓?fù)鋬?yōu)化方法實(shí)現(xiàn)材料漸變分布與結(jié)構(gòu)拓?fù)鋮f(xié)同設(shè)計(jì)。采用基于NURBS的建模、分析和設(shè)計(jì)參數(shù)化方法來實(shí)現(xiàn)梯度材料結(jié)構(gòu)的一體化優(yōu)化設(shè)計(jì)。等幾何方法中NURBS基函數(shù)用于梯度材料分布拓?fù)鋬?yōu)化使得設(shè)計(jì)結(jié)果可以自然避免拓?fù)鋬?yōu)化的常見數(shù)值不穩(wěn)定問題，高階連續(xù)性可以促進(jìn)材料分布的光滑性以及設(shè)計(jì)優(yōu)化的靈敏度分析，因此建立了合理的設(shè)計(jì)變量與材料分布的關(guān)系。數(shù)值算例表明，本文方法具有較好的收斂性，在不需要顯著增加計(jì)算單元和設(shè)計(jì)參數(shù)的情況下能夠有效獲得具有清晰細(xì)節(jié)特征的結(jié)構(gòu)拓?fù)渑c材料漸變的協(xié)同優(yōu)化設(shè)計(jì)。通過與已有研究結(jié)果對比，本文優(yōu)化設(shè)計(jì)方法展現(xiàn)了有效性和優(yōu)越性。