充氫工藝對復(fù)合材料儲氫氣瓶殘余應(yīng)力的影響

江 勇，許 明，危書濤

(合肥通用機械研究院有限公司，合肥，安徽 230031)

由于可再生、零排放、高能量密度等優(yōu)點，氫能被認為是21 世紀綠色能源，吸引了世界各國的關(guān)注。儲氫是關(guān)系到氫能源能否廣泛應(yīng)用的關(guān)鍵技術(shù)，常用的儲氫技術(shù)包括高壓儲氫、低溫液態(tài)儲氫和金屬氫化物儲氫，其中高壓儲氫易于實現(xiàn)且成本低，成為了研究的焦點。由于具有高比強度和比剛度的優(yōu)點，復(fù)合材料氣瓶被認為是最有可能成為氫燃料電池汽車用儲氫載體。

為了實現(xiàn)快速充氫，充氫時間通常為3 min～5 min。由于負焦耳-湯普森系數(shù)、動能轉(zhuǎn)化為熱能和氫氣壓縮效應(yīng)，充氫過程溫度會上升。當溫度高于85 ℃，會導致復(fù)合材料性能下降，甚至汽車零部件損壞，引起安全事故。因此，國際標準ISO 15869 和行業(yè)標準SAE J2579?2018 要求充氫過程中，氣瓶溫度低于85 ℃。

為了確保復(fù)合材料儲氫氣瓶的安全性，國內(nèi)外很多工作研究復(fù)合材料氣瓶充氫溫升問題。Yang[1]建立了熱動力分析模型，研究儲氫瓶充氫過程中的溫升。Cristina 等[2]采用Computational Fluid Dynamics (CFD)模型，研究IV 型氣瓶不同工況下充氫到70 MPa 后的溫度分布。Melideo 等[3]研究采用預(yù)冷方式充氫和泄氫的溫度變化。一些研究者[4]通過試驗研究了充氫過程中的溫度變化。文獻[5]建立了氫氣泄放模型，該模型考慮了真實氣體效應(yīng)。浙江大學鄭津洋[6]研究了III 型復(fù)合材料儲氫氣瓶快速充氫的溫升問題，發(fā)現(xiàn)氣瓶尾部溫度較高。文獻[7]研究了氣瓶初始氫氣質(zhì)量和充氫時間對氣瓶溫度分布的影響。文獻[8]研究了不同體積和周圍溫度以及汽車初始溫度對溫升的影響。文獻[9]研究了幾何尺寸對溫度升高的影響，發(fā)現(xiàn)小長徑比會加劇溫升效應(yīng)，大充氫口徑會降低溫升效應(yīng)。文獻[10]研究了充氫流量、氣瓶初始壓力和環(huán)境溫度對最高溫度的影響，結(jié)果表明隨著充氫流量和周圍環(huán)境溫度的增加，最高溫度會提高，隨著初始壓力的提高，最高溫度會下降。文獻[11]報道了一種控制系統(tǒng)，該系統(tǒng)通過控制充氫速率來降低充氫成本和充氫時間。文獻[12 ? 14]研究了不同的參數(shù)和因素對氣瓶溫升的影響規(guī)律。文獻[15]研究了車用III 型和IV 儲氫氣瓶的溫升問題，比較了兩種氣瓶溫升效應(yīng)。文獻[16]研究了在火源熱輻射作用下，氣瓶的溫升效應(yīng)，并提出了降低溫度的防控措施。文獻[17]發(fā)展了一種簡化的解析模型，研究復(fù)合材料圓柱體的溫度特性，基于該解析模型，建立了三維熱-機械耦合有限元模型研究復(fù)合材料圓柱體的熱機械耦合行為。然而，該模型過于簡單，且做了很多與實際快速充氫過程不符的假設(shè)。

綜上所述，學者們已經(jīng)開展了不少關(guān)于復(fù)合材料氣瓶充氫過程溫升效應(yīng)的研究。但是，溫度分布不均勻，且復(fù)合材料和內(nèi)膽材料熱膨脹系數(shù)的差別導致熱機械耦合效應(yīng)存在，熱機械耦合效應(yīng)的存在必然會影響結(jié)構(gòu)的應(yīng)力、應(yīng)變分布[18?19]，目前這方面的研究工作較少。因此，有必要研究復(fù)合材料儲氫氣瓶快速充氫過程中的熱機械耦合行為。本文提出流-熱-固耦合分析方法，研究充氫速率、充氫位置和氣瓶尺寸對氣瓶應(yīng)力分布的影響。本文分成以下幾部分：第一部分介紹流-熱-固分析方法；第二部分研究充氫速率、充氫位置和幾何尺寸對氣瓶應(yīng)力分布的影響；第三部分是總結(jié)。

