球形封頭開孔接管結構屈曲特性研究

董乃逸，賀小華

(南京工業大學 機械與動力工程學院，南京 211816)

0 引言

球殼封頭結構簡單，受力均勻，廣泛應用于高壓容器、深海耐壓潛水器等承壓結構[1]。因工藝需要，耐壓球殼需開設各類接管結構，開孔接管破壞了球殼結構的整體連續性，在外壓作用下，導致球殼承壓能力下降，易產生屈曲失穩破壞。對球形封頭開孔接管結構屈曲特性及臨界失穩壓力進行研究意義重大。

對于外壓球殼的屈曲行為研究由來已久。早期針對理想球殼，ZOELLY[2]基于小變形理論，推導了外壓完整球殼的臨界失穩壓力公式。由于Zoelly公式計算得到的結果與試驗值有較大差距，KARMAN[3]基于外壓球殼的軸對稱失穩問題，提出了大變形理論公式以計算外壓球殼的臨界失穩壓力，結果與試驗值相對較為接近。20世紀60年代，KRENZKE等[4]對200多個具有不同缺陷的耐壓球殼進行試驗，提出了臨界失穩壓力的計算公式。目前，國內外標準規范：GB/T 150—2011《壓力容器》、ASME Ⅷ-1和ASME Ⅷ-2，以及EN 13445-3:2014分別基于上述小變形理論及修正公式，提出了臨界失穩壓力的計算公式。文獻[5]分析描述了球殼在均勻外壓下的前屈曲狀態，使用有限元計算球殼的臨界失穩壓力和屈曲模式，顯示了球殼的后屈曲行為；文獻[6]研究了外壓半球封頭的彈塑性屈曲，提出了預測半球封頭臨界失穩壓力的解析公式，并與數值結果進行對比。

由于球形封頭開孔接管結構在其開孔處的應力狀態較為復雜，難以通過解析方法求解球形封頭開孔接管結構臨界失穩壓力。文獻[7]基于外壓球殼開孔區域的軸對稱性和等面積補強法，對大深度潛水器載人球殼開孔加強形式進行幾何變換，完善了球殼開孔加強理論計算方法，并對球殼進行外壓試驗，驗證了計算方法的有效性和適用性。文獻[8]運用有限元軟件ANSYS，分析了側壁加強、墊板加強、肘板加強這三種耐壓球殼開孔加強方式的極限強度并進行對比，結果表明側壁加強的極限強度最高；同時考慮到結構性能和工藝性方面，提出了梯形加強這一新型的開孔加強形式。文獻[9]運用有限元軟件ANSYS對潛水器開孔耐壓殼半球封頭承載能力進行非線性分析，討論了開孔個數、孔徑大小、圍壁厚度3個參數對耐壓殼半球封頭承載力的影響，結果表明，開孔越多、孔徑越大、圍壁厚度越小都會使球形封頭承載能力下降。

本文系統討論開孔率di/Di、接管與封頭厚度比δet/δe及徑厚比Di/δe這3個無量綱結構參數對球形封頭開孔接管結構臨界失穩壓力Pcr及屈曲特性的影響規律，為外壓球形封頭開孔接管結構穩定性研究及工程設計提供參考。

1 有限元模擬計算

1.1 材料性能參數

球形封頭材料采用Q345R，接管材料采用16Mn，常溫下的材料性能參數如表1所示。

表1 常溫下材料性能參數

1.2 有限元計算方案

有限元分析結構為球形封頭徑向開孔接管結構，如圖1所示。球形封頭厚度取定值δe=10 mm，取筒體直邊段長度L1=25 mm，接管外伸長度L2=200 mm，通過改變開孔率di/Di、接管與封頭厚度比δet/δe和封頭徑厚比Di/δe的值進行有限元模擬計算。di/Di分別取值0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8；δet/δe分別取值0.2,0.4,0.6,0.8,1.0,1.25,1.5,1.75,2.0；Di/δe分別取值50,100,150,200。有限元分析模型共計288組。

1.3 材料本構關系

文獻[10]指出本構關系對殼體結構的臨界失穩壓力影響甚少，因此,本文選取雙線性材料本構模型進行有限元分析，塑性段斜率取E/100，其本構關系如圖2所示。

