為應對氣候變化，全球各國制定了一系列政策條款要求在2030年實現碳中和，尋求和發展清潔可再生能源，完成從化石能源向清潔能源的快速過渡。在可再生能源中，氫能作為一種幾乎零CO2和溫室氣體排放的能源[1~3]，引起了各國的重點關注。在氫能的眾多運輸方式中，管道運輸具有體量大、距離遠、能耗損失低等優點，可增大氫能資源的市場覆蓋度，在能源結構優化調整中具有很高的經濟性[4, 5]，但新建管道成本較高制約了氫能的長距離輸送環節發展[6]。目前，利用現有天然氣管道摻氫輸送可實現氫能長距離、大規模、網絡化運輸需求，國內外也已經開展大量天然氣和氫氣摻混技術研究[7~10]。然而，現有天然氣管道長期服役后內外均會出現腐蝕缺陷，在此區域會有應力集中，并且氫在此區域容易發生聚集。隨后應力集中和氫富集會誘發材料產生微裂紋，最終導致材料發生脆性斷裂，造成安全風險。

在役天然氣管道在運行過程中，不可避免地受到自然災害、外部壓力、環境影響或第三方破壞而產生腐蝕、凹陷、裂紋等各種宏觀結構缺陷[11]，同時金屬管材內部存在大量位錯、晶界、析出相、非金屬夾雜物等微觀組織缺陷[12]。目前，對位錯、晶界和非金屬夾雜物[10,13]等各種微觀組織缺陷的氫捕獲行為進行了廣泛研究，而關于氫原子在管道腐蝕、凹陷和裂紋等宏觀結構缺陷上的擴散和分布研究不足。郭詩雯等[14]針對雙腐蝕缺陷管道上的氫濃度分布規律研究表明，缺陷長度會顯著影響其周圍氫分布狀態。劉韋辰等[15]基于有限元的數值分析模型，研究了拉伸應變作用下X80管道腐蝕缺陷處的氫原子分布，研究表明隨著缺陷長度的減小和缺陷深度的增大，氫原子更集中在缺陷的中心處。由此不難看出，內部壓力、缺陷位置、缺陷長度和深度等因素會對缺陷處氫原子富集狀態產生顯著影響。然而，管材內局部富集的氫濃度無法通過現有技術手段測量，因此明確氫原子在腐蝕缺陷處擴散富集行為及濃度分布狀態，對現有老化管道輸氫適用性評估具有重要意義。

本文以含腐蝕缺陷的X52管道為研究對象，開展腐蝕缺陷處氫擴散和分布特征研究。通過有限元模擬，耦合應力場和氫擴散場，研究不同內壓下管道缺陷處的應力和氫濃度分布，分析缺陷長度和深度對最大氫濃度的影響。本文重點分析氫與腐蝕缺陷的相互作用，對老化管道的風險評估和氫管道的完整性管理具有重要意義。

氫在金屬中擴散是物質交換過程，它可以由不同因素驅動，直至達到熱力學平衡[13]。其驅動力包括濃度、應力、溫度和其他化學勢梯度，它們相互關聯并相互作用，導致管道中的氫原子發生擴散[16]，最終氫原子均勻分布或聚集在晶格和氫陷阱內。此時，氫通量J可以表示為[17]：

(1)

式中，D是氫擴散系數，4.32 × 10-11 m2·s-1 [18]；R是氣體常數，8.41 J·mol-1·K-1；T是絕對溫度，273 K；C是氫濃度，mol·m-3；μ是化學勢。

當不考慮應力狀態時，化學勢可以寫成[19]：

(2)

式中，μ0是參考溫度和壓力下的化學勢。

然而，如果系統受到外力，化學勢μ應改寫為[20]：

(3)

式中，μσ是μ的與應力相關部分，描述為[21]：

(4)

