基于ABAQUS的載人潛水器觀察窗結(jié)構(gòu)蠕變行為分析

杜青海，江海港，胡曉康

1 上海海洋大學 深淵科學工程技術(shù)研究中心, 上海 201306

2 上海海洋大學 海洋科學學院, 上海 201306

3 上海海洋大學 工程學院, 上海 201306

0 引 言

目前，陸上資源的日益匱乏讓全世界的目光開始轉(zhuǎn)向海洋資源。隨著對海洋資源的深入研究，各國對海洋的關(guān)注重點開始從近淺海逐漸向深遠海發(fā)展[1]，而載人潛水器（HOV）在其中將發(fā)揮關(guān)鍵性作用。觀察窗是載人潛水器的關(guān)鍵結(jié)構(gòu)之一，其安全性對于人員生命安全以及潛水器性能具有十分重要的作用。相對于金屬艙體結(jié)構(gòu)而言，透明玻璃材質(zhì)觀察窗結(jié)構(gòu)的設(shè)計更具復雜性和不確定性[2]，故國際上的通用規(guī)范[3]采取提高安全系數(shù)(4～16倍)的方法來確保結(jié)構(gòu)安全。

國內(nèi)外學者對于觀察窗的研究始于上世紀中期，其中最為著名的是美國的Stachiw[4]，其畢生都在進行該方面的研究，成果已被寫入設(shè)計規(guī)范[3]并受到國際上相關(guān)行業(yè)的認可。值得注意的是，Stachiw的研究主要集中在7 000 m潛深范圍以內(nèi)，對于更大的深度，其采取的是增大安全系數(shù)按比例進行保守設(shè)計的方法。近10年來，我國載人深潛裝備隨著7 000 m級“蛟龍”號的成功研制，也得到了蓬勃發(fā)展，4 500 m級“深海勇士”號以及萬米級的“奮斗者”號已相繼完成研制。針對載人潛水器的觀察窗結(jié)構(gòu)，我國科研人員在結(jié)構(gòu)應力[5]、蠕變特性[6-7]和結(jié)構(gòu)協(xié)調(diào)性[8]等方面開展了研究。對于全海深大深度觀察窗結(jié)構(gòu)，筆者所在課題組通過數(shù)值方法，針對觀察窗結(jié)構(gòu)規(guī)范方法開展了探索性研究[9]，指出該設(shè)計方法過于保守；隨后，又針對觀察窗結(jié)構(gòu)進行了優(yōu)化計算，以尋找最優(yōu)方案[10]。但對于觀察窗的長時間周期和蠕變壽命等問題，尚未開展相關(guān)研究。

因此，本文擬采用ABAQUS軟件構(gòu)建錐臺型觀察窗的有限元計算模型，然后基于窗玻璃材料（丙烯酸塑料材料）的基本性能確定蠕變模型參數(shù)，進而分析、驗證觀察窗結(jié)構(gòu)的蠕變特性，最后，基于有限元蠕變模型計算分析不同加載速度下觀察窗結(jié)構(gòu)的蠕變特性及結(jié)構(gòu)性能，為觀察窗結(jié)構(gòu)的疲勞壽命研究奠定基礎(chǔ)。

1 錐臺型觀察窗

1.1 觀察窗幾何結(jié)構(gòu)參數(shù)

錐臺型觀察窗的基本結(jié)構(gòu)參數(shù)如下：承受的外載荷P、厚度δ、玻璃低壓面（即內(nèi)側(cè)圓面）直徑Di、窗座內(nèi)口徑Df、錐角α，以及厚?徑比δ/Di和窗座的懸伸比Di/Df。如圖1所示，取錐臺型觀察窗的內(nèi)側(cè)圓面為AB面，外側(cè)高壓圓面為EF面，選擇EF面的圓心O點作為柱坐標系rOz的原點，觀察窗軸對稱旋轉(zhuǎn)方向，也即錐臺觀察窗圓面的圓周方向記為θ，AB面的圓心為C點，OF長記為r。