圖1 球形封頭徑向開孔接管結構

圖2 雙線性材料本構模型

1.4 單元類型及邊界條件

分析結構選擇Shell 181殼單元，該單元計算精度已在文獻[10-11]中證實。在滿足網格無關性要求下進行網格劃分，分析結構的邊界條件為封頭直邊段施加固支約束，開孔接管端面施加環向約束以保證接管端面的圓度，封頭與接管均受外壓作用。球形封頭開孔接管有限元模型見圖3。

圖3 封頭開孔接管有限元模型

2 有限元模型驗證

2.1 完整球殼臨界失穩壓力值比較

對球殼屈曲行為進行的研究歷程是從最初的不考慮缺陷，到后來發現初始缺陷對其屈曲耐壓特性有巨大的影響。實際工程中半球封頭是沖壓成型的，由于模具和工藝不可避免地會造成一定的缺陷，從而使半球封頭的實際形狀和完美的半球封頭存在一定的差距。本文選取文獻[12]中的4個試驗模型，模型參數分別為球殼直徑D=117.54,117.98,117.54,117.5 mm，計算厚度分別為0.422,0.423,0.406,0.415 mm，彈性模量E=1.93×105MPa，泊松比μ=0.28，屈服強度σs=205 MPa。分別從理論公式、試驗結果和有限元數值模擬三方面，對外壓球殼的臨界失穩壓力進行計算對比。

小變形理論公式是基于完美幾何形狀且無缺陷的球殼，其臨界失穩壓力的計算公式[2]：

Pcr=1.21E(δe/R)2

(1)

式中，E為材料彈性模量,MPa；δe為球殼計算厚度,mm；R為球殼半徑,mm。

大變形理論公式基于結構的非線性因素影響，其臨界失穩壓力的計算公式[3]：

Pcr=0.36516E(δe/R)2

(2)

本文有限元數值模擬時施加的缺陷方法有兩種，一種是采用“一致缺陷模態法”[13]設置初始幾何缺陷(一致缺陷模態法是基于最低階屈曲模態，從而模擬結構的初始缺陷分布),通過ANSYS中的UPGEOM命令里的FACTOR值對模型統一施加10%的幾何缺陷，然后進行非線性求解；另一種是施加局部厚度缺陷。GB/T 150.4—2011《壓力容器 第4部分：制造、檢驗和驗收》中6.4.2規定球形封頭的最大形狀偏差外凸不超過1.25%Di，內凹不超過0.625%Di。文獻[14-15]表明內凹0.625%Di下臨界失穩壓力低于外凸1.25%Di下臨界失穩壓力，因此本文模型采用形狀偏差內凹0.625%Di的球形封頭模型，并在封頭頂部最大壁厚減薄位置施加單個缺陷厚度減薄15%進行非線性分析。最終求解提取出位移最大點的載荷-位移圖，并根據零曲率準則得到臨界失穩壓力值Pcr。

將文獻[12]試驗結果、公式(1)解、公式(2)解及本文有限元特征值、有限元非線性解1(一致缺陷模態法)以及有限元非線性解2(施加局部厚度缺陷)列于表2，并將其數據作成圖4，以便于直觀比較。

表2 球形封頭臨界失穩壓力值Pcr比較

圖4 球形封頭臨界失穩壓力值Pcr比較

由表2和圖4可以看出，有限元特征值與小變形公式(1)解結果較為接近，由于這兩種方法均未考慮幾何非線性以及缺陷的存在，導致其結果遠大于其他解；大變形公式(2)解略大于試驗值，計算結果偏于保守；有限元非線性解1及有限元非線性解2兩者均考慮了結構的幾何大變形和初始幾何缺陷，因此其值均與試驗結果較為接近。考慮到一致缺陷模態法結果的有效性及操作的方便性，后續有限元模型均采用此法設置初始幾何缺陷。

2.2 開孔球殼臨界失穩壓力結果比較