式中，σh是靜水應力，等于3個主應力的平均值，σh =$\frac{\mathrm{1}}{\mathrm{3}}{\sum }_{i=\mathrm{1}}^{\mathrm{3}}{\sigma }_{ii}$；Vh是氫的偏摩爾體積，2 × 10-6 m3·mol-1 [22]。將 式(3)和(4)代入 式(1)，得到下式：

(5)

考慮一個任意體積V的材料，其邊界是表面S，質量守恒要求V內部總氫的變化率等于通過表面S的通量[17]：

(6)

式中，∂/∂t是對時間的偏導數；n是向外的單位法向量。將 式(6)中的J替換為 式(5)中所給出的可以得到：

(7)

通過應用散度定理， 式(7)可以表示為：

(8)

該方程是Fick第二定律的修正形式。在有限元軟件COMSOL Multiphysics中，稀物質傳遞的質量守恒方程如下[23]：

(9)

式中，Ra是物質的反應速率表達式。因此，在COMSOL中將Ra定義為-$\nabla$·[(DCVh/RT)$\nabla$σh]模擬管道中的氫擴散。

在COMSOL Multiphysics 5.6軟件中，建立了一個三維(3-D)有限元模型，用于模擬和分析內壓為6、8和10 MPa下管道內外腐蝕缺陷處的氫擴散和分布，如圖1所示。圖1a中，模擬管段的幾何參數為管長L = 2000 mm，外徑D = 508 mm，壁厚t = 9.525 mm[24]。腐蝕缺陷寬度為10 mm，缺陷長度分別為10、12、14和16 mm，缺陷深度分別為壁厚的20%、30%、40%和50%，即1.10、2.86、3.81和4.76 mm，以確定圖1b中腐蝕缺陷幾何尺寸的影響。為了節省計算時間，在模擬中使用了四分之一的模型以及其中的內部和外部腐蝕缺陷，如圖1c和d所示。

假設模型中初始氫濃度為10 mol/m3 (≈ 1.29 × 10-3 mg/g[25])，即管道內表面的邊界氫濃度。在模擬過程中，沒有氫通量通過模型邊界。遠離缺陷的管段一端為固定約束，其余管段邊界沿壁厚方向為對稱約束。值得注意的是，當建模管道的長度超過4.9(Dt/2)1/2時，其中D是管道的直徑(mm)，t是管道壁厚(mm)，邊界條件不會影響建模結果[26]。在本研究中，4.9(Dt/2)1/2 = 241.02 mm，2000 mm的管道長度足以消除邊界條件對有限元模型的影響。

為了平衡計算精度和計算時間，如圖2所示，該模型由兩個正方形分成3個區域。其中，區域1是邊長l1 = 30 mm的第一個正方形，包含腐蝕缺陷；區域2是邊長為l2 = 100 mm的第二個正方形，但不包含區域1；余下部分為區域3。區域2和3采用粗六面體網格，壁厚方向均勻分為4層。區域1用更精細的四面體網格劃分，腐蝕缺陷附近最小網格尺寸為1~5 mm。

圖3為網格敏感性分析結果。可以看出，當最小單元尺寸小于3 mm時，最大應力和最大氫濃度基本不變，因此最小單元尺寸選擇為1 mm，最大單元增長速率為1.1，完整網格由53949個單元組成。采用MUMPS (multi-frontal massively parallel sparse)求解器進行求解。模擬的擴散過程持續5 × 106 s。

模擬中使用的Young′s模量E = 207 GPa，Poisson比ε = 0.3。描述鋼彈塑性行為的硬化函數σhard如下[27]：

(10)

式中，σexp為根據X52鋼工程應力-應變曲線實測的應力函數；εp是塑性應變；σe是有效應力(MPa)；σe/E是彈性應變；σy是屈服強度(MPa)；X52鋼的屈服強度和極限抗拉強度分別為359和455 MPa。