圖1 錐臺型觀察窗的基本結(jié)構(gòu)參數(shù)Fig.1 Basic structural parameters of conical observation window

1.2 材料的蠕變特性

由于潛水器觀察窗大多數(shù)使用丙烯酸塑料（PMMA）[4]，這種材料具有比較明顯的蠕變特性，能使應力保持在較低水平，且不會導致嚴重的應力集中現(xiàn)象，有力地推動了觀察窗的飛躍式發(fā)展。本文錐臺型觀察窗主要結(jié)構(gòu)材料的機械物理屬性[10]如表1所示。

表1 丙烯酸塑料的材料力學參數(shù)[10]Table 1 Mechanical parameters of PMMA[10]

圖2給出了丙烯酸塑料試件在不同外載荷作用下拉伸后的應變ε隨時間t的變化曲線[11]，揭示了此窗玻璃材料固有的蠕變特性。

圖2 PMMA拉伸試件試驗應變隨時間的變化曲線[11]Fig.2 Variation of strain level of PMMA with time in tensile specimen experimental test[11]

2 錐臺型觀察窗蠕變有限元理論

2.1 蠕變理論模型

由于采用傳統(tǒng)的ABAQUS線彈性觀察窗模型得出的結(jié)果忽略了丙烯酸塑料蠕變這一特性，導致計算結(jié)果與真實的試驗數(shù)據(jù)間存在較大差異，無法作為普遍的計算方法，因此需要尋找一種能體現(xiàn)蠕變特性的有限元仿真計算方法。為此，本文將在進一步研究觀察窗蠕變特性的基礎(chǔ)上，采用單步靜力?線彈性與多步?非線性相結(jié)合的分析步控制方法，同時引入加載和保壓的時間，重點得出觀察窗的軸向位移、應力和應變等隨時間的變化規(guī)律。

ABAQUS軟件中有3種蠕變模型：雙曲正弦模型、應變硬化模型和時效硬化模型，可用來仿真模擬材料的蠕變等特性。其中，恒溫與定負載條件下的蠕變特性選擇時效硬化模型或應變硬化模型，變溫條件下的蠕變特性則使用雙曲正弦法則模型。應變硬化的基本思想是：在蠕變發(fā)生的全過程中，起到強化作用的主要原因是蠕變變形與時間無任何關(guān)系，故其適用于短時間的試驗。然而，載人潛水器在工作深度下長時間作業(yè)的觀察窗結(jié)構(gòu)蠕變特性[6-7]與時間有很大的關(guān)系，故觀察窗結(jié)構(gòu)蠕變分析可選擇時效硬化模型。

時效硬化理論[12]的公式可以表述為：當外界溫度一定時，材料的應變ε、應力σ和時間t這3個參數(shù)滿足的狀態(tài)函數(shù)Φ為：

式中，ε˙為應變隨時間變化的應變率。

由文獻[12]可知，在蠕變發(fā)生階段，影響蠕變的是時效、擴散及恢復這幾個參數(shù)，且因蠕變曲線在第1、第2階段具有幾何相似的性質(zhì)，故

其中，應力函數(shù)f1(σ)符合下面的冪函數(shù)形式，即

式中，n為應力冪函數(shù)的指數(shù)參數(shù)。

而時間響應函數(shù)f2(t)則用Ω(t)函數(shù)表示，即

式中，A，m分別為時間冪函數(shù)的系數(shù)和指數(shù)參數(shù)。

因此，將式(3)和式(4)代入式(2)，可得到時效硬化模型下應變的數(shù)學表達式，即

2.2 蠕變模型參數(shù)

