復合材料層合殼體熱彈耦合簡化模型

鐘軼峰，陳 磊，余文斌，張亮亮

（1．重慶大學a．土木工程學院；b．山地城鎮建設與新技術教育部重點實驗室，重慶400045；2．美國猶他州立大學 機械與航空航天工程系，洛根84322）

近20年來，先進復合材料結構因其高強度、高模量、可設計性等優點已廣泛應用于航天航空、機械、土木等領域。由于復合材料不同成分、不同纖維方向的熱膨脹系數相差很大，相比于各向同性材料對溫度的改變更加敏感，復合材料板／殼結構的熱響應已成為近年來復合材料力學熱點課題之一［1－2］，出現了各種復合材料板／殼理論。Ng等［3］應用古典層合理論（CLT）分析了復合材料圓柱殼的熱殘余應力和自由振動，由于忽略了橫向剪切效應，該理論僅對薄殼有效；為克服CLT的缺陷，Pradyumna等［4］基于一階剪切變形理論（FOSDT）分析了復合材料層合殼體在熱環境下非線性動態穩定性，但模擬彎曲時數值計算結果嚴重失真，產生“橫剪自鎖”現象，其原因在于總應變能中包含橫向剪切應變能的項在量級上不正確，沿厚度方向變化的剪切應變線性假設與上下表面剪應力為零自相矛盾，需人為地引入與幾何形狀和材料相關的剪切修正因子來解決這一矛盾；為提高沿厚向的應力／應變預測精度，各國學者紛紛把研究重點轉向了高階剪切變形理論，Li等［5］用高階殼理論分析了復合材料圓柱殼在熱環境下的非線性屈曲和后屈曲性能，Nosier等［6］基于層合理論確定濕熱環境下正交鋪層／角鋪層殼體的局部位移函數和層間應力，但仍無法準確預測面內的不連續性和沿厚度方向的橫向位移分量。

由上述分析可知，現有板殼熱彈性分析理論大多基于位移場假設，無法較好地反映層合板殼的三維效應和鋪層之間的相互作用。筆者嘗試采用變分漸近法［7－8］分析層合圓柱殼的單向耦合熱彈性問題。通過對降維模型近似能量中變分項的漸近修正，得到與原三維能量盡可能接近的近似能量，從而構建一種無需先驗性假設且便于工程應用的簡化殼體模型，為復合材料層合殼體在工程中的實際應用提供有價值的參考。

1 復合材料層合殼體三維能量表達式

基于Hamilton擴展原理，三維結構的彈性動力性能可表述為

式中：t1和t2是任意兩固定時間；K和U分別是動能和Helmholtz自由能；δˉW是外力所做虛功，上劃線表明不需要對該項精確變分。

圖1 層合殼體幾何構型及坐標系

由于僅考慮單向耦合熱彈性問題，因殼體變形引起的溫度改變忽略不計，得到無二次項的Helmholtz自由能泛函［9］：

式中：D為6×6階三維材料矩陣；α為6×1階三維熱膨脹系數列陣；T為相對于零應力狀態的溫差；為殼體幾何修正系數，其中kαβ為殼體面外曲率，因選用坐標系為圖1所示曲率線，k12＝k21＝0。

Γ為三維應變場［10］，可定義為

式中：為相應的算子；

其中：εαβ，καβ統稱二維廣義應變，其階數分別用n，n／h表示；w i為未知翹曲函數，若將參考面定義為殼體中面，則w i須滿足如下3個約束

殼體參考面上質點的動能和外力所做虛功可分別表示為

式中：μ為質量密度；U i為沿坐標軸x i的位移，其上的點表示對時間求導；＝－δU3；α為虛擬旋轉；廣義力f i和力矩mα可定義為

式中：p i，τi，βi分別為體力、殼體頂／底面表面力。

由式（1）、（6）和式（7），殼體參考面的動力性能可改寫為

由于式（9）中存在未知翹曲函數w i，若直接求解該式，將會遇到與求解原三維問題同樣的困難。常規做法是對翹曲場作先驗性假設，由于復合材料的各向異性和非均質性，這種假設可能會產生較大的誤差。變分漸近法可利用殼體固有的小參數對式（9）的變分項進行漸近分析得到w i，從而構建與原三維模型盡可能接近的二維殼體模型。

2 降維方法及近似能量推導