1 流-熱-固耦合分析方法

本文提出的流-熱-固耦合分析方法分成兩步：第一步，建立二維軸對稱CFD 模型，模擬復(fù)合材料氣瓶充氫過程溫度升高；第二步，建立三維熱機械耦合模型，以第一步CFD 模型計算獲得的溫度場作為輸入條件，研究氣瓶應(yīng)力分布。由于CFD 模型和三維熱機械耦合模型的網(wǎng)格差異，需要通過線性插值將溫度場施加到熱機械耦合模型中。

1.1 CFD 模型

本文是以實驗室設(shè)計、纏繞和固化而成的鋁合金內(nèi)膽碳纖維纏繞增強復(fù)合材料氣瓶為研究對象，復(fù)合材料氣瓶實物如圖1(a)所示。該氣瓶的設(shè)計工作壓力為70 MPa，內(nèi)膽采用鋁合金T6061旋壓成型，復(fù)合材料層采用T700 碳纖維增強環(huán)氧樹脂纏繞獲得，分為環(huán)向纏繞和螺旋纏繞，其中環(huán)向纏繞角度為90°，螺旋纏繞角度為±16°，鋪層順序為[902/±16]21/902，筒身纏繞的復(fù)合材料單層厚度0.31 mm，鋁合金內(nèi)膽筒身厚度8 mm，氣瓶總體長度為1950 mm，體積為150 L，內(nèi)膽外徑為350 mm，如圖1(b)所示?？紤]到溫升效應(yīng)會導致氣瓶壓強升高，按照標準要求，充氫時的目標壓強為工作壓力的1.25 倍，即87.5 MPa。

圖 1 本文采用的復(fù)合材料氣瓶及內(nèi)膽尺寸Fig. 1 The composite overwrapped pressure vessel and the dimension of the liner in this work

為了節(jié)約計算時間，同時氣瓶具有軸對稱性，本文采用商用軟件ANSYS Fluent 17.0 建立二維軸對稱CFD 分析模型(如圖2 所示)，分析快速充氫過程中的溫升效應(yīng)，這一方法也為其他工作所采用。整個模型包括三部分：內(nèi)部氫氣、內(nèi)襯和外層復(fù)合材料。模型通過有限體積法保守形式求解Navier-Stokes 方程，采用二階隱式。利用共軛傳熱模型，模擬固體材料與流體之間熱傳導與溫度變化的耦合關(guān)系。復(fù)合材料與外界環(huán)境的導熱系數(shù)為 4.5 W·m?2·K?1。氣瓶和內(nèi)部初始溫度為20 ℃，周圍環(huán)境溫度在整個充氫過程保持不變，為20 ℃。氫氣源的溫度為?40 ℃，氣瓶內(nèi)部初始壓力為0.1 MPa，為了準確模擬氣瓶內(nèi)部壓強從0.1 MPa 變到87.5 MPa，模型中采用Redlich-Kwong真實氣體方程，k-ω 湍流模型。模型中鋁合金T6061 和復(fù)合材料的物理屬性如表1 所示。由于該模型中熱量從內(nèi)襯往復(fù)合材料傳導，即主要沿著徑向方向傳導，該方向熱傳導系數(shù)主要由樹脂決定，因此本文中復(fù)合材料的熱傳導系數(shù)認為各向同性，這與論文[20]的假設(shè)一致。

圖 2 CFD 模型網(wǎng)格圖Fig. 2 The mesh of the CFD model

表 1 材料物理性能Table 1 Physical properties of different materials

CFD 的計算結(jié)果將作為溫度場輸入到熱機械耦合模型，影響熱機械應(yīng)力的準確性，因此CFD模型的準確性至關(guān)重要。為了驗證CFD 模型的準確性，將其計算結(jié)果與解析解[21]進行對比，解析解見式(1)：