根據(jù)對觀察窗結(jié)構(gòu)進行的線彈性數(shù)值模擬仿真分析與優(yōu)化[4,9-10]，可知規(guī)范規(guī)定的錐臺型觀察窗結(jié)構(gòu)所對應的最大工作壓力作用下的整體應力水平基本集中在30～45 MPa之間。鑒此，為得到式(5)中的A，m，n這3個待定系數(shù)，根據(jù)圖2所示的試驗數(shù)據(jù)[11]，盡量選擇高應力條件下的試驗值，即選取30～40 MPa間的材料拉伸蠕變曲線。并且，選擇在拉伸試驗曲線的兩端取點，對拉伸試件試驗值與計算值進行對比，如圖3(a)和圖3(b)所示。由圖可見，在這2種取點情況下所得的計算值與試驗值相差較大，故應在試驗值的中間部分取點，并需滿足二者之間誤差很小，在此取點情況下的試驗值與計算值對比如圖3(c)所示。由圖可見，擬合出的計算值與試驗值基本吻合，證明該組參數(shù)有效。具體取點方式如表2所示。

表2 確定參數(shù)的取點方式及其解Table 2 Point selecting schemes and solution of determined parameters

圖3 PMMA拉伸試件應變試驗值與計算值對比Fig.3 Comparison of strain level between expreimental test and calculation of PMMA tensile specimen

2.3 有限元模型

根據(jù)錐臺型觀察窗結(jié)構(gòu)和均勻靜水壓力載荷的軸對稱性特點，對觀察窗結(jié)構(gòu)選取軸對稱形式的二維平面模型進行仿真模擬計算。根據(jù)模型尺寸，將選取大小合適的CAX4R單元對模型進行網(wǎng)格劃分，并設(shè)置模型部件之間的接觸形式、窗體與窗座之間的摩擦系數(shù)和罰函數(shù)因子來描述接觸[10]。觀察窗仿真計算模型的載荷邊界約束分布及網(wǎng)格單元劃分如圖4所示。

圖4 ABAQUS有限元模型的載荷邊界約束與網(wǎng)格劃分Fig.4 Boundary constraints and grid division of finite element model on ABAQUS

2.4 有限元計算方法驗證

現(xiàn)有的行業(yè)標準ASME PVHO是建立在Stachiw[4]的系列試驗模型研究基礎(chǔ)之上的，故選取其典型模型尺寸，即主尺度參數(shù)。具體如下：Di=25.4 mm，δ/Di= 0.75，1.0，1.25，α=90°。根據(jù)真實試驗情況將計算壓力加載到規(guī)范規(guī)定[3]的界限壓力68.9 MPa，然后保壓100 h。

在觀察窗結(jié)構(gòu)數(shù)值模擬仿真研究[9-10]中，已針對上述尺度的模型驗證了線彈性計算環(huán)境下，摩擦系數(shù)μ=0.05時有限元計算結(jié)果隨網(wǎng)格粗細變化的準確性及收斂性，得出在該尺度下觀察窗和窗座分別采用1 和3 mm的網(wǎng)格比較合適，故此處不再贅述。

運用ABAQUS軟件對上述尺度的模型建立幾何模型并進行模擬仿真分析，獲得加壓和保壓階段的應變與變形等，然后將數(shù)值模擬仿真所得內(nèi)側(cè)圓面中心軸向位移fc及其與試驗結(jié)果的相對誤差Δ進行對比分析，結(jié)果如圖5所示。由圖5可見，在前100 h，有限元計算結(jié)果與試驗結(jié)果[4]吻合，誤差絕對量在10%以內(nèi)；在100～1 000 h間，有限元計算結(jié)果稍大于試驗結(jié)果，二者之間的誤差保持在11%左右，這說明在此階段采用數(shù)值模擬仿真方法具有保守性?？梢姡疚姆抡嬗嬎憬Y(jié)果與試驗結(jié)果較相符。而實際工程應用中保壓的時間一般都小于100 h，觀察窗體結(jié)構(gòu)因其玻璃材料的黏彈特性，在保壓的起始階段(0～2 h)其結(jié)構(gòu)性能主要由材料的黏彈屬性決定[4]，故本文將著重選擇10～100 h這一階段進行研究。