圖4展示了在不同內壓下，長度為16 mm，深度為4.76 mm的內部和外部腐蝕缺陷處von Mises應力分布。可以看出，內壓的存在使腐蝕缺陷處von Mises應力分布不均勻，發生局部應力集中，該區域應力遠高于管段其他部位。6 MPa變化到8 MPa壓力時，von Mises應力較大區域變大、應力變大。當內部壓力為10 MPa時，沿缺陷邊緣von Mises應力更大，并且僅存少部分應力較小區域。內、外腐蝕的von Mises應力分布基本相同，應力值隨內壓的增加而增大。

為了分析應力和氫濃度沿著腐蝕缺陷邊緣的線性分布，對缺陷邊緣進行標記，如圖1e和f所示，其中內腐蝕的左邊緣和右邊緣分別為1-2-3-4和1-2-5-6，外腐蝕的左邊緣和右邊緣分別為1-2-3-4-5和1-2-6-7-8。

圖5a和b為沿內部腐蝕缺陷左邊緣和右邊緣von Mises應力分布。左邊緣von Mises應力從最大值減小到最小值，并在小幅增加后保持穩定。右邊緣的應力也從最大值減小到最小值，后上升到一個局部峰值，并在略微減小后趨于穩定。不同內壓下的von Mises應力大小順序為10 MPa > 8 MPa > 6 MPa。圖5c和d為沿外部腐蝕缺陷的左邊緣和右邊緣的von Mises應力分布。左邊緣的von Mises應力從最小值到最大值增加并略有下降，隨后迅速下降，最后達到穩定狀態，此時6 MPa的應力等于開始時的最小值。右邊緣的應力首先迅速上升，然后迅速下降到最小值，并再次迅速增加到一個小的局部峰值。在此峰值之后，它會再降低，隨后緩慢地上升到穩定值。總的來看，內壓越大，局部缺陷處的von Mises應力越大。

圖6為不同內壓下，內部和外部腐蝕缺陷(長度為16 mm，深度為4.76 mm)處的靜水應力分布。圖中的彩色圖例表示靜水應力(MPa)的大小。在內壓作用下，腐蝕缺陷處會產生靜水應力集中。在6 MPa下，缺陷處的較大靜水應力在應力條上為黃色，較小的為淺藍色。在8 MPa的內壓下，較大的應力位于紅色區域，面積擴大，但藍色區域略有縮小。當內壓達到10 MPa時，出現暗紅色和深藍色區域，并且紅色區域略微偏離缺陷邊緣。總的來看，內、外腐蝕缺陷處的靜水應力分布相似，并隨內壓的增加而增大。

圖7a和b為沿內部腐蝕缺陷的左邊緣和右邊緣的流體靜壓力分布。在10 MPa下，沿左右邊緣的靜水應力先略有上升，然后迅速下降到最小值，而在6和8 MPa下，均從初始最大值下降至最小值。3者共同增加后，左邊界應力趨于穩定，右邊界應力出現一個小峰，然后趨于穩定。圖7c和d為沿外部腐蝕缺陷的左邊緣和右邊緣的流體靜壓力分布。左邊緣的靜水應力先上升到最大值，隨后下降到最小值，之后經過一個波動后趨于穩定。右邊緣的靜水應力先迅速增大，隨后迅速減小到0附近，再增大至穩定。并且10 MPa的初始值最小，為負值。

圖8為10 MPa下外部腐蝕缺陷(長度16 mm，深度4.76 mm)處氫擴散過程，其中圖例條為氫濃度(mol/m3)的大小。從腐蝕開始到1 × 106 s，氫的擴散是均勻的，只在腐蝕缺陷處氫的擴散速度較快。1.5 × 106 s后氫在腐蝕缺陷處積累，濃度條上的暗紅色表示氫濃度變大，直到2.5 × 106 s擴散達到穩定狀態，氫分布保持不變。并且由于內部壓力，缺陷處最終氫濃度超過初始邊界值10 mol/m3。