選用文獻[16]中按1∶10比例自制的不銹鋼縮比試驗模型，模型結構如圖5(a)所示。該模型參數為：球殼外半徑R0=100 mm，球殼厚度t=6.3 mm，接管內徑d1=50 mm，接管外徑d=77 mm，接管厚度t1=13.5 mm，內伸高度h1=32.2 mm，外伸高度h2=32.2 mm，彈性模量E=1.95×105MPa，泊松比μ=0.3，屈服強度σs=300 MPa。試驗測得的縮比模型臨界失穩壓力為20.5 MPa，采用一致缺陷模態法，施加10%的初始缺陷進行有限元非線性求解，有限元模型如圖5(b)所示，得到臨界失穩壓力為19.223 MPa，與試驗值比較，誤差為6.64%。兩者誤差基本在工程允許范圍內，由此也說明本文非線性屈曲分析方法的有效性。

(a)開孔球殼圍壁結構示意

(b)開孔球殼圍壁有限元模型

3 球形封頭開孔接管結構有限元非線性屈曲分析

3.1 臨界失穩壓力削弱系數K

為了分析討論開孔接管結構對球形封頭臨界失穩壓力的影響，通過引入削弱系數K來反映開孔接管的作用，K值為球形封頭開孔接管結構的臨界失穩壓力Pcr與球形封頭未開孔結構的臨界失穩壓力Pcrs(有限元非線性模擬值)之比，即K=Pcr/Pcrs。

3.1.1 接管與封頭厚度比δet/δe與K值關系

圖6示出了不同徑厚比Di/δe及不同開孔率di/Di下，臨界失穩壓力削弱系數K與δet/δe的分布圖。可以看出，不同徑厚比及不同開孔率下，K值與δet/δe變化規律基本相同；同一δet/δe下，K值隨di/Di的增大而減小；當δet/δe較小時，K值明顯偏低，隨著δet/δe增大，K值顯著增大，當δet/δe≥1.0后，K值隨δet/δe增大變化趨于平緩。

(a)Di/δe=50

(b)Di/δe=100

(c)Di/δe=150

(d)Di/δe=200

3.1.2 開孔率di/Di與K值關系

圖7示出了不同徑厚比Di/δe及不同接管與封頭厚度比δet/δe下，臨界失穩壓力削弱系數K與開孔率di/Di的分布圖。可以看出，Di/δe與δet/δe一定時，K值隨著di/Di增大而減小；當δet/δe較小時，K值隨di/Di增大顯著降低；而當δet/δe較大時，K值隨di/Di增大、降低幅度較小；當δet/δe≥1.0時，同一di/Di下，δet/δe對臨界失穩壓力削弱系數K的影響較小。

(a)Di/δe=50

(b)Di/δe=100

(c)Di/δe=150

(d)Di/δe=200

3.2 計算模型失穩模態分析

由圖6可以看出，臨界失穩壓力削弱系數K與接管與封頭厚度比δet/δe變化有個轉折點，當δet/δe≥1.0時，K值隨δet/δe增大變化不顯著。同樣由圖7可以看出，當δet/δe≥1.0時，同一開孔率下，不同計算模型臨界失穩壓力削弱系數K值較為接近。為進一步討論球形封頭開孔接管結構屈曲特性，圖8示出徑厚比Di/δe=150時部分計算模型的失穩模態圖。

圖8 球形封頭中心開孔接管結構失穩模態圖(Di/δe=150)

由圖8(a)可以看出，較小開孔率(di/Di=0.4)，當接管較薄(δet/δe=0.2)時，計算模型失穩位置位于接管；當δet/δe增大到0.6時，由圖8(b)可知，計算模型失穩位置由接管轉移到封頭-接管連接處；當δet/δe增至2.0時，由圖8(c)可知，計算模型失穩位置位于封頭處。

由圖8(d)～(f)可看出，當接管較薄(δet/δe=0.2)時，計算模型失穩位置位于接管，隨著接管壁厚增大，計算模型失穩位置逐漸由接管轉移到封頭。由此，解釋了圖6中臨界失穩壓力削弱系數K與δet/δe分布圖中轉折點的由來。當接管與封頭厚度比δet/δe低于圖6中轉折點時，接管較薄，承壓能力較弱，計算模型失穩位置位于接管，此時分析結構臨界失穩壓力較小；當δet/δe增大到一定值時，接管承壓能力超過封頭，計算模型失穩位置轉到封頭，此時分析結構臨界失穩壓力較大，進一步由圖6可以看出，隨著開孔率di/Di增大，失穩位置由接管轉到封頭所需接管與封頭厚度比越大。本文計算參數范圍內，當δet/δe≥1.0時，計算模型失穩位置均位于封頭上。