圖5 有限元計算結(jié)果與試驗結(jié)果對比Fig.5 Comparison of finite element calculation and experimental test results

3 全海深觀察窗結(jié)構(gòu)蠕變分析

3.1 高壓作用下的蠕變特性

現(xiàn)有行業(yè)規(guī)范[3]中的設(shè)計方法是通過安全系數(shù)法則來確保結(jié)構(gòu)疲勞壽命安全，但卻無法揭示長周期載荷下觀察窗所具有的蠕變現(xiàn)象，導致計算結(jié)果與實測結(jié)果存在較大誤差[7,9,10]，且不利于結(jié)構(gòu)的優(yōu)化。因此，本文將基于時效硬化模型和試驗確定參數(shù)的蠕變計算方法深入研究觀察窗結(jié)構(gòu)在高壓作用下的蠕變等特性?！膀札垺碧栍^察窗的結(jié)構(gòu)尺寸參數(shù)和最大工作壓力恰好是規(guī)范[3]中方法——曲線設(shè)計參數(shù)分界[9-10]的典型結(jié)構(gòu)尺寸參數(shù)和壓力范圍，其基本幾何結(jié)構(gòu)參數(shù)如表3所示。

表3 “蛟龍?zhí)枴庇^察窗主尺度參數(shù)[5,7]Table 3 Main dimensions of conical observation window of Jiaolong ship[5,7]

利用第2節(jié)的數(shù)值模擬仿真方法對觀察窗結(jié)構(gòu)進行有限元建模和模擬仿真，并在計算過程中分別選擇線彈性和時效硬化的蠕變模型。在數(shù)值計算過程中，同時采用單步靜力?線彈性與多步?非線性相結(jié)合的分析步控制的方法加入時間非線性響應。最終，將作用在觀察窗模型上的壓力加載至71.6 MPa并保壓10 h，觀察窗結(jié)構(gòu)軸向位移的計算結(jié)果，如圖6所示。圖中，Uz為結(jié)構(gòu)剖面的軸向位移。

圖6 保壓10 h后觀察窗的軸向位移云圖Fig.6 Axial displacement contours of conical observation window after holding pressure for 10 h

圖6中給出了觀察窗結(jié)構(gòu)在保壓10 h后線彈性和蠕變這2種模型下的軸向位移云圖，而圖7則給出了這2種模型在加載和保壓過程中觀察窗內(nèi)側(cè)圓面中心點軸向位移fc隨時間的變化曲線。由圖7可見，采用線彈性模型計算出的軸向位移在保壓階段無任何改變，而采用蠕變模型計算出的軸向位移fc在保壓階段則是隨著時間的增長逐步增大，這直接反映出了觀察窗結(jié)構(gòu)蠕變這一特性，與實際結(jié)果比較契合。圖8進一步展示了最大應變εmax隨時間的變化趨勢。由圖8可見，在保壓階段，采用蠕變模型計算的最大應變隨著時間的推移逐漸增大，而采用線彈性模型計算的最大應變在保壓階段則無任何變化，可見觀察窗玻璃材料的固有蠕變特性直接導致了觀察窗結(jié)構(gòu)在高壓下的蠕變特性。

圖7 觀察窗軸向位移隨時間的變化曲線Fig.7 Variation of the axial displacement of conical observation window with time

圖8 觀察窗最大應變隨時間的變化曲線Fig.8 Variation of the maximum strain level of conical observation window with time

根據(jù)對深海觀察窗結(jié)構(gòu)強度優(yōu)化的研究[10]，發(fā)現(xiàn)在主尺度參數(shù)下，錐角對錐臺型觀察窗的結(jié)構(gòu)優(yōu)化與安全評估具有重要意義。因此，本文選取與表3所示觀察窗主尺度參數(shù)一致的4組結(jié)構(gòu)模型進行分析，即Di=220 mm，δ/Di=1.0，α=70°，80°，90°和100°，將作用壓力加載至71.6 MPa并保壓100 h，計算得到4組模型結(jié)構(gòu)的基本力學特性，如圖9～圖12所示。圖中， σ為 最大等效應力，τmax為最大剪應力。