式中：p2、T2分別為氣瓶充氫后的壓力和平均溫度；cv2為充氫后氫氣的等體積比熱容；p1、T1和cv1分別為氣瓶初始壓力、溫度和氫氣的等體積比熱容；p0、T0和cp0分別為氫氣源的壓力、溫度和等壓比熱容；α=1.9155×10?6為常數(shù)。

CFD 的計算結(jié)果與解析解對比如圖3 所示，從圖中可以看出，CFD 計算結(jié)果與解析解吻合度較高，最大誤差低于3K，證明本文采用的CFD 模型有效。

圖 3 氣瓶內(nèi)平均溫度隨著壓強變化而變化曲線Fig. 3 The variation of the average temperature in COPV during fast filling with the pressure

1.2 熱-機械耦合模型

圖 4 熱機械耦合模型及應(yīng)變測試點位置Fig. 4 Thermal-mechanical model and the location of the strain measurement point

表 2 復(fù)合材料力學和熱物理性能Table 2 Mechanical property and thermal physical property of composite materials

表 3 鋁合金T6061 力學和熱物理性能Table 3 Mechanical property and thermal physical property of aluminum alloy T6061

圖 5 線性插值方法獲得熱機械耦合模型溫度場示意圖Fig. 5 Comparison between the experimental result and the thermal-mechanical model result

為了驗證熱機械耦合模型的正確性，本文開展試驗測試充氫過程中復(fù)合材料氣瓶封頭最外層16°螺旋纏繞層上一點在纖維方向應(yīng)變，其位置如圖4 應(yīng)變測試點所示，并與熱機械耦合模型計算結(jié)果進行對比，如圖6 所示，從圖中可以看出熱機械耦合模型給出的結(jié)果與試驗測量結(jié)果吻合度較高，證明了模型的正確性。

圖 6 數(shù)值模擬和試驗測量的纖維方向應(yīng)變隨著內(nèi)壓變化曲線Fig. 6 The variation of micro-strain at the measurement point with the inner pressure by experiment and numerical method

2 各因素的影響

從其他研究工作可以看出，充氫速率、充氫位置和幾何尺寸會對復(fù)合材料氣瓶溫升產(chǎn)生影響。由于熱機械耦合效應(yīng)，溫度場將影響應(yīng)力場。因此，本文將研究這些因素對復(fù)合材料氣瓶的應(yīng)力分布的影響。

2.1 充氫速率

通常充氫過程將在3 min～5 min 內(nèi)完成，因此本文分別研究了3 min 和5 min 的充氫時間進行對比，充氫結(jié)束后，兩者的溫度場如圖7 所示。內(nèi)襯內(nèi)表面溫度沿著軸向分布如圖8 所示，從圖中可以看出，5 min 充氫時間更長，因此溫度要略微高點，特別是氣瓶右側(cè)封頭處。

圖 7 不同充氫速率下的復(fù)合材料氣瓶云圖Fig. 7 Temperature contour of COPV at different rates

圖 8 不同充氫速率下內(nèi)襯溫度沿軸向分布Fig. 8 Temperature distribution in the liner along the axis direction at different rates

以上兩種情況的內(nèi)襯Mises 等效應(yīng)力云圖如圖9 所示，內(nèi)襯的應(yīng)力沿著軸向分布如圖10 所示，圖中“機械應(yīng)力”為不考慮熱應(yīng)力，內(nèi)壓作用產(chǎn)生的機械應(yīng)力。從圖中可以看出，在氣瓶右封頭3 min 和5 min 充氫內(nèi)襯的應(yīng)力都將比機械應(yīng)力低，最大差別約70 MPa，同時3 min 完成充氫過程產(chǎn)生的熱機械應(yīng)力要比5 min 產(chǎn)生的熱應(yīng)力要高，但是兩者差別較小，氣瓶其他部分兩者幾乎一致。這是由于氣瓶底部溫度較高，將導致內(nèi)襯膨脹，但是外層復(fù)合材料將阻止內(nèi)襯膨脹，因此，內(nèi)襯部分產(chǎn)生的熱應(yīng)力為壓應(yīng)力，而復(fù)合材料層產(chǎn)生的熱應(yīng)力為拉應(yīng)力。最外層16°螺旋纏繞層纖維方向應(yīng)力如圖11 所示，從圖中可以看出在封頭處復(fù)合材料層應(yīng)力波動幅度較大，這是由于封頭處本身應(yīng)力變化幅度較大，同時復(fù)合材料層厚度和纖維方向都發(fā)生著變化，加劇了波動幅度，文獻[22]也觀察到類似的變化。兩種充氫條件導致的復(fù)合材料層熱機械應(yīng)力，在氣瓶右封頭處將比機械應(yīng)力高，而在氣瓶其他部位兩者相近，同時5 min 充氫的應(yīng)力要比3 min 充氫應(yīng)力略高。如前所述，由于右封頭處溫度較高，復(fù)合材料承受的熱應(yīng)力為拉應(yīng)力，因此熱機械耦應(yīng)力將提高。