對殼體結構，可選擇h／l?1，h／R?1作為變分漸近計算所需的小參數，由此可估計各荷載階數為

式中μ為彈性材料常量的階數。

利用h／l，h／R可將三維能量漸近擴展為如下形式的系列二維近似能量泛函

為能處理體多層殼結構，并與二維有限元求解器相銜接，可將三維翹曲場離散為一維有限單元形式，

式中：S為形函數；V為沿橫法線方向的翹曲場節點值。

將式（11）代入式（9），得到離散形式能量泛函為

式中：L為荷載相關項；新引入的與幾何形狀和材料屬性有關的變量為

式（5）翹曲約束的離散形式可表示為

式中：H＝ ［STS］；ψ為零初始曲率E0的正交化核心矩陣，ψTHψ＝I。這樣，未知翹曲函數的求解問題轉化為式（14）約束下式（12）最小化問題。

2．1 零階近似

應用變分漸近法，需根據不同階數找到泛函的主導項。式（12）中與未知翹曲函數有關的零階近似主導項為

式中：E0，D hε0，Dεε0，αh0，αε0分別由式（14）中φ＝1（無幾何修正）定義的相關矩陣。

相應的零階翹曲函數為

將式（16）代入式（15）得到漸近修正到O（1）階的能量泛函為

這與如下形式的熱彈性古典層合模型相一致，

式中：A，NT分別為二維剛度矩陣和溫度產生的應力合力，其計算式為

2．2 一階近似

零階近似可用于分析薄殼的全局性能和面內分量，更高階近似可用來分析對中厚殼失效十分重要的面外應力和應變（σi3，Γi3），從而構建更精確殼體模型。

首先，將能量泛函漸近修正到O（h／R）階以考慮初始曲率效應，O（h／R）階翹曲因對能量無貢獻，可不必計算。得到

式中：

其中帶星矩陣由φ－1代替式（14）的φ定義。

其次，將能量漸近修正到O（h2／l2）以考慮橫向剪切變形。為此，將零階翹曲函數攝動為

將式（23）代回式（17），可得到一階近似的總能量泛函主導項為

式中：

相應的一階翹曲可求解為

最后得到漸近修正到O（h2／l2），O（h／R）的總能量泛函為

式中：

2．3 近似能量轉換及三維場重構

盡管式（27）能量泛函漸近修正到O（h2／l2），O（h／R）階，但因含有廣義應變的導數ε；α，難以直接應用。為得到實用的能量范函，可將式（27）轉換為工程中常用的Reissner－Mindlin模型形式。

在Reissner－Mindlin模型中有兩個附加橫向剪切應變Reissner－Mindlin模型應變量R與ε的關系可表示為

式中：

將式 （29）代 入 式 （27），可得到 由 Reissner－Mindlin應變量表示的能量泛函為

而實際應用的廣義Reissner－Mindlin模型形式為

為得到與式（32）等效的 Reissner－Mindlin模型，可通過以下2個彎矩平衡方程消除所有二維應變量的偏導數

由式（33），可將式（31）改寫為

式中：

若對任何R，U＊都趨于零，則可得到漸近修正Reissner－Mindlin殼體模型。對于一般各向異性殼體，該項往往并不為零，可通過最小二乘法等優化技術最小化U＊，等效Reissner－Mindlin模型的精確性取決于U＊趨近于零的程度。

由于降維模型的可靠性最終取決于其對三維場預測的精確度，因此還需提供重構關系以完善降維模型。由式（3），重構的三維應變場可表示為

三維應力場可使用材料本構關系得到，

3 算例

柱形彎曲問題已作為評估新提出的3D／2D分析模型精確性的基準。本章將所述理論和方法編制成變分漸近板／殼分析程序VAPAS，對圖2所示4層簡支復合材料層合殼體在熱／載荷下的柱形彎曲問題進行分析。

圖2 層合殼體結構示意圖

3．1 模型參數

殼體 各 層 傾 角 為 ［90°／0°／90°／0°］；殼 厚h＝1 mm，半徑R＝10 cm，夾角φ＝π／3；采用的坐標系為x1∈ ［0，φ］，x2∈ ［0，∞］，x3∈ ［－h／2，h／2］；材料為石墨／環氧復合材料，材料屬性為

殼體承受的溫度變化和正弦面荷載分別為

3．2 數值分析與討論

