海水環(huán)境下金屬腐蝕過程數(shù)值模擬*

馮麗, 蔡琦（.海軍工程大學(xué)核能科學(xué)與工程系，武漢430033；2.海軍上海監(jiān)修室，上海 20036）

0 引言

我國南海海域有著“高溫、高濕、高鹽”的“三高”氣候特點(diǎn)，其年平均氣溫、平均濕度和高溫持續(xù)時(shí)間高于東海、黃海和北海海域，同時(shí)由于高溫蒸發(fā)等因素的影響，其表面海水鹽度也高于其它海域。金屬材料在南海“三高”環(huán)境下的腐蝕現(xiàn)象尤為突出[1]。材料的腐蝕形貌是評(píng)價(jià)材料腐蝕性能的重要特征。當(dāng)前腐蝕研究中主要采用概率統(tǒng)計(jì)方法，即首先通過試驗(yàn)方法獲得材料腐蝕數(shù)據(jù)，再利用數(shù)字圖像處理技術(shù)得到腐蝕形貌圖像，進(jìn)而研究試樣腐蝕形態(tài)和數(shù)據(jù)間的關(guān)系。然而腐蝕是材料和環(huán)境作用造成的緩慢變化過程，變化速度往往以年甚至更長時(shí)間為單位，導(dǎo)致試驗(yàn)過程需要更多的時(shí)間和精力。采用數(shù)值方法模擬三維表面即正向建模方法，對(duì)試樣腐蝕損傷過程進(jìn)行數(shù)值模擬，能夠避免繁雜耗時(shí)的試驗(yàn)過程，便于模擬不同環(huán)境因素對(duì)材料腐蝕過程的影響。目前有少數(shù)二維或三維腐蝕模型[2-3]采用正向建模方法，然而建立的三維腐蝕模型，均未考慮材料腐蝕前表面缺陷或粗糙度的影響，將材料表面作平面處理，得到的數(shù)值模擬結(jié)果勢必與真實(shí)表面腐蝕結(jié)果存在差異。

本文首先基于分形理論用MATLAB編程生成三維W-M分形表面，然后利用三維元胞自動(dòng)機(jī)方法在模擬的分形表面基礎(chǔ)上，對(duì)海水濃度變化情況下金屬腐蝕過程進(jìn)行數(shù)值模擬，研究不同海水濃度下，金屬材料腐蝕前后表面形貌特征及其變化規(guī)律。

1 分形表面模擬

MANDELBROT最早將W-M函數(shù)引入分形領(lǐng)域，用于模擬二維粗糙表面輪廓高度，Komvopoulos[4]通過函數(shù)變換推導(dǎo)出直角坐標(biāo)系下的三維表面分形函數(shù)：

由上式對(duì)隨機(jī)粗糙表面進(jìn)行模擬，如圖1所示。

2 腐蝕模型

元胞自動(dòng)機(jī)（CA）是一種用來仿真局部規(guī)則和局部聯(lián)系的方法。元胞自動(dòng)機(jī)已被應(yīng)用于物理模擬、生物模擬等領(lǐng)域。許多問題通過簡化都可以轉(zhuǎn)換成二維問題，因此二維元胞自動(dòng)機(jī)方法是應(yīng)用范圍最廣的元胞自動(dòng)機(jī)方法。由于三維元胞自動(dòng)機(jī)方法能更真實(shí)地模擬實(shí)際情況，近年來，隨著研究問題的深入及計(jì)算機(jī)能力的加強(qiáng)，研究者們?cè)絹碓蕉嗟亻_始采用三維元胞自動(dòng)機(jī)方法來模擬各種工程實(shí)際問題[5-14]。

圖1 W-M分形表面模擬（G=1.36×10－2nm,M=10,γ=1.5，D=2.4）

本文在模擬的三維分形表面基礎(chǔ)上采用三維元胞自動(dòng)機(jī)方法，對(duì)海水濃度變化對(duì)金屬表面腐蝕形貌的影響進(jìn)行數(shù)值模擬。

2.1 數(shù)值模擬過程

數(shù)值模擬過程如下：

1）首先利用MATLAB編程，在x和y方向1 nm×1 nm平面內(nèi)（其中0≤x≤1,0≤y≤1），由W-M分形函數(shù)在100×100平面網(wǎng)格上生成三維分形表面，此表面作為腐蝕前即t=0時(shí)刻的表面形貌。

