金屬表面點蝕的穩態過程數值研究

崔靜，楊帆，楊霆浩，楊廣峰

中國民航大學機場學院，天津300300

0 引 言

現代艦船裝備結構的系統化和復雜化程度越來越高，在艦船上使用的新材料也越來越廣泛，因此艦船的腐蝕防護工作面臨著巨大挑戰。金屬腐蝕會導致艦船裝備的損壞，造成經濟損失以及腐蝕污染。預測腐蝕的損傷過程，對于保證艦船裝備的安全和正常運行起著重要的作用。艦船所處的海洋環境復雜，運用實際測量和經驗估算等方法對艦船材料表面腐蝕進行檢測時，評估標準不統一，不能實時觀測腐蝕發展過程，而采用計算機技術進行數值模擬，在一定程度上可彌補傳統方法的不足。

腐蝕過程是一個包含了多組分之間的質量傳輸、化學反應或電化學反應、界面動力學等多重因素影響的復雜過程，而格子Boltzmann方法（Lattice Boltzmann Method，LBM）作為一種介觀方法（Mesoscopic technique），近年來發展迅速，成功地應用于多種微觀相互作用的流體系統問題中。LBM既能兼顧實際的計算成本，又能保證模型或方法的合理性。國內外運用LBM進行腐蝕模擬的研究較少，更多學者采用的是元胞自動機模型（Cellular Automata Model，CAM）進行腐蝕模擬，采用CAM進行腐蝕模擬對于多組分之間的質量傳輸和腐蝕化學反應只進行了簡單的元胞轉化，并沒有理論支撐。He等［1］首次將化學反應與LBM結合，解決了流體—壁面表面反應與擴散問題的耦合。Kang等［2-3］對溶解沉淀問題中單相流體和多相流體的傳輸及流動進行了研究，將LBM與結點體積法（Volumn of Pixel，VOP）結合來模擬物質的傳輸過程。Chen等［4］為了解決多相流動中的相變和溶解沉淀問題，將LBM和VOP方法在孔隙尺度下相結合，解決了多相流中的相分離、化學反應、質量傳輸和溶解沉淀等問題。張婷等［5］使用MRT-LBM來模擬腐蝕溶液的流動，并分析了不同參數對腐蝕的影響。Liu等［6］使用D3Q19模型模擬了巖石的腐蝕情況。Zhang等［7］在腐蝕模型中加入了電場并與流場成功耦合，更好地模擬了電化學腐蝕。但是，現有的格子Boltzmann（LB）腐蝕模型中，并未對具體的金屬材料進行模擬，更多是使用化學腐蝕的方法進行研究，僅考慮了陽極固體腐蝕，很少涉及陰極反應的發生。

本文將針對現有LB腐蝕模型的不足，采用陰陽極共同作用的方式研究電化學腐蝕，將使用改進后的模型對金屬表面點蝕現象的穩態過程進行模擬，得到點蝕坑的形貌變化特征以及各組分濃度的變化情況。通過數值模擬，分析腐蝕化學反應速率、腐蝕溶液擴散系數、腐蝕產物擴散系數對腐蝕程度的影響。

1 格子Boltzmann腐蝕模型

腐蝕過程包含多相多組分的流動與擴散、電化學反應、金屬的溶解以及腐蝕產物沉淀。根據以上問題建立多個子模型，并將子模型耦合形成LB腐蝕模型。其中，子模型包括電化學反應模型、多組分模型、質量傳輸模型、固體邊界遷移模型和結點鈍化模型。

1.1 電化學反應模型

模擬浸沒于液體腐蝕環境的金屬（以不銹鋼為例）表面的點蝕現象，其腐蝕過程如圖1所示。

圖1 物理模型示意圖Fig.1 Schematic illustration of corrosion model

為了便于數值建模，本文將復雜的電化學反應過程進行如下簡化。

陽極：

陰極：

沉淀反應：

沉淀氧化：

亞鐵離子氧化：

鐵離子水解：

自催化反應：

假設腐蝕環境已達到蝕孔可發生穩態點蝕時（腐蝕深度達到一定值，腐蝕產物溶度達到其飽和濃度的 60%～80%［8］），如圖1（a）所示。蝕孔內金屬基體處于活化溶解狀態，作為陽極，蝕孔外鈍化膜作為陰極，孔內外形成了“大陰極，小陽極”的電化學腐蝕體系，陰陽極發生式（1）和式（2）的反應，如圖1（b）所示；由于孔外發生式（2）反應，孔口處OH-逐漸增多，孔內的金屬離子向孔外擴散，發生式（3）和式（4）的反應，在孔口附近產生沉淀；隨著孔內金屬不斷腐蝕溶解，金屬離子逐漸增多，繼而發生式（5）的水解，孔內H+增多，呈酸化環境，進一步加快式（6）的腐蝕，如圖1（c）所示。