圖9顯示了在5 × 106 s時不同內壓下內部和外部腐蝕缺陷(長度16 mm，深度4.76 mm)處的氫分布，其中圖例條為氫濃度(mol/m3)的大小。在內壓作用下，氫在腐蝕缺陷處局部積累，積累氫濃度隨著內壓的增加而增加。在5 × 106 s時，內、外腐蝕的最大濃度均高于初始邊界條件10 mol/m3。但值得注意的是，外腐蝕缺陷處的氫濃度明顯高于內腐蝕缺陷處。

圖10a和b為在5 × 106 s時氫濃度沿著內部腐蝕缺陷的左邊緣和右邊緣的分布。沿左右邊緣的氫氣濃度均從10 mol/m3的初始濃度短暫上升到最大值，然后下降到最小值。缺陷附近的氫濃度從最小值再次增加到10 mol/m3。并在遠離缺陷處保持在10 mol/m3。其中，在10 MPa下的最大值高于6和8 MPa，最小值低于6和8 MPa，并且10 MPa內壓條件下氫在缺陷處分布最不均勻。圖10c和d為在5 × 106 s時氫濃度沿著外部腐蝕缺陷的左邊緣和右邊緣的分布。氫濃度均從10 mol/m3的初始濃度上升到最大值，然后下降，并且右邊緣濃度下降到10 mol/m3以下，然后上升到10 mol/m3并且在遠離缺陷處穩定。隨后，在左邊緣的濃度在未達到10 mol/m3時略有反彈，然后再次降低，并最終在遠離缺陷處穩定在10 mol/m3。總的來看，內部壓力越大，最大氫濃度越高，但最小氫濃度的差異不大。

通過COMSOL中的體積最大值和表面最大值計算方法，在不同內壓下，5 × 106 s時最終氫濃度的最大值(Cmax)和位置坐標列于表1中，并在圖11中進行了展示說明。由此進一步表明，內腐蝕缺陷處的Cmax比外腐蝕處的高，表明在相同的內壓下，外腐蝕缺陷比內腐蝕缺陷更容易促進氫的積累，并且Cmax位置會隨內壓而變化。

圖12為在5 × 106 s時10 MPa下具有不同長度(深度為4.76 mm)的內部和外部腐蝕缺陷處的氫分布。顏色圖例為氫濃度(mol/m3)的大小。從圖12a~d可以看出，當長度為10 mm時，氫在內腐蝕缺陷處有明顯聚集，濃度條上以暗紅色表示，且缺陷處暗紅色濃度僅隨長度增加而略有增加。當缺陷長度為12 mm時，缺陷處的氫濃度在濃度條上為淺紅色，表示為較高氫濃度。如圖12e~h所示，不同長度的外腐蝕缺陷處的氫積累也很明顯，當缺陷長度從10 mm增加到16 mm時，濃度條上的氫濃度由紅色變為暗紅色，表示為濃度增加。由此可得，隨著缺陷長度的增加，氫的聚集區域擴大。不同缺陷長度下的最大氫濃度和位置坐標列于表2中。

由圖13可得，在5 × 106 s時在10 MPa下具有不同深度(長度為16 mm)的內部和外部腐蝕缺陷處的氫分布。顏色圖例為氫濃度(mol/m3)的大小。如圖13a~d所示，在1.91 mm深度處，腐蝕缺陷處的氫濃度為深紅色略高于其他部位。當深度增加到2.86 mm及以上時，缺陷處的氫濃度在濃度棒上呈暗紅色表示為濃度最高，并且氫積累區域變得大于1.91 mm處。從圖13e~h中可以看出，當缺陷深度為1.91 mm時，氫在缺陷處不會明顯地積聚，如濃度條上的黃色所示。由此可得，隨著缺陷深度的逐漸增加，外腐蝕缺陷處的局部氫濃度逐漸升高，依次出現橙、紅、暗紅3種顏色。表3中列出了不同缺陷深度處的最大氫濃度和位置坐標。