2）然后采用三維元胞自動(dòng)機(jī)方法在1 nm×1 nm×1 nm立方體內(nèi)（其中-0.5≤z≤0.5），將金屬和海水系統(tǒng)劃分為100×100×100網(wǎng)格空間，使金屬材料和海水系統(tǒng)離散劃分成有序的元胞，本文模擬系統(tǒng)中考慮三種類型元胞：金屬元胞（金屬原子）、腐蝕元胞（氯離子）和非腐蝕元胞（水分子）。

圖2 c=0.400時(shí)元胞自動(dòng)機(jī)模型網(wǎng)格圖

圖3 c=0.035時(shí)z=0截面腐蝕后原子分布圖

圖4 c=0.100時(shí)z=0截面腐蝕后原子分布圖

圖5 c=0.400時(shí)z=0截面腐蝕后原子分布圖

圖6 海水濃度c=0.400時(shí)分形表面腐蝕前表面形貌

圖7 海水濃度c=0.400時(shí)分形表面腐蝕后表面形貌

3）設(shè)定邊界條件和各相關(guān)參數(shù)，用MATLAB編程對(duì)每個(gè)時(shí)間步長進(jìn)行迭代計(jì)算，得到不同海水濃度下金屬材料腐蝕前后表面形貌數(shù)據(jù)。為了提高運(yùn)算速度、簡化模擬程序，文中腐蝕模型的算法采用Than和Buttgenbach等[15]提出的觀點(diǎn)，在模擬程序運(yùn)行時(shí)，某一時(shí)刻只考慮金屬材料與海水接觸面上相關(guān)的元胞。

2.2 計(jì)算實(shí)例

程序輸入中各參數(shù)取下列值：γ=1.5、M=10、D=2.8、G=1.36×10-2nm、△γ=0.25 J/m2、海水濃度c為0.035、0.070、0.100、0.200、0.400。圖2顯示了c=0.400時(shí)元胞自動(dòng)機(jī)模型網(wǎng)格圖，圖3~圖5顯示了c為0.035、0.100、0.400時(shí)相同時(shí)間步長下z=0截面腐蝕后原子分布圖（其中，黃色代表金屬原子，紅色代表氯離子，藍(lán)色代表水分子），圖6、圖7顯示了海水濃度c=0.400時(shí)金屬表面腐蝕前后表面形貌。

2.3 結(jié)果分析

任何表面在微觀上都是由許多不規(guī)則的凸峰和凹谷組成的。文中在模擬的W-M分形表面基礎(chǔ)上，由元胞自動(dòng)機(jī)方法生成的腐蝕前模型部分網(wǎng)格圖如圖2所示，凹凸不平的粗糙表面金屬原子與海水中的水分子和氯離子接觸，氯離子隨機(jī)分布在海水中。

數(shù)值迭代后，由圖3~圖5可以看出，隨著海水濃度的增大，氯離子數(shù)量隨之增加，加大了金屬與海水接觸界面上氯離子與金屬原子接觸的可能性。當(dāng)氯離子與表面金屬原子接觸并向金屬原子的方向運(yùn)動(dòng)時(shí)將發(fā)生腐蝕現(xiàn)象，氯離子將替換金屬原子，金屬原子則被腐蝕消耗掉。

海水濃度c=0.400時(shí)分形金屬表面腐蝕前后表面形貌如圖6、圖7所示，可見腐蝕現(xiàn)象發(fā)生后，更多的表面金屬原子被水分子或氯離子替換而被腐蝕消耗掉。

3 結(jié)論

1）在W-M分形函數(shù)基礎(chǔ)上對(duì)海水濃度變化情況下金屬腐蝕過程進(jìn)行數(shù)值模擬，為考慮表面缺陷或表面粗糙度效應(yīng)的腐蝕規(guī)律和腐蝕機(jī)理的進(jìn)一步研究提供了思路。

2）由圖3~圖5可以看出，隨著海水濃度的增加，金屬表面腐蝕原子數(shù)量增多，表面腐蝕進(jìn)程加快，越來越多的金屬原子被氯離子腐蝕消耗掉，而被水分子或氯離子填充，因而表面金屬原子逐漸減少，水分子和氯離子數(shù)量逐漸增多，符合客觀規(guī)律，該模型的建立有其合理性。

3）由圖6~圖7可以看出，腐蝕前后金屬表面形貌發(fā)生明顯變化，表面金屬原子逐漸被氯離子腐蝕消耗掉。