1.2 多組分模型

使用Shan-Chen模型［2-4］來處理多組分介質間的流動，其演化方程為

式中：ωα為權系數，ω0=4/9，ω1-4=1/9，ω5-8=1/36；宏觀量密度ρ和速度u可由密度分布函數求得：

在Shan-Chen模型中，為了降低誤差和保證計算結果的穩定性，使用了分子作用力對速度進行修正，分子間的作用力計算如下：

式中：Fσσ為液體分子間的作用力；Gσσ為作用強度；ψσ為有效密度函數［9-12］；為權系數，其中

式中：Fσw為液體分子與固相分子之間的作用力；wσw為液體與固體壁面的作用強度，可通過接觸角的變化進行調整；s為固體壁面密度［13-15］，與通過分布函數算出的密度不是同一概念。

將上述計算得到的分子間作用力代入下式進行速度修正：

式中：Fσ為該結點所受到的合力；u′為式（11）所求的速度；式（14）計算的u需代入式（9）進行平衡分布函數的計算。

1.3 質量傳輸模型

采用質量傳輸模型來解決溶液組分之間的相互傳輸問題，本文涉及到的液體組分為腐蝕溶液組分和腐蝕產物組分。LB質量傳輸模型的演化方程為

式中：Ck為第k組分的宏觀量濃度；u由式（11）計算得到；系數Jk，α為

其中，Jk，0的取值范圍在0～1之間。

濃度場中的擴散系數Dk為

式中，Jk，0取值為0.9。

1.4 固體邊界遷移模型

由于發生腐蝕溶解反應，被腐蝕的金屬基體表面結點的體積會隨之減小，并且當新生成的腐蝕產物濃度達到其飽和濃度時，在固體表面沉淀結點的體積會隨之增大，因此，固體邊界在腐蝕與沉淀的影響下不斷更新。為了能更好地實時追蹤新的固體邊界，采用VOP方法予以解決。

根據式（1）反應，金屬基體表面結點的體積變化［19］如下：

根據式（2）反應，新生成的沉淀結點在固體表面的體積變化如下：

然后我就發現，身邊很多女孩跟你一樣，也跟當初的我一樣，總是對自己的某一個缺點耿耿于懷。見人就說：我太胖了，我腿太粗了，我眼睛太小了，我牙齒太丑了……

金屬與腐蝕溶液在邊界處發生式（1）反應的腐蝕，腐蝕邊界處各組分濃度變化［7］為：

式中：CO2(aq)為腐蝕邊界上O2的消耗濃度；為腐蝕邊界上Fe3+的生成濃度；n為反應固體表面的外法線方向；為O2擴散系數；為Fe3+擴散系數。

1.5 結點鈍化模型

模擬過程中的腐蝕邊界，每個結點都有一定的概率被鈍化，導致該處的腐蝕速率變慢，使腐蝕坑的形狀呈現不對稱性，這也符合真實的腐蝕情況。在VOP方法中，腐蝕邊界處結點發生體積的變化，如果腐蝕邊界某一結點發生鈍化，該結點的體積減小速率會變慢。結點發生鈍化的概率稱為鈍化概率P，鈍化概率系數ε為

式中，R為0～1之間的任意隨機實數。因此，需對式（22）進行修正：

在VOP方法中引入鈍化概率函數，可以避免腐蝕界面上的陽極金屬受到的腐蝕作用相同，以使腐蝕坑的形貌隨機變化且不對稱，更加符合真實的腐蝕情況。

2 數值模擬

2.1 物理模型

基于上文構建的LB腐蝕模型，對金屬表面點蝕穩態過程進行模擬。計算區域選取100×160的網格。整個計算域內，初始的初生蝕孔假設為一部分圓，其圓心位于（71，80），半徑為10，鈍化膜的厚度為5，如圖1（a）所示。

2.2 邊界條件

對于液體四周的邊界采用了非平衡外推格式進行處理，而溶液與金屬之間、溶液與鈍化膜以及溶液與腐蝕產物沉淀之間的流動碰撞采用反彈邊界進行處理。

2.3 模擬工況

本文模擬金屬表面點蝕穩態過程的形貌的演化過程，參數如表1所示，表中v為運動粘度。采用變參數的數值模擬得到腐蝕程度的影響因素，各工況參數如表2所示［3，14］。所有擴散系數均由各濃度場的松弛時間決定，如式（14）所示。為了滿足運行過程中的數值穩定性，選取數值需大于0.5。