在本研究中，鋼管只承受由內壓引起的應力。通常情況下，應力水平會隨著內壓的增加而增加。但由于腐蝕缺陷的存在導致缺陷處應力集中，其原因為腐蝕使缺陷處的管壁厚度比未腐蝕處薄，導致該處在內壓作用下應力較大。從圖4和6可以看出，缺陷處的von Mises應力和靜水應力的集中分布是不一樣的。圖4中的von Mises應力集中沿著一個缺陷邊緣向內或外腐蝕缺陷的曲面擴散，而圖6中的靜水應力沿缺陷的一側向壁厚方向分布。當內壓為6 MPa時，這種現象并不明顯，但在內壓為8和10 MPa時卻很容易觀察到，特別是在10 MPa時的暗紅色區域。總的來看，在不同的內壓作用下，遠離腐蝕缺陷處的管道上應力分布相對均勻。

從圖5中可以看出，在3種內壓作用下，沿缺陷邊緣的von Mises應力最大值分別為260、290和325 MPa，均未超過X52鋼的屈服應力，因此缺陷處于彈性狀態，沒有局部塑性變形。在弧長超過0.04 mm的缺陷邊緣上，均勻分布的von Mises應力分別為140、170和220 MPa左右。此外，圖7示出了缺陷邊緣上的最大靜水應力為160、210和260 MPa，遠離缺陷的穩定靜水應力約為70、90和120 MPa。值得注意的是，缺陷處的最小靜水應力可以小于0，例如在10 MPa下可以達到-70~-80 MPa。根據 式(3)，驅動氫擴散的化學勢與應力有關，并且 式(4)也進一步表明靜水應力是影響氫擴散的直接因素。

圖8以存在外腐蝕缺陷的管道為例，模擬了10 MPa內壓下管道中氫的擴散過程。在這種情況下，濃度梯度和靜水應力梯度一起作為氫擴散的驅動力。在未被腐蝕的部分，氫氣沿壁厚方向從管內壁向外壁均勻擴散。在5 × 106 s后，氫擴散到腐蝕缺陷附近，由于缺陷的存在，氫在腐蝕缺陷處的擴散速度要快于其他地方。在1.5 × 106 s時，缺陷處出現氫積累，且積累面積和氫濃度隨時間而增加。當擴散進行到2.5 × 106 s時，氫的分布達到穩定狀態，并且缺陷處的最大氫濃度高于管道上為10 mol/m3的均勻氫濃度。

從圖9可以看出，當缺陷長度和深度不變時，無論是內腐蝕還是外腐蝕，內壓越大，缺陷處的氫積累面積和濃度就會越大。并且氫的分布與圖6中的靜水應力分布相似。通過比較圖10中的氫濃度沿著缺陷邊緣的線性分布和圖7中的靜水應力更容易看出這一點。從表1中可以看出，隨著內壓從6 MPa增加到10 MPa，內腐蝕的Cmax從10.2 mol/m3增加到10.6 mol/m3，外腐蝕的Cmax從10.94 mol/m3增加到11.53 mol/m3。由此可得，內壓使局部氫濃度過大，容易達到氫致開裂的門檻值。在這種情況下，更高的氫濃度意味著管道上氫致開裂的風險更高。

此外，Cmax的x坐標均為0，表明它們位于y-z平面中，也就是圖9中模型的左界面。值得注意的是，在圖10中，盡管這兩個濃度看起來很接近，最大值卻不一定在缺陷邊緣。此外，從表1的數據來看，不同內壓下Cmax的y、z坐標發生了變化，說明內壓可能對最大氫濃度產生影響，進而影響氫致開裂的起始位置，但此現象在本研究中并不顯著，因此無法定性或定量地解釋這種影響。