由圖9可見，在壓力載荷與摩擦系數(shù)都相同的情況下，相同厚?徑比觀察窗內(nèi)側(cè)圓中心處的軸向位移fc與錐角α成反比；由圖10可見，觀察窗同一時間處的最大等效應力 σ隨著錐角α的增大呈現(xiàn)先增大后減小的趨勢；由圖11可見，觀察窗同一時間的最大剪應力τmax隨著錐角α的增大呈增大的趨勢；而圖12則揭示了觀察窗同一時刻的最大應變εmax是隨著錐角α的增大呈先增大后減小的規(guī)律。綜合上述分析，可得錐角α=70°時觀察窗盡管在軸向位移方面其數(shù)值是4個模型中最大的，但其最大等效應力、最大剪應力和最大應變卻是最小的，且應力集中現(xiàn)象不明顯。因此，α=70°時的觀察窗模型是相對優(yōu)化的結(jié)果，這與線彈性數(shù)值模擬優(yōu)化研究結(jié)論[10]一致。

圖9 軸向位移隨時間的變化趨勢Fig.9 The changing trend of axial displacement with time

圖10 最大等效應力隨時間的變化趨勢Fig.10 The changing trend of maximum equivalent stress with time

圖11 最大剪應力隨時間的變化趨勢Fig.11 The changing trend of maximum shear stress with time

圖12 最大應變隨時間的變化趨勢Fig.12 The changing trend of maximum strain with time

3.2 不同加載速度下的蠕變特性

經(jīng)過大量的仿真計算，發(fā)現(xiàn)不同的加載（下潛）速度對觀察窗的位移、應變及最大等效應力等輸出值均有影響，會影響結(jié)構(gòu)的蠕變特性。

為便于研究，針對3.1節(jié)中優(yōu)化結(jié)果，設(shè)Di為固定值220 mm，α=70°，δ/Di=1.0，摩擦系數(shù)μ=0.05，采取3種不同的加載速度，即0.67，2和6 MPa/min，然后將觀察窗模型的作用壓力P加載至71.6 MPa并保壓100 h，最后研究采取不同加載速度時對觀察窗結(jié)構(gòu)基本力學性能的影響。

圖13和圖14給出了加載和保壓過程中觀察窗結(jié)構(gòu)的位移和最大應變特性。由圖可見，在將作用壓力加載至71.6 MPa的過程中，3種加載速度下的位移和應變值相差無幾，只是隨著加載時間的縮短，二者略有降低；而在保壓階段二者差距卻十分明顯，盡管保壓時間相同，但隨著加載時間的增加，位移和最大應變在保壓階段增長的幅度反而降低了。圖15顯示出在整個加載過程中不同加載速度下的最大等效應力值幾乎一致，其變化主要體現(xiàn)在保壓階段，即隨著加載時間的增大，最大等效應力的降幅在減小。由圖15可見，在相同的長時間保壓作用下，延長加載時間相對縮短加載時間而言，與正常加載速度（例如2 MPa/min）相比，觀察窗模型結(jié)構(gòu)計算的等效應力變化幅度受到的影響較大。綜上所述，壓力加載速度不宜過快，否則不利于高壓保壓階段的結(jié)構(gòu)蠕變應變及位移釋放。

圖13 不同加載速度下位移的變化曲線Fig.13 Variation of displacement under different loading speeds

圖14 不同加載速度下最大應變的變化曲線Fig.14 Variation of maximum strain under different loading speeds

圖15 不同加載速度下最大等效應力的變化規(guī)律Fig.15 Variation law of the maximum equivalent stress under different loading speeds

3.3 觀察窗蠕變特性下的結(jié)構(gòu)強度