表1 模擬中使用的參數Table1 Simulation parameters

表2 模擬工況Table 2 Simulated condition

上述工況中所有參數均為無量綱參數，均通過表1中的轉化關系得到。由于化學反應速率受溫度、反應接觸面積、催化劑等影響，且腐蝕介質擴散系數受濃度、溫度等影響，這兩種參數可選取值范圍較廣［4，15］，因此選取參數時，僅在表1中的模擬數值附近選取。

3 數值模擬結果與分析

3.1 腐蝕演化過程分析

通過模擬得到不銹鋼材料表面點蝕穩態過程的形貌變化以及各組分濃度變化情況。模擬工況為k1=8.14×10-5，圖2所示為不銹鋼材料點蝕坑形貌隨時間變化的情況。從圖中可以看到，腐蝕破損區域逐漸變大，并且逐漸向金屬基體縱深發展。當蝕孔出現Fe（OH）3沉淀，隨著時間的推移沉淀增多，逐漸覆蓋了蝕孔（圖2（c））；除了初生的蝕孔發生腐蝕，在蝕孔底部出現次生蝕孔繼續向金屬基體深處腐蝕（圖2（d）～圖2（f））。

圖2 腐蝕坑形狀變化Fig.2 Shape evolution of corrosion pit

圖3 腐蝕溶液濃度變化Fig.3 Change in the concentration of corrosion solution

圖3所示為含氧腐蝕溶液的濃度隨時間變化的情況。隨著時間的延續，腐蝕的區域變大，并向金屬基體縱深發展，腐蝕溶液向蝕坑內滲入，腐蝕溶液逐漸被消耗，在蝕孔處出現沉淀，阻礙了腐蝕溶液向蝕孔內擴散，導致蝕孔內腐蝕溶液的濃度逐漸降低。

圖4所示為蝕孔內陽極金屬發生腐蝕溶解時，腐蝕產物Fe3+的濃度變化情況。隨著時間的延續，腐蝕產物在腐蝕邊界不斷生成，并且向蝕坑外以及液體環境擴散；當Fe3+擴散至蝕孔與孔外大量的OH-結合產生Fe（OH）3，Fe（OH）3達到其飽和濃度繼而沉淀。蝕孔產生沉淀，阻礙了Fe3+向孔外擴散，孔內的Fe3+升高，導致Fe3+進一步發生水解反應。

圖4 腐蝕產物Fe3+濃度變化Fig.4 Changes in the concentration of corrosion products Fe3+

圖5所示為蝕孔外鈍化膜作為陰極及其表面產生的OH-濃度變化及OH-擴散的情況。隨著時間的延續，蝕孔內陽極金屬發生失電子的腐蝕溶解反應，電子通過金屬基體轉移至電位較高的鈍化膜，O2在鈍化膜表面得到電子，發生陰極反應，在表面產生大量的OH-，并向液體環境擴散；OH-擴散至蝕孔處與Fe3+生成Fe（OH）3發生沉淀反應，孔口沉淀的增多阻礙了OH-向蝕孔內的擴散，導致孔內OH-濃度降低。

圖5 OH-濃度變化Fig.5 Changes in the concentration of OH-

圖6所示為H+濃度隨時間的變化情況。由于蝕孔處的沉淀不斷增多，阻礙了Fe3+向孔外擴散以及OH-向孔內擴散，導致孔內大量Fe3+聚集，繼而發生水解反應，產生大量的H+，初生蝕孔底部H+濃度變大，進一步腐蝕金屬基體，金屬加速溶解，產生次生孔，推動點蝕的穩定發展，蝕孔逐步向深處生長。

圖6 H+濃度變化Fig.6 Changes in the concentration of H+

圖7所示為整個腐蝕過程中被腐蝕的金屬體積隨時間的變化關系。當t=0～20 000時，金屬腐蝕僅發生初生蝕孔的吸氧腐蝕，初生蝕孔被腐蝕的體積V逐漸變大，蝕孔的深度h逐漸增大；當t＞20 000時，由于孔內閉塞及酸化的腐蝕環境，導致金屬自催化腐蝕反應的發生，腐蝕速率大于先前單獨發生的吸氧腐蝕的速率，并且V隨t的增加而增大，h也逐漸增大。

圖7 腐蝕程度隨時間的變化Fig.7 Variation of the degree of corrosion with respect to time

3.2 影響因素分析

3.2.1 腐蝕化學反應速率