圖 9 不同充氫速率下的內(nèi)襯應(yīng)力云圖Fig. 9 Stress contour of liner at the end of filling at different rates

圖 10 不同充氫速率下內(nèi)襯應(yīng)力沿著軸向分布Fig. 10 Stress distribution along the axis direction of the liner at different rates

圖 11 不同充氫速率下最外層16°螺旋纏繞層纖維應(yīng)力沿著軸向分布Fig. 11 The fiber stress versus the location in the outer 16° layer at different rates

2.2 充氫位置

由于充氫口的位置將影響氣體的流動和溫度分布，圖12(a)中充氫位置在氣瓶左封頭附近，以下簡稱a 情況，圖12(b)充氫位置在氣瓶中間，以下簡稱b 情況。從圖12(a)中可以看出，右封頭的溫度較高，而12(b)中左封頭溫度較高。內(nèi)襯溫度沿著軸向分布如圖13 所示。從圖12(c)中的流線圖可以看出渦流在左封頭，使得該處熱對流充分，溫度較為均勻，而底部沒有渦流，溫度較高。圖12(d)流線圖可以看出，氣瓶右封頭渦流較為明顯，該處傳熱效果較好，溫度較低。

兩種充氫位置下內(nèi)襯的應(yīng)力沿著軸向分布如圖14 所示，從圖中可以看出，a 情況下內(nèi)襯Mises應(yīng)力與內(nèi)壓作用下產(chǎn)生的機械應(yīng)力幾乎一致，但是在右封頭處偏離了機械應(yīng)力，最大偏差能到70 MPa；對于b 情況下，左封頭處熱機械耦合應(yīng)力小于機械應(yīng)力，最大應(yīng)力差能到150 MPa，其他部位兩者差別不大。原因在于氣瓶的溫度場分布，對于a 情況2.1 部分已做分析；對于b 情況，高溫在氣瓶左封頭處，并且b 情況下最高溫度比a 情況大，因此由于熱應(yīng)力導致的應(yīng)力下降要比b 情況要比a 情況大。

圖 12 不同充氫位置下復(fù)合材料氣瓶溫度云圖和流線圖Fig. 12 Temperature contour and streamline of COPV with different inlet locations

圖 13 不同充氫位置下溫度沿軸向分布Fig. 13 Temperature distribution in the liner along the axis with different inlet locations

圖 14 不同充氫位置內(nèi)襯應(yīng)力沿著軸向分布Fig. 14 Stress distribution along the axis with different inlet locations

最外層16°螺旋纏繞層纖維方向的應(yīng)力如圖15所示，對于a 情況，熱機械耦合應(yīng)力在氣瓶左側(cè)和機械應(yīng)力一致，但是當來到氣瓶右封頭處，熱機械耦合應(yīng)力要比機械應(yīng)力大；對于b 情況，熱機械耦合應(yīng)力在氣瓶左側(cè)比機械應(yīng)力大，但是慢慢的降到機械應(yīng)力。對于熱機械耦合應(yīng)力和機械應(yīng)力的差，b 情況下要比a 情況下大，原因在于兩種情況最高溫度存在差別。

圖 15 不同充氫位置下最外層16°螺旋纏繞層纖維方向應(yīng)力沿著軸向分布Fig. 15 The fiber stress versus the location in the outer 16°layer with different inlet locations

從以上分析可以看出，充氫口位置將影響溫度分布，進而影響熱機械耦合應(yīng)力分布。對于充氫口位于左封頭處，高溫將降低右側(cè)封頭內(nèi)襯應(yīng)力而提高外層復(fù)合材料壓力；對于充氣口位于氣瓶中間位置，高溫將降低左側(cè)封頭處內(nèi)襯應(yīng)力，提高該處的復(fù)合材料應(yīng)力。

2.3 幾何尺寸