為了更好地研究觀察窗材料蠕變特性下其強度的力學結(jié)構(gòu)屬性，需進一步分析觀察窗在長時間保壓作用下的結(jié)構(gòu)強度特性，這對觀察窗的疲勞壽命分析來說意義重大。因此，選取前述優(yōu)化的觀察窗結(jié)構(gòu)主尺度參數(shù)，即Di=220 mm，δ/Di=1.0，α=70°，再按2 MPa/min加載速度將作用壓力加載至71.6 MPa，并保壓100 h。

圖16(a)所示為保壓后觀察窗上表面（EF面）在3個方向的應變曲線。由圖可見，在EF面中心附近存在高應變集中區(qū)域，數(shù)值上呈現(xiàn)從中心向四周遞減的趨勢，圖中虛線范圍所表示的區(qū)域隨著外部加載的不斷增大，出現(xiàn)了“凹坑”情況[4]。圖16(b)所示為內(nèi)側(cè)圓面（AB面）在3個方向（即圖1所示柱坐標系下）的應變曲線。由圖可見，內(nèi)側(cè)圓面在邊緣處易出現(xiàn)因幾何結(jié)構(gòu)和應變突變所導致的裂紋等。圖中 ，εr，εθ和 εz分別為3個坐標軸方向上的主應變。

圖16 觀察窗外/內(nèi)表面的應變曲線Fig.16 The strain distribution of outside/inside surface of conical observation window

圖17所示為觀察窗保壓后的3/4立體應變云圖。從圖中可以直觀地看出，EF面中心O點附近3個應變分量的絕對值明顯高于四周，這也進一步揭示了在高壓下觀察窗上表面會出現(xiàn)“凹坑”的破損失效現(xiàn)象。

圖17 觀察窗結(jié)構(gòu)的應變云圖Fig.17 The strain contours of conical observation window structure

圖18給出了在保壓過程中，在0，10，100 h這3個時刻觀察窗中心軸線OC上沿z軸3個方向的應變變化情況，圖中縱坐標為軸線OC位置坐標z與厚度δ的比值，橫坐標為應變值大小。從圖中可以清晰地看出，虛線區(qū)域內(nèi)的應變量會隨著保壓時間的增加而持續(xù)增加，相比其他中心軸線OC上區(qū)域，該區(qū)域受到的應變突變更嚴重，也更容易出現(xiàn)損傷破壞。

圖18 觀察窗中心軸線OC的應變曲線Fig.18 The strain distribution along OC at the center axis of conical observation window

4 結(jié) 論

本文基于丙烯酸塑料材料的蠕變特性，采用ABAQUS軟件對錐臺型觀察窗結(jié)構(gòu)進行了數(shù)值模擬仿真計算，并通過時效硬化模型對材料和觀察窗結(jié)構(gòu)的蠕變特性展開了分析，主要得到如下結(jié)論：

1) 基于ABAQUS軟件，構(gòu)建了用于評估錐臺型觀察窗在長周期載荷下材料與結(jié)構(gòu)蠕變特性的可靠數(shù)值方法。

2) 觀察窗結(jié)構(gòu)有限元線彈性模型與蠕變模型的比較分析，以及筆者所在課題組之前對蠕變模型與典型參數(shù)模型試驗的對比分析表明，蠕變模型下的計算結(jié)果更符合試驗情況。通過多參數(shù)優(yōu)化計算分析觀察窗結(jié)構(gòu)的蠕變特性，并在特定主尺度參數(shù)下優(yōu)化對比錐角對結(jié)構(gòu)性能的影響，結(jié)果顯示錐角α=70°時觀察窗結(jié)構(gòu)的蠕變特性更優(yōu)。

3) 不同加載（下潛）速度下觀察窗的蠕變強度研究顯示，為避免蠕變特性下的結(jié)構(gòu)損傷，不宜采取過快的加載速度。

本文針對觀察窗結(jié)構(gòu)的蠕變特性分析可為今后觀察窗結(jié)構(gòu)的優(yōu)化設(shè)計提供參考。