模擬過程中，通過改變式（1）反應的化學反應速率，得到某一時刻（t=50 000）不同腐蝕化學反應速率下點蝕坑的形貌（圖8）。從圖中可以看出，在相同腐蝕時間下，隨著腐蝕化學反應速率的降低，腐蝕程度減小。

通過對圖8中V和h進行統計測量，得到t=50 000時V和h與k1的關系曲線（圖9）。從圖9可以看出，隨著k1增加，V越大，h就越大。由于式（1）中化學反應速率的增加，導致了吸氧腐蝕加快，腐蝕產物增多，蝕孔內Fe3+濃度增大，加速了Fe3+水解過程，產生更多的H+，使孔內陽極金屬進一步加速腐蝕，因此呈現出以上結果，并且該結果與真實情況相符。

圖8 不同腐蝕化學反應速率下的點蝕坑形貌Fig.8 The corrosion pit appearance of different chemical reaction rate

圖9 腐蝕程度隨腐蝕反應速率的變化Fig.9 The degree of corrosion varies with corrosion reaction rate

綜合分析圖8與圖9可知，由于k1增加，單位時間內V變大，h增加。

3.2.2 腐蝕溶液擴散系數

通過改變模擬過程中濃度場含氧溶液的擴散系數，研究腐蝕反應中介質間流動與擴散對腐蝕的影響。圖10給出了某一時刻（t=50 000），不同腐蝕溶液擴散系數下腐蝕坑的形貌。從圖中可以看出，相同的腐蝕時間下，隨著的增大，腐蝕程度加重。

圖10 不同腐蝕溶液擴散系數下的點蝕坑形貌Fig.10 The corrosion pit appearance of different solution diffusion coefficient

通過對圖10中V和h進行統計測量，得到t=50 000時V和h與DO2的關系曲線（圖11）。從圖11可以看出，隨著的增大，V越大，h越大。由于的增加，導致腐蝕溶液向蝕孔內擴散速率加快，金屬腐蝕邊界的濃度相對較高，加快了腐蝕反應，因此，呈現出以上模擬結果。

圖11 腐蝕程度隨O2擴散系數的變化Fig.11 Variation of the degree of corrosion with respect to diffusion coefficient of O2

3.2.3 腐蝕產物擴散系數

在模擬中改變濃度場的腐蝕產物Fe3+的擴散系數，研究腐蝕產物擴散系數對腐蝕程度的影響。圖12所示為某一時刻（t=50 000）不同腐蝕產物Fe3+擴散系數下腐蝕孔的形貌特征。從圖中可以看出，在相同的腐蝕時間下，隨著的減小，腐蝕程度加重。

圖12 不同腐蝕產物擴散系數下的腐蝕坑形貌Fig.12 The corrosion pit appearance of different corrosion products diffusion coefficient

通過對圖12中V和h進行統計測量，得到t=50 000時V和h與的關系曲線（圖13）。從圖13可以看出，隨著的降低，V越大，h越大。由于的降低，導致蝕孔內Fe3+的大量聚集，加速了Fe3+的水解過程，產生了更多的氫離子，使孔內的酸度增加，金屬基體的腐蝕溶解速率增加。因此，呈現出以上模擬結果。

圖13 腐蝕程度隨Fe3+擴散系數的變化Fig.13 Variation of the degree of corrosion varies with respect to diffusion coefficient of Fe3+

4 結 論

本文基于改進后的LB腐蝕模型，將鈍化概率函數引入VOP方法中，對腐蝕化學反應速率進行修正，建立了LB腐蝕模型。使用該模型，模擬了浸沒于液體腐蝕環境中的金屬表面單坑點蝕穩態過程，通過模擬得到以下結論：

1）LBM算法可以明確得到溶液中不同組分的流動與擴散，由此得到其濃度分布，通過流場和濃度場之間的耦合模擬真實的液體腐蝕環境；引入鈍化概率函數使得模擬的腐蝕坑的形貌隨機變化且不對稱，更加符合真實的腐蝕坑形貌；陰、陽極共同作用使電化學腐蝕的模擬更加完善；VOP方法可以控制腐蝕的快慢，而不是簡單的固、液相隨機轉換。

2）腐蝕反應速率越大，在相同的腐蝕時間下，在金屬表面點蝕的穩態過程中，金屬被腐蝕的體積越大，腐蝕的深度就越大；腐蝕含氧溶液擴散系數越大，在相同的腐蝕時間下，金屬被腐蝕的體積越大，腐蝕的深度就越大；腐蝕產物Fe3+擴散系數越大，在相同的腐蝕時間內，金屬被腐蝕的體積越小，腐蝕的深度就越小。