3.3 外壓球形封頭開孔接管結構臨界接管與封頭厚度比(δet/δe)c

由圖6可以看出，臨界失穩壓力削弱系數K與δet/δe變化存在轉折點。當δet/δe低于轉折點對應值時，計算模型失穩位置位于接管；而當δet/δe高于轉折點對應值時，計算模型失穩位置位于封頭。本文將失穩位置由接管轉向封頭的最小接管與封頭厚度比稱為臨界接管與封頭厚度比 (δet/δe)c。表3列出本文計算參數范圍下的臨界接管與封頭厚度比(δet/δe)c。

表3 臨界接管與封頭厚度比(δet/δe)c

由表3可以看出，臨界接管與封頭厚度比(δet/δe)c是3組定值(0.6,0.8,1.0)，其與徑厚比及開孔率之間的關聯式如下：

(3)

表3中列出的臨界接管與封頭厚度比(δet/δe)c基本對應于圖6中曲線的轉折點。由圖6可以看出，當δet/δe<(δet/δe)c時，臨界失穩壓力削弱系數K明顯較小，結構承載能力較弱；當δet/δe≥(δet/δe)c時，K值相對較大，且K值隨δet/δe變化不顯著。實際工程設計中，建議參考表3中給出的臨界接管與封頭厚度比，取δet/δe=(δet/δe)c，可以在滿足結構輕量化設計基礎上，獲得較優的結構承載能力。

3.4 臨界失穩壓力削弱系數K值與結構參數的經驗關系式

基于臨界失穩壓力削弱系數K隨Di/δe，di/Di，δet/δe的變化曲線，進一步擬合出K值與3個結構參數的經驗關系式。由于所研究參數δet/δe變化范圍(0.2～2.0)較大，且K值隨δet/δe的變化可分為兩個階段。為了保證擬合精度，采用分段擬合方法。當δet/δe位于[0.2,0.6]區間內，采用以下關系式進行擬合：

(4)

式中，a1～a7為待定系數。

基于有限元模擬所得數據點，最終擬合結果為：

(5)

當δet/δe位于(0.6,2.0]區間內，基于有限元模擬所得數據點，最終擬合結果為：

(6)

為了驗證所得關系式的可靠性，進一步給出當開孔率di/Di=0.2,0.5,0.8時，利用以上經驗關系式所得擬合K值與實際有限元K值數據點對比圖，如圖9所示。圖9中補充計算了δet/δe=0.7時相應有限元模型的臨界壓力削弱系數K，并與第二區間段公式擬合值進行比較。由圖9可以看出，不同區間內擬合所得K值與實際有限元K值兩者較為接近，且擬合曲線能夠反映出K值隨不同結構參數的變化趨勢。

(a)Di/δe=50

(b)Di/δe=100

(c)Di/δe=150

(d)Di/δe=200

計算表明,除了δet/δe=0.2時個別點誤差較大，其他點公式擬合值和有限元模擬值誤差較小，所有分析計算點平均誤差小于10%。

4 結論

分析討論3個結構參數：徑厚比Di/δe、開孔率di/Di、接管與封頭厚度比δet/δe對外壓球形封頭開孔接管結構臨界失穩壓力Pcr的影響。研究結果如下。

(1)球形封頭臨界失穩壓力削弱系數K值隨著開孔率di/Di增大而減小。

(2)接管與封頭厚度比較小時，臨界失穩壓力削弱系數K值越小，計算模型失穩位置位于接管，當接管與封頭厚度比增至一定值時，K值顯著增大，計算模型失穩位置位于封頭。

(3)開孔率越大，失穩位置由接管轉到封頭所需接管與封頭厚度比越大。本文計算參數范圍內，當δet/δe≥1.0時，計算模型失穩位置均位于封頭上。

(4)提出了臨界接管與封頭厚度比(δet/δe)c，實際工程設計中，當接管與封頭厚度比取臨界接管與封頭厚度比時，可以在滿足結構輕量化設計基礎上，獲得較優的結構承載能力。

(5)建立了臨界失穩壓力削弱系數K值與不同結構參數的經驗關系式，對比發現所得關系式預測精度較高。