在圖9中，相同的內壓下內腐蝕缺陷處積累的氫濃度明顯低于外腐蝕。在濃度條上，前者用橙色和紅色表示，后者用暗紅色表示，盡管它們在圖6中的靜水應力顏色是相同的。比較圖10的a，b和c，d，外腐蝕的最大邊緣濃度也明顯高于內腐蝕。這表明，管道腐蝕部位影響缺陷處的氫積累程度，而且內壓對外腐蝕處氫積累的影響更大。同樣，表1中的Cmax值也證明了這一點。外腐蝕Cmax超過邊界條件10 mol/m3以上，并且當內壓從6 MPa變化到10 MPa時，內腐蝕Cmax增加0.4 mol/m3，外腐蝕Cmax增加0.59 mol/m3，變化幅度更大。以上結果表明，在相同條件下，外腐蝕管道相比內腐蝕管道更容易發生氫致開裂。

此外，內、外腐蝕缺陷的Cmax位置變化也不同。從圖11中可以看出，內腐蝕的Cmax位置隨內壓的變化比外腐蝕的大，但都遠離缺陷弧面。因此，相同條件下的內、外腐蝕缺陷也可能具有不同的氫致開裂起始部位。

圖12a~d中不同長度的內腐蝕缺陷處氫積累的顏色變化不明顯，均為暗紅色，但積累區域隨著長度的增加而擴大，因為缺陷長度增加擴大了應力集中區，進而增大氫積累區。如圖12e~h所示，外腐蝕缺陷處的濃度顏色隨缺陷長度的增加而加深，在14和16 mm深度處出現暗紅色，氫積累面積增大。在表2中，對于內部腐蝕，Cmax值都在10.5 mol/m3左右，并且只有非常輕微的變化。只有10 mm深度處的Cmax位置接近缺損，其他三點幾乎重合。而且在各種長度下外腐蝕的Cmax位置沒有變化，濃度變化也很小，僅為0.33 mol/m3。由此可見，缺陷長度影響氫的積累范圍，但對最大氫濃度影響不大。

圖13中，缺陷處的氫濃度隨著缺陷的深度而明顯地從淺紅色變為深紅色，并且積累區域也逐漸擴大。隨著缺陷深度的增加，缺陷處的壁厚變薄，導致更大程度的應力集中和更大的累積面積。然而，從表3中的數據來看，當缺陷深度從1.91 mm增加到4.76 mm時，內部腐蝕Cmax增加了0.27 mol/m3，而外部增加了0.91 mol/m3，遠遠大于內部腐蝕。因此，缺陷深度對外腐蝕缺陷處氫濃度分布的影響較大。并且深度增加導致最大值發生在更靠近管壁的位置，即內部腐蝕的Cmax移動到外壁，外部腐蝕的Cmax移動到內壁。

(1) 缺陷的存在加速了缺陷附近處氫的擴散，并且由于內壓的存在導致腐蝕缺陷處應力集中、氫分布不均勻。在相同內壓下，von Mises應力和靜水應力分布并不相同。但是內壓的增加使兩種應力均升高，并使缺陷處的應力集中程度增大，從而導致氫在腐蝕缺陷處積聚，并使局部氫濃度增大。而且不同內壓下，氫氣濃度最大值的位置會發生移動，但移動變化不大。

(2) 在相同內壓下，外腐蝕缺陷處的氫濃度最大值大于內腐蝕缺陷處，且外腐蝕缺陷處的氫濃度隨內壓的變化較大。此外，內部和外部腐蝕的最大氫濃度的位置是不同的，但它們通常會遠離缺陷處。

(3) 缺陷長度的增加使外腐蝕氫濃度最大值略有增大，而對內腐蝕氫濃度影響不大。在不同的長度上，外腐蝕最大濃度位置保持不變，而應力集中區和氫積累區均隨缺陷長度的增加而擴大。隨著缺陷深度的增加，氫濃度隨之增大，內腐蝕和外腐蝕的氫積累面積也增大，同時最大氫濃度產生于更靠近管壁的位置。