圖3（a）～（f）繪出了重構的沿厚度方向應力分布，并與一階剪切變形理論（FOSDT）、古典層合理論（CLT）和精確解［11－12］進行對比。應變和位移的變化趨勢和精度與應力相同，限于篇幅未在此繪出。由于應力分量σαβ，σ33分別是正弦和余弦函數，應力分布分別繪于x1＝π／6和x1＝π／3處。為方便比較，圖中縱橫坐標分別正則化為

由圖3可看出：VAPAS、FOSD和CLT都能較準確地預測面內應力分量σαβ分布，由于VAPAS能得到優化剪切剛度矩陣G，其結果比FOSDT和CLT更精確；對于橫向應力分量σα3，因特殊的鋪層設計（正交鋪設），CLT無法計算該值，而FOSDT因對位移場所做先驗性假設，誤差較大，VAPAS與精確解相一致，且重構的應力分布在不同層分界面處連續，這與大多數位移基殼理論完全不同。

圖3 4層簡支復合材料圓柱殼體在熱／載荷下沿厚度方向的應力分布

4 結論

1）基于變分漸近法構建了復合材料圓柱殼熱彈性簡化Reissner－Mindlin模型和重構關系，可考慮熱／載荷單向耦合效應，并使用4次多項式表示沿厚度方向的任意溫度分布，相較于假設溫度沿厚度方向線性分布（單層板理論）或沿層線性分布（層合理論）更符合實際情況。

2）在推導降維模型過程中，使用變分漸近法求解未知翹曲函數，不需任何動力學假設和剪切修正因子，也不同于傳統的板／殼理論假設翹曲場的一般形式，用高階翹曲作參數來求解假設函數中的未知參數的方法。

3）通過算例驗證：重構的沿厚度方向應力分量與三維精確解吻合很好，且簡化模型為等效單層殼模型，計算量與一階剪切變形理論相當，若殼體層數增加，其高效性更加顯著。

［1］Shen H S．Postbuckling of nanotube－reinforced composite cylindrical shells in thermal environments， Part II：pressure－loaded shells［J］．Composite Structures，2011，93（10）：2496－2503．

［2］Li Z M，Yang D Q．Thermal postbucking analysis of 3D braided composite cylindrical shells［J］．Journal of Mechanics，2010，26（2）：113－122．

［3］Ng T Y，Lam K Y，Reddy J N．Dynamic stability of cross－ply laminated composite cylindrical shells ［J］．International Journal of Mechanical Sciences，1998，40（8）：805－823．

［4］Pradyumna S，Gupta A．Nonlinear dynamic stability of laminated composite shells integrated with piezoelectric layers in thermal environment［J］．Acta Mechanica，2011，218（3／4）：295－308．

［5］Li Z M，Lin Z Q，Chen G L．Nonlinear buckling and postbuckling behavior of 3D braided composite cylindrical shells under external pressure loads in thermal environments［J］．Journal of Pressure Vessel Technology，2009，131（4）：258－269．

［6］Nosier A，Miri A K．Boundary－layer hygrothermal stresses in laminated，composite，circular，cylindrical shell panels［J］．Archive of Applied Mechanics，2010，80（4）：413－440．

［7］Pinhas B Y．New variational－asymptotic formulations for interlaminar stress analysis in laminated plates［J］．Journal of Applied Mathematics and Physics，1986，37（8）：305－321．

［8］Lee C Y．Dynamic variational asymptotic procedure for laminated composite shells：part I：low－frequency vibration analysis［J］．Journal of Applied Mechanics，2009，76（9）：110－122．

［9］Reddy J N．Mechanics of laminated composite plates and shells，theory and analysis［M］．2 ed．CRC Press，2004．

［10］鐘軼峰，余文斌．用變分漸近法進行復合材料層合板仿真及三維場重構［J］．復合材料學報，2010，27（4）：174－179．ZHONG Yifeng，YU Wenbin．Simulation and 3D field recovery of composite laminated plates by use of variational asymptotic method ［J］． Acta Materiae Compositae Sinica，2010，27（4）：174－179．

［11］Ren J G．Exact solutions for laminated cylindrical shells in cylindrical bending ［J］．Composite Science and Technology，1987，29（3）：169－187．

［12］Jing H S，Tzeng K G．Elasticity solution for laminated anisotropic cylindrical panels in cylindrical bending［J］．Composite Structures，1995，30（4）：307－317．