以兩只內(nèi)徑分別為350 mm 和450 mm，但體積均為150 L 的復(fù)合材料氣瓶為研究對象，開展了復(fù)合材料氣瓶的幾何尺寸的影響研究。充氫時間為3 min，充氫口位于氣瓶左封頭處。氣瓶溫度云圖和流線如圖16 所示，沿著軸向的溫度場分布如圖17 所示。由于兩者幾何尺寸上的差別，橫坐標為軸向長度正則化坐標，定義為實際坐標值與氣瓶長度比值。

圖 16 不同內(nèi)徑復(fù)合材料氣瓶溫度云圖和流線圖Fig. 16 Temperature contour and streamline of COPV with different inner radius

圖 17 不同內(nèi)徑氣瓶溫度沿軸向分布Fig. 17 Temperature distribution in the liner along the axis with different inner radius

從對比中可以看出，內(nèi)膽外徑350 mm 的氣瓶最高溫度在右側(cè)，但內(nèi)膽外徑450 mm 氣瓶的最高溫度在中間。這可由圖中的流線給出解釋，在內(nèi)徑450 mm 氣瓶渦流大，幾乎占滿了整個瓶子，因此熱對流充分，溫度分布較為均勻。由于兩只瓶子幾何尺寸的差別，只研究了兩者熱應(yīng)力分布，內(nèi)襯應(yīng)力沿著軸向分布如圖18 所示。從圖中可以看出，內(nèi)徑450 mm 瓶子的應(yīng)力分布更均勻，在氣瓶右封頭處，內(nèi)徑350 mm 瓶子的熱應(yīng)力比450 mm氣瓶大，在氣瓶其他部位，內(nèi)膽外徑350 mm 氣瓶的熱應(yīng)力比450 mm 氣瓶小。

兩只氣瓶最外層16°螺旋纏繞層纖維方向應(yīng)力沿著軸向分布如圖19 所示。從圖中可以看出，從左側(cè)到右側(cè)，應(yīng)力從壓應(yīng)力逐漸變?yōu)槔瓚?yīng)力。內(nèi)膽外徑450 mm 氣瓶的應(yīng)力分布較為均勻，而內(nèi)徑350 mm 氣瓶應(yīng)力變化幅度較大。這是由于兩者溫度場分布造成，內(nèi)徑350 mm 氣瓶高溫存在于右側(cè)封頭處，而內(nèi)徑450 mm 氣瓶溫度分布更為均勻。

圖 18 不同內(nèi)徑氣瓶內(nèi)襯應(yīng)力沿軸向分布Fig. 18 Stress distribution along the axis with different inner radius

圖 19 不同內(nèi)徑氣瓶最外層16°螺旋纏繞層纖維方向應(yīng)力分布Fig. 19 The fiber stress versus the location in the outer 16°layer with different inner radius

從以上分析可以看出，與小內(nèi)徑氣瓶相比，大內(nèi)徑氣瓶熱傳導較為充分，溫度場分布較為均勻，因此熱應(yīng)力較小。

3 結(jié)論

本文提出了復(fù)合材料氣瓶快充過程的流-熱-固耦合分析方法。該方法中首先建立模擬氣瓶快充的CFD 模型；然后，在CFD 計算結(jié)果的基礎(chǔ)上，建立熱-機械耦合模型，并通過理論解和試驗驗證方法的準確性。通過該方法研究不同充氫速率、充氫位置和幾何尺寸對熱應(yīng)力的影響，為復(fù)合材料氣瓶的設(shè)計和充氫工藝提供技術(shù)支持。本文獲得以下結(jié)論：

(1) 充氫速率會影響復(fù)合材料氣瓶的溫度和應(yīng)力場分布，降低充氫速率可以為熱傳導提供更多時間，因此溫度場和應(yīng)力場將不同，但是區(qū)別不明顯。

(2) 不同充氫位置將導致不同的溫度場分布，當充氫口位于左側(cè)瓶口處，高溫將位于氣瓶右側(cè)封頭，充氫口位于氣瓶中間位置，高溫將位于氣瓶左側(cè)封頭處，高溫降低該處內(nèi)襯應(yīng)力，提高外層復(fù)合材料應(yīng)力。

(3) 大內(nèi)徑氣瓶溫度場較為均勻，因此熱應(yīng)力相對于小內(nèi)徑氣瓶較低。