金屬腐蝕的多尺度計算模擬研究進展

魏 薪，董超芳，徐奧妮，李 妮，肖 葵，李曉剛

(北京科技大學(xué) 腐蝕與防護中心 教育部腐蝕與防護重點實驗室，北京 100083)

1 前 言

金屬的腐蝕受到材料因素、環(huán)境因素和界面狀態(tài)的影響，任何單一因素或耦合條件的改變都可能造成腐蝕規(guī)律的變化[1]，影響因素的復(fù)雜程度造成了腐蝕行為的多樣性，其機理研究的難度也正源于此。伴隨著計算機的發(fā)展，“計算和模擬”已經(jīng)成為了繼實驗、形式理論之后的第3條發(fā)現(xiàn)新概念和機理的途徑[2, 3]。在通過實驗建立的腐蝕模型基礎(chǔ)上，可以實現(xiàn)對點蝕[4-6]、晶間腐蝕[7]、應(yīng)力腐蝕[8]、沖蝕[9, 10]、擴散傳質(zhì)[11, 12]和電極過程動力學(xué)[13]等的理論研究。近年來，隨著試驗觀測技術(shù)的提高，掃描電化學(xué)顯微鏡(SECM)[14]、原子力顯微鏡(AFM)[15]、掃描開爾文力顯微鏡(SKPFM)[16]等微觀觀測儀器已經(jīng)應(yīng)用到了腐蝕研究中。這些儀器的觀測尺度一般在1 mm到1 μm的范圍內(nèi)，可以細致到觀測單個蝕坑的形貌，研究單個蝕坑周圍的電位，這對揭示腐蝕的深層機理有很大的促進。觀測手段的發(fā)展，為介觀尺度下的腐蝕過程的建模和仿真提供了試驗數(shù)據(jù)基礎(chǔ)。

另一方面，受反應(yīng)時間和尺度限制，大多數(shù)腐蝕實驗研究方法，并不能直接觀察材料的腐蝕萌生行為。研究者只能通過實驗獲得的數(shù)據(jù)進行分析進而推斷反應(yīng)發(fā)生過程，這一方法容易發(fā)生錯誤致使人們不能正確獲得真正的腐蝕機理，從而不能正確指導(dǎo)材料的設(shè)計與使用。原位測試技術(shù)的發(fā)展在某種程度上解決了部分問題，但由于其對試樣和實驗條件要求極高并不能廣泛的應(yīng)用。近年來，理論計算用軟件得到了迅速發(fā)展，使得腐蝕機理的研究取得了突破性的進展，對材料的設(shè)計和使用具有積極的指導(dǎo)作用。目前，腐蝕領(lǐng)域中應(yīng)用較廣泛的計算方法包括基于密度泛函的第一性原理和分子動力學(xué)方法、基于牛頓力學(xué)的分子動力學(xué)計算、基于電化學(xué)的熱力學(xué)和動力學(xué)計算以及基于斷裂力學(xué)的有限元方法等。

計算方法的選擇與所研究對象的尺度有關(guān)，已有從原子到宏觀尺度的方法包括：第一性原理計算、Monte Carlo方法、分子動力學(xué)、元胞自動機、有限元、邊界元等，其中第一性原理計算方法受到較大關(guān)注。Di Stefano等[17]利用基于密度泛函的第一性原理計算方法研究了∑3(111)[-110]晶界和∑5(210)[001]晶界與H的相互作用，結(jié)果表明：這兩種晶界與H原子的相互作用行為截然不同。密排的∑3晶界既不能捕獲H原子也不能促進其在晶界的擴散，在材料內(nèi)部充當二維的擴散阻擋層；而∑5晶界能夠提供大量的H溶解位，且相比于塊體材料，H原子沿∑5晶界的擴散系數(shù)提高了兩個數(shù)量級。Kulkoua等[18]利用第一性原理計算方法研究了H與B2結(jié)構(gòu)TiFe合金的相互作用，通過計算H原子晶界和自由表面的吸附能和偏聚能表明，H有助于晶間開裂的發(fā)生。

本文將對不同尺度下金屬腐蝕研究的計算方法進行總結(jié)，并以大氣環(huán)境下鋁合金的腐蝕與氫致開裂為計算對象，利用第一性原理和內(nèi)聚有限元的跨尺度計算方法，研究了氫原子偏聚導(dǎo)致的晶間裂紋的萌生與擴展行為，將第一性原理計算獲得的晶界結(jié)合能輸入到內(nèi)聚有限元計算中，實現(xiàn)了從原子尺度到宏觀尺度的跨尺度計算。

2 計算方法在金屬腐蝕研究中的應(yīng)用

2.1 第一性原理計算方法及其在腐蝕研究中的應(yīng)用

隨著量子理論的建立和計算機技術(shù)的發(fā)展，人們希望能夠借助計算機對微觀體系的量子力學(xué)方程進行數(shù)值求解，然而量子力學(xué)的基本方程—薛定諤方程的求解是極其復(fù)雜的。然而，電子密度泛函理論(DFT)的確立[19, 20]克服了這種復(fù)雜性。電子密度泛函理論是20世紀60年代在Thomas-Fermi理論的基礎(chǔ)上發(fā)展起來的量子理論的一種表述方。傳統(tǒng)的量子理論將波函數(shù)作為體系的基本物理量，而密度泛函理論則通過粒子密度來描述體系基態(tài)的物理性質(zhì)。密度泛函理論是一種完全基于量子力學(xué)的從頭算理論，但是為了與其他的量子化學(xué)從頭算方法區(qū)分，人們通常把基于密度泛函理論的計算叫做第一性原理。

在密度泛函理論中，所有的近似都被集中到稱為交換相關(guān)能的一項上，所以密度泛函理論的精度直接由交換相關(guān)能量泛函的近似形式?jīng)Q定。尋找更好的交換相關(guān)近似就成為密度泛函理論體系發(fā)展的一條主線。交換相關(guān)能量泛函的一個最初的簡單近似方法是局域密度近似(LDA)，即用具有相同密度的均勻電子氣的交換相關(guān)泛函作為對應(yīng)的非均勻系統(tǒng)的近似值，這樣一個簡單近似往往能給出很好的結(jié)果。如果進一步分別考慮不同自旋分量的電子密度，可得到自旋極化的局域密度近似L(S)DA。盡管L(S)DA獲得了巨大的成功，但是也有許多不足之處，比如系統(tǒng)地高估結(jié)合能。在L(S)DA基礎(chǔ)上的改進有廣義梯度近似(GGA)。在GGA近似下，交換相關(guān)能是電子密度及其梯度的泛函。構(gòu)造GGA交換相關(guān)泛函的方法分為兩個流派。一個是以Becke 為首的一派，認為“一切都是合法的”，人們可以任何理由選擇任何可能的泛函形式，而這種形式的好壞由實際計算來決定。通常，這類泛函的參數(shù)由擬合大量的計算數(shù)據(jù)得到。另外一個流派以Perdew 為首，他們認為發(fā)展交換相關(guān)泛函必須以一定的物理規(guī)律為基礎(chǔ)，這些規(guī)律包括標度關(guān)系、漸進行為等?；谶@種理念構(gòu)造的一個著名的GGA泛函是PBE泛函[21]，也是現(xiàn)在用得最廣泛的GGA泛函之一。電荷密度泛函理論在局域電荷密度近似或是梯度修正的版本，這是由Perdew和Wang所發(fā)展的GGA。DFT所描述的電子氣體交互作用被認為是對大部分的狀況都是精確的，并且他是唯一能實際有效分析周期性系統(tǒng)的理論方法。梯度修正的方法在研究表面的過程，小分子的性質(zhì)，氫鍵晶體以及內(nèi)部空間的晶體是比較準確的。眾所皆知，LDA會低估分子的鍵長或者鍵能，以及晶體的晶格參數(shù)，而GGA會補救這些缺點。而GGA在里斯晶體會過度修正LDA結(jié)果。當LDA與實驗符合非常好的時候，GGA會高估晶格長度。

利用第一性原理計算方法研究金屬材料在原子尺度的腐蝕行為對探究腐蝕機理有著十分重要的意義?；诿芏确汉挠嬎惴椒蓮脑又g相互作用的角度計算金屬材料的固有性質(zhì)，或描述金屬基體或鈍化膜表面與環(huán)境中分子或離子等之間的化學(xué)反應(yīng)。

2.1.1 表面吸附導(dǎo)致鈍化膜的形成與破壞

金屬暴露于各類環(huán)境下，與氣體或溶液中的氣體和離子發(fā)生相互作用，具體宏觀表現(xiàn)為金屬氧化膜的形成與破壞以及金屬表面與腐蝕性氣體和離子的相互作用。由此可見，吸附對腐蝕行為影響顯著，金屬表面對氣體或離子的吸附分為物理吸附和化學(xué)吸附。判斷被吸附物與表面是否發(fā)生化學(xué)反應(yīng)需要計算體系的電子性質(zhì)，如態(tài)密度、電荷密度、電子局域函數(shù)、布局分析、Bader電荷分析等。

本課題組[22]利用密度泛函方法計算了H2O和SO2在面心立方金屬Cu(111) 表面的單吸附和共吸附行為。通過對幾何構(gòu)型，能量，電荷轉(zhuǎn)移及態(tài)密度的分析可得到以下結(jié)論。單分子H2O和SO2均不能與化學(xué)鍵的形式吸附在表面。對于共吸附構(gòu)型，當覆蓋度為0.25 ML時，H2O和SO2分子分別吸附在Cu(111)表面；當覆蓋度增大到0.5 ML時，H2O分子解離為OH和H，其中OH吸附在表面，而H與SO2鍵合遠離表面。

Guo等[23]利用第一性原理計算表明，H2O分子在Al(111)純表面的解離能為248.32 kJ/mol，然而當存在預(yù)吸附O時，解離能降低到128.53 kJ/mol，說明O原子的存在大大降低了H2O在鋁表面的解離能?；诿芏确汉挠嬎隳軌蛴行У难芯縃2O與金屬的相互作用，但是在計算中卻存在一個嚴重的問題，那就是標準的密度泛函方法無法描述H2O分子與金屬表面的弱相互作用[24]。目前，已有大量研究[25,26]證明在考慮H2O分子與金屬表面相互作用時應(yīng)考慮范德華色散力的作用。本課題組[27,28]采用范德華密度泛函計算方法研究了O2和H2O分子在Al(111)表面的單吸附和共吸附行為。結(jié)果表明，O2分子可自發(fā)解離，以O(shè)原子的形式吸附在表面，伴隨著較大的吸附能；而H2O分子只能以分子的形式弱吸附在表面。同樣地，預(yù)吸附的O原子僅能促進H2O分子的變形，大大降低其解離能壘，不能使其解離吸附。當O2分子存在時，且兩個O原子同時與H2O分子中的H原子產(chǎn)生相互作用時，H2O分子解離為OH和H，吸附在金屬Al(111)表面。通過對O2和H2O分子的共吸附行為的研究可知，Al表面鈍化膜的形成機理為H2O分子的存在促進了吸附的原子O由表面向內(nèi)部的滲透而增厚，而不是金屬原子向外部的遷移。

Digne等[29]利用第一性原理計算的方法建立了γ-Al2O3表面模型。利用密度泛函方法，作者詳細研究了γ-Al2O3三種表面(100)、(110)和(111)的酸基性質(zhì)，并考慮了與溫度相關(guān)的羥基化表面的覆蓋度。模擬與大量實驗結(jié)合解決了在紅外光譜中出現(xiàn)的OH-伸縮頻率的問題。研究還確定了酸化表面位的精確本性。酸化的強度可以通過模擬目標分子，如CO，和計算表面電子性質(zhì)來量化。這些結(jié)果嚴重挑戰(zhàn)了以前的無缺陷尖晶石γ-Al2O3的表面模型并建立了γ-Al2O3酸化性質(zhì)更為嚴格的描述。本課題組的研究成果[30]利用考慮強關(guān)聯(lián)效應(yīng)的密度泛函方法研究了SO2分子在金屬鎳表面鈍化膜NiO(111)表面的吸附行為，結(jié)果表明，表面O空位的存在大大加強了表面對SO2分子的吸附，進一步說明表面鈍化膜的完整性對于材料耐蝕性的關(guān)鍵作用。這一系列的研究充分說明第一性原理可用于鋁合金與大氣環(huán)境中氣體相互作用的計算。

2.1.2 氫原子在界面處的偏聚

材料內(nèi)部存在大量的晶界以及第二相與基體的界面，這些界面屬于高能區(qū)。在腐蝕環(huán)境下，腐蝕優(yōu)先發(fā)生在界面處，且氫原子在界面的聚集導(dǎo)致氫脆的發(fā)生。因此，計算界面的結(jié)合能以及氫原子聚集對界面結(jié)合能的影響對研究應(yīng)力腐蝕、氫致開裂等局部腐蝕機理具有重要的意義。Fernandea等[31]利用密度泛函與基于過渡態(tài)理論和熱力學(xué)的靜態(tài)模型相結(jié)合的方法研究了氫原子與金屬鎢的相互作用。該模型得到的與溫度相關(guān)的數(shù)據(jù)有助于理解宏觀實驗獲得的結(jié)果。研究結(jié)果有效的解釋了實驗測得的H擴散系數(shù)與密度泛函計算獲得的結(jié)果存在偏差是由于空位濃度與溫度的差異。Kristinsdóttir等[32]利用密度泛函方法研究了氫原子在23種過渡態(tài)金屬密排表面的擴散行為。作者選擇了第4、5和6周期，從Sc到Au，除了Mn、Tc和Hf的d-金屬。氫原子在這些金屬表面的勢能面被建立，且從一個極小值到另一個極小值的能壘與微動彈性帶計算相比較。結(jié)果表明，所有的擴散激活能相對低或從Pt/0.04 eV到Y(jié)和Zr/0.28 eV。作者還估計了隧道效應(yīng)開始發(fā)生的溫度。Ferrina等[33]利用周期自洽DFT-GGA(PW91)方法研究氫原子與17個過渡族金屬表面，即面心立方金屬的(100)和(111)面，密排六方金屬的(0001)和體心立方金屬的(100)和(110)的相互作用。計算得到與實驗結(jié)果一致的幾何結(jié)構(gòu)和結(jié)合能。對于同一金屬的密排面和較開放的表面，計算得到的結(jié)合能有很大的不同。同一金屬不同亞表面的氫原子幾何結(jié)構(gòu)基本一樣，但是結(jié)合能有很大不同。除了Ag和Au表面H的吸附相對于H2時放熱的。對于所有研究的金屬，H在亞表面層的優(yōu)先吸附位沒有在表面的優(yōu)先吸附位穩(wěn)定。H從表面到第一層亞表面的擴散激活能由兩個表面熱力學(xué)穩(wěn)定性的差異控制。擴散從第一亞表面擴散進入下一層通常沒有很大的熱力學(xué)能壘，但有適當?shù)膭恿W(xué)能壘。對于同一種金屬，H在不同表面向內(nèi)擴散的激活能壘也有很大的不同。Mao等[34]利用第一性原理和分子動力學(xué)方法研究了間隙原子H在Er2O3內(nèi)的擴散，計算了H在塊體材料，(001)表面和 Σ13 (4-3-1)/[111] 扭轉(zhuǎn)晶界的擴散率和擴散能壘。結(jié)果表明，H在(001)表面的擴散比在塊體材料和晶界處快。由DFT和MD方法估計的H擴散能壘高于由實驗得到的沿晶界的能壘。這表明H在Er2O3內(nèi)的擴散由晶界控制。

2.2 有限元方法及其在腐蝕研究中的應(yīng)用

有限單元法是隨著電子計算機的發(fā)展而迅速發(fā)展起來的一種現(xiàn)代計算方法。它是50年代首先在連續(xù)體力學(xué)領(lǐng)域-飛機結(jié)構(gòu)靜、動態(tài)特性分析中應(yīng)用的一種有效的數(shù)值分析方法，隨后很快廣泛的應(yīng)用于求解熱傳導(dǎo)、電磁場、流體力學(xué)等連續(xù)性問題。有限元法是一種為求得偏微分方程邊值問題近似解的數(shù)值技術(shù)。它通過變分方法，使得誤差函數(shù)達到最小值并產(chǎn)生穩(wěn)定解。

ABAQUS是學(xué)術(shù)研究中使用最為廣泛的有限元分析軟件之一。ABAQUS基于其豐富的單元庫，可用于模擬各種復(fù)雜的幾何形狀，并且擁有豐富的材料模型庫，可用于模擬絕大多數(shù)常見的工程材料，如金屬、聚合物、復(fù)合材料、橡膠、可壓縮的彈性泡沫、鋼筋混凝土以及各種地質(zhì)材料。隨著實驗室實驗方法得到很大發(fā)展的同時，有限元模擬方法也取得了很大的進展，除一般的線性和非線性有限元模擬，還出現(xiàn)一些新的有限元模擬技術(shù)。本文主要介紹內(nèi)聚力單元(Cohesive elements)的有限元模擬方法。

ABAQUS提供一種內(nèi)聚力單元，用以模擬兩個部分之間的粘性連接，一般來說，它要求粘結(jié)材料尺寸和強度都小于粘結(jié)部分(比如多層復(fù)合材料的膠粘層)，進而可以利用內(nèi)聚力單元模擬材料的斷裂。

這個單元似乎與C3D8實體單元的結(jié)構(gòu)一致，但最大的區(qū)別在于橫縱向尺寸的比例，C3D8單元在模擬這種大橫縱比結(jié)構(gòu)的時候，已經(jīng)無法給出精確的解答了。這就是內(nèi)聚力單元存在的意義，此外，從內(nèi)聚力單元的變形機理來看，中面雖然能夠承受拉伸和剪切應(yīng)變，但并不能產(chǎn)生任何應(yīng)力。因此，內(nèi)聚力單元只能支持垂直于上下表面的牽引-分離破壞準則。當需要模擬諸如扁平結(jié)構(gòu)破壞或者某結(jié)構(gòu)沿一個平面撕裂這種工況的時候，內(nèi)聚力單元提供了一種無需細化網(wǎng)格的解決方案，可以有效的降低網(wǎng)格數(shù)量，提高運算效率。

基于力學(xué)的有限元模擬在腐蝕中應(yīng)用的研究并不是很多，大多數(shù)研究集中在材料表面和內(nèi)部裂紋萌生和擴展問題，尤其是氫敏感體系。Simonovski等[35]提出了利用內(nèi)聚單元法在晶粒尺度模擬晶間開裂行為的方法。該方法將0厚度的內(nèi)聚力單元嵌入晶粒之中代表晶界。利用traction-seperation法則進行應(yīng)力和晶界損傷的計算，結(jié)果表明，該方法能夠有效的評價多晶體的晶間損傷和開裂行為。Scheider等[36]用內(nèi)聚力有限元模型研究了氫擴散對裂紋擴展的影響。Pouillier等[37]利用內(nèi)聚力模型模擬鋁合金的氫致晶間開裂行為，所得結(jié)果與電子背散射技術(shù)測得的結(jié)果一致。

綜上所述，計算模擬方法在腐蝕領(lǐng)域得到了廣泛的應(yīng)用，并取得了突破性的進展。計算模擬可為實驗結(jié)果提供理論指導(dǎo)，也可為數(shù)據(jù)分析提供依據(jù)。因此，計算模擬方法是研究腐蝕行為與機理的有效手段。

2.3 跨尺度計算與模擬方法研究進展

金屬結(jié)構(gòu)材料在環(huán)境中腐蝕的過程從表面對腐蝕性離子的吸附、點蝕與微裂紋的萌生到點蝕坑與裂紋的形成，整個過程從肉眼難以察覺的原子尺度到可觀察的宏觀尺度。通常來講，受到反應(yīng)時間和尺度的限制，原子尺度的反應(yīng)很難利用實驗手段進行觀測，這也導(dǎo)致原子尺度的計算與模擬結(jié)果很難得到驗證。因此，建立原子尺度到宏觀尺度的跨尺度的計算方法，將原子尺度的輸出結(jié)果輸入到宏觀尺度的計算之中，這樣既解決了小尺度實驗難以實現(xiàn)的問題，又可對宏觀尺度的結(jié)果進行實驗驗證進而間接驗證原子尺度的計算結(jié)果。除此之外，腐蝕環(huán)境復(fù)雜多變，利用跨尺度的計算與模擬方法能夠?qū)我灰蛩氐挠绊戇M行系統(tǒng)研究。

近年來，跨尺度計算已經(jīng)出現(xiàn)在不同領(lǐng)域的研究之中。多尺度的計算方法一般可分為尺度分離和尺度間耦合兩種，前者著眼于采用不同尺度來分析對象的不同部分，后者注重尋找尺度之間的關(guān)聯(lián)[38]。通常來講，尺度之間的劃分為宏觀(>10-4m)、介觀(介于10-6m與10-4m之間)、微觀(介于10-7m到10-6m)和納觀(<10-7m)[39]。對于腐蝕計算而言，既需要利用不同尺度的計算方法研究腐蝕發(fā)生的不同階段，也需要探究不同尺度之間的耦合關(guān)系，通過選擇合適的計算方法對腐蝕的行為與機理進行研究。目前，常見的跨尺度計算模型大多集中在對材料的塑性形變和斷裂過程的研究，如本課題組[40]利用第一性原理與內(nèi)聚有限元相結(jié)合的跨尺度計算方法研究了鋁合金的晶間開裂行為。首先，利用第一性原理研究了H原子在Al內(nèi)部和晶界處的吸附、溶解和擴散行為，結(jié)果表明H原子優(yōu)先溶解在晶界間隙位且晶界為H原子的快速擴散通道；其次，將不同濃度的H原子引入到晶界間隙位并計算其晶界結(jié)合能，結(jié)果表明晶界結(jié)合能，即晶界強度隨H原子濃度的增加而線性下降；最后，將第一性原理計算獲得的晶界結(jié)合能作為斷裂能輸入到內(nèi)聚有限元計算之中，模擬鋁合金的晶間開裂行為。

3 鋁合金腐蝕與氫致開裂的模擬計算

3.1 鋁鈍化膜點蝕萌生機理

鋁表面覆蓋的氧化物在溶液中或在潮濕空氣中，會吸附環(huán)境中的H2O分子形成羥基化的水合表面。當溶液中存在Cl-時，它會與溶液中的OH-在鈍化膜表面發(fā)生競爭吸附形成Me-Cl-或是MeO(H)-Cl-的產(chǎn)物，產(chǎn)物的生成導(dǎo)致鈍化膜中的金屬離子與表面的結(jié)合力降低從而易于離開表面溶于溶液之中。隨著吸附反應(yīng)的發(fā)生，更多的金屬陽離子遷移到溶液中導(dǎo)致鈍化膜快速減薄，降低其耐蝕性能。關(guān)于鈍化膜的局部減薄造成點蝕的機理已被廣泛討論過，但尚未達成一致的說法，尤其是在原子尺度的機理研究還未見報道。本部分利用基于密度泛函的第一性原理計算研究Cl-在金屬鋁鈍化膜表面的吸附和向亞表面滲透的過程。通過對取代和滲透反應(yīng)的構(gòu)型、能量和能帶結(jié)構(gòu)的分析來探究Cl-與鈍化膜反應(yīng)的本質(zhì)及其導(dǎo)致鈍化膜破裂的機理。

圖1a為Cl-取代羥基化γ-Al2O3(010)表面OH-的構(gòu)型圖。隨著數(shù)值OH-/Cl-的減小，表面發(fā)生不同程度的重構(gòu)。對于OH∶Cl>8∶8，Cl-取代的OH-來自解離的H2O，對氧化膜本身的結(jié)構(gòu)未產(chǎn)生明顯的影響。對于OH∶Cl=2∶14，有14/16的表面OH-吸附位被取代，由構(gòu)型圖可知表面第一層Al原子均被拉出表面，發(fā)生明顯的重構(gòu)現(xiàn)象。對于OH∶Cl=0∶16來說，表面OH-均被Cl-取代，第一層表面原子與Cl-形成腐蝕產(chǎn)物。圖1b為隨表面Cl-濃度的增加Cl-取代OH-的能量變化，當表面Cl-濃度較低時(小于OH∶Cl=12∶4中所示濃度)，隨著表面Cl-濃度的增大取代能負移，說明Cl-濃度的增加促進其在表面的吸附；當表面Cl-濃度介于OH∶Cl=12∶4和OH∶Cl=8∶8之間時，隨著Cl-濃度的增加取代能正移，說明當表面存在一定吸附Cl-時，后續(xù)的Cl-吸附反應(yīng)難度增加；當表面Cl濃度大于OH∶Cl=8∶8所示，Cl開始取代表面OsH，取代的難度進一步加大。造成能量正移的原因為表面吸附Cl-之間的靜電斥力。當表面全部被Cl-取代時，取代能略有降低，但仍明顯高于低Cl-濃度時。這是因為全部取代使得表面處于對稱結(jié)構(gòu)，降低了靜電斥力間的不均勻，使得吸附結(jié)構(gòu)變得穩(wěn)定。綜上所述，吸附的Cl-和OH-對后續(xù)Cl-的取代反應(yīng)均有阻礙作用，控制吸附反應(yīng)發(fā)生的因素為表面吸附離子之間的靜電斥力。

Cl-取代表面OH-的反應(yīng)為放熱反應(yīng)，隨著取代Cl-濃度的增加，取代能正移，反應(yīng)驅(qū)動力減小。若Cl-取代羥基化γ-Al2O3(010)表面的吸附的來自H2O解離的OH-，則對鈍化膜的結(jié)構(gòu)無明顯影響；若取代鈍化膜本身與H+結(jié)合成的OH-時，則發(fā)生明顯結(jié)構(gòu)變化，隨著取代Cl-濃度的增大，表面Al與Cl-結(jié)合生成無序氯化物層。

圖1 Cl-在水合γ-Al2O3(010)表面發(fā)生取代反應(yīng)：(a)原子構(gòu)型圖; (b)能量隨表面Cl-濃度的變化Fig.1 Substitution reaction between chloride ions and hydrous γ-Al2O3(010) surface: (a)atomic configuration diagram; (b)variation of energy as chloride ions’ concentration

3.2 氫原子的吸附、溶解和擴散

H原子能夠在鋁表面吸附并進入材料內(nèi)部，溶解在間隙位。Al (111)表面共有4類對稱吸附位，分別為fcc位、hcp位、top位和bridge位。通過計算表明單一H原子僅能吸附在表面top位和fcc位，當H原子在hcp和bridge位時會自動遷移到最近鄰fcc位。對于H原子在材料內(nèi)部的溶解，本節(jié)共研究了Al(111)層狀表面的5種間隙溶解位，如圖2所示，分別為表面與亞表面之間的兩類四面體間隙位和一類八面體間隙位；塊體材料內(nèi)部的四面體和八面體間隙位。

圖3a和3b為H原子在晶界∑5(100)/36.87°不同溶解間隙位的遷移路徑和Al原子空位的位置，其中I3GB為靠近晶界的晶內(nèi)溶解位，I4-1GB和I4-2GB為等效最優(yōu)溶解位。本節(jié)選取了兩條擴散路徑，分別為I1GB→I4-1GB和I5GB→I4-1GB。由圖3c和3d，對于互相垂直的兩條遷移路徑I1GB→I4-1GB和I5GB→I4-1GB，Al原子空位的存在并未降低其擴散能壘，而是大幅度的提高了擴散能壘。當不存在Al原子空位時，I1GB→I4-1GB的遷移能壘為12.96 kJ/mol。Al原子空位的存在使其遷移能壘提高到76.41 kJ/mol，能壘大幅度升高63.45 kJ/mol。對于I5GB→I4-1GB的遷移路徑，無空位存在時的能壘為9.31 kJ/mol，空位的出現(xiàn)使其能壘降低到3.55 kJ/mol，說明Al空位促進H原子沿該路徑的遷移。由以上研究可知，Al原子空位的存在對H原子在晶界的遷移具有很大的影響，Al原子空位能夠促進晶內(nèi)的H原子遷移進入晶界；Al原子空位對晶界內(nèi)H原子遷移的作用具有方向性，在某些方向能夠促進H原子的遷移，而在某些方向卻大大抑制H原子的遷移。

圖2 H原子沿Al(111)表面吸附和向內(nèi)部溶解的間隙位：(a)TIS1-1，表面與亞表面之間的第一類四面體間隙；(b)TIS1-2，表面與亞表面之間的第二類四面體間隙位；(c)OIS1，表面與亞表面之間八面體間隙位；(d)TIS2，金屬Al內(nèi)部四面體間隙位；(d)OIS2，金屬Al內(nèi)部八面體間隙位[40]Fig.2 Interstitial sites of H atom adsorbed and dissolved along Al(111) surface: (a) TIS1-1, the first type of tetrahedral interstitial site between surface and subsurface; (b) TIS1-2, the second type of tetrahedral interstitial site between surface and subsurface; (c) OIS1, octahedral interstitial site between surface and subsurface; (d) OIS2, octahedral interstitial site in Al bulk[40]

3.3 氫致鋁晶間開裂行為的跨尺度計算

利用Neper軟件產(chǎn)生多晶體結(jié)構(gòu)模型，如圖4所示。晶粒之間的界面利用FORTRAN語言編寫的程序產(chǎn)生。多晶體模型包含300個晶粒，1387個晶界面。線性四面體單元(ABAQUS 類型 C3D4) 被用于劃分晶內(nèi)部分網(wǎng)格, 線性三角單元(ABAQUS類型 COH3D6) 被用于劃分晶界部分網(wǎng)格。 多晶體晶粒的取向為正交各向異性。沿x方向的10-3mm的位移邊界條件被施加在右側(cè)表面，左側(cè)表面的節(jié)點沿x、y、z三個方向固定。

圖3 H原子在晶界(100)/36.87°，∑5的溶解位示意圖：(a)側(cè)視圖；(b)俯視圖(其中箭頭表示為不同的遷移路徑，黃色四邊形標記處為Al原子空位的位置)。H原子在鋁晶界間隙位擴散的路徑：(c)原子在有無空位情況下沿I1GB→I4-1GB的遷移路徑；(d)H原子在有無空位情況下沿I5GB→I4-1GB的遷移路徑(綠色箭頭表示的數(shù)值為有無空位情況下遷移能壘的差值)Fig.3 Dissolved sites of H atom in (100)/36.87°，∑5 interface: (a) side view; (b) top view; (c) migration path of H atom along I1GB→I4-1GB with and without vacancies; (d) migration path of H atom along I5GB→I4-1GB with and without vacancies. Green arrows represent the difference values between migration energy barrier with vacancies and without vacancies

圖4 多晶體模型圖Fig.4 Polycrystal model diagram

利用第一性計算原理方法研究晶界偏聚對晶界強度的影響，將不同濃度H原子置于Al∑5 [100]/36.87°晶界超胞穩(wěn)定吸附位，計算其晶界結(jié)合能，表明晶界結(jié)合能隨H濃度升高而線性降低，如圖5所示，即晶界強度線性降低。H原子對晶界的弱化效應(yīng)歸因于原子H的存在弱化了晶界處Al原子之間的金屬鍵。

圖5 Al∑5 [100]/36.87°晶界結(jié)合能隨H濃度增加的變化Fig.5 The relationship between H concentration and cohesive energy of Al GB ∑5 [100]/36.87°

當一定量的H原子從環(huán)境中進入材料內(nèi)部，H原子能夠偏聚在晶界，導(dǎo)致原子之間結(jié)合力降低，這是鋁合金失效的主要原因之一；雖然局部晶界應(yīng)變已能夠通過實驗測試獲得，但應(yīng)力的分布仍然難以測得,內(nèi)聚力有限元方法能夠?qū)Ы绲膽?yīng)力應(yīng)變分布和損傷行為進行細致的計算與模擬，其中最重要的輸入?yún)?shù)—晶界結(jié)合能，能夠通過第一性原理計算獲得。通過內(nèi)聚有限單元法對鋁多晶體模型的晶間開裂行為進行模擬，其中利用第一性原理計算得到的晶界結(jié)合能將作為斷裂能這一參數(shù)輸入到內(nèi)聚有限元計算中，實現(xiàn)原子尺度與宏觀尺度之間的銜接，從而實現(xiàn)跨尺度計算。圖6為晶界單元的損傷云圖，可以看出，隨著晶界結(jié)合能的降低，失效單元增多，隨著失效單元的累積使得微觀晶間出現(xiàn)微觀裂紋。

圖6 不同晶界結(jié)合能下內(nèi)聚有限單元法計算得到的SDEG損傷分布云圖[40]Fig.6 SDEG damage distribution diagrams of polycrystal model with different binding energy of grain boundary via cohesion finite element method[40]

4 結(jié) 語

本文對計算模擬方法在金屬腐蝕研究中的應(yīng)用與進展進行了綜述，將第一性原理與有限元計算相結(jié)合，建立了從原子尺度到宏觀尺度的跨尺度計算方法。本文以鋁及鋁合金為例，利用第一性原理計算了Cl與表面鈍化膜的相互作用來研究鈍化膜的點蝕機理。本文還建立了Al-H體系，研究了H原子在Al表面和內(nèi)部的吸附、溶解，擴散與偏聚，并計算H原子偏聚對晶界強度的影響。最后將原子尺度第一性原理計算獲得的能量輸入到宏觀尺度的有限元計算之中模擬開裂行為。本研究基本完成鋁大氣腐蝕行為與機理的跨尺度計算與模擬框架的搭建，此方法適用于服役于大氣環(huán)境中具有Cl-點蝕敏感性和氫脆敏感性的結(jié)構(gòu)材料，對工程應(yīng)用具有指導(dǎo)作用。

本課題組在后續(xù)的工作中還針對金屬鈍化膜的性能開展了大量的研究，利用第一性原理與電化學(xué)測試相結(jié)合的方法研究鈍化膜的穩(wěn)定性[41,42]及其半導(dǎo)體性質(zhì)對其耐蝕性能的影響。

References

[1] Xiao Jimei(肖紀美), Cao Chunan(曹楚南).PrinciplesofMaterialCorrosion(材料腐蝕學(xué)原理)[M]. Beijing: Chemical Industry Press, 2002.

[2] Zhu X, Zi G.Construction&BuildingMaterials[J], 2017, 137: 330-344.

[3] Maulana A, Dipoyono H K, Khairurrijal K.IndonesianJournalofPhysics[J], 2015, 18(1): 5-9.

[4] Wang H, Han E H.ElectrochimicaActa[J], 2013, 90(5): 128-134.

[5] Laycock N J, Krouse D P, Hendy S C,etal.ElectrochemicalSocietyInterface[J], 2014, 23(4): 65-71.

[6] Wang H T, Han E H.Materials&Corrosion[J], 2015, 66(9): 925-930.

[7] Harrison T J, Crawford B R, Brandt M,etal.ComputationalMaterialsScience[J], 2014, 84(1): 74-82.

[8] Li Q K, Zhang Y, Shi S Q,etal.MaterialsLetters[J], 2012, 56(6): 927-932.

[10] Zhang G A, Zeng L, Huang H L,etal.CorrosionScience[J], 2013, 77: 334-341.

[11] Liu Zhixiao(劉智驍).TransactionsofNonferrousMetalsSocietyofChina(中國有色金屬學(xué)報)[J], 2013, 23(4): 1160-1167.

[12] Pan Changchang(潘昌昌). Dissertation for Master[D]. Lanzhou：University of Lanzhou Technology, 2016.

[13] Venkatraman M S, Cole I S, Emmanuel B.ElectrochimicaActa[J], 2011, 56(20): 7171-7179.

[14] Indira K, Nishimura T.JournalofMaterialsEngineering&Performance[J], 2016, 25(10): 1-14.

[15] Giner I, Keller A.CorrosionScience[J], 2015, 100: 496-503.

[16] Ceylan S, Sergiy B, Martin S,etal.CorrosionScience[J], 2012, 58(5): 307-314.

[17] Di Stefano D, Mrovec M, Els?sser C.ActaMaterialia[J], 2015, 98: 306-312.

[18] Kulkova S E, Kulkov S S, Bakulin A V,etal.InternationalJournalofHydrogenEnergy[J], 2012, 37(8): 6666-6673.

[19] Hohenberg P, Kohn W.PhysicalReview[J], 1964, 136：B864-B871.

[20] Kohn W, Sham L J.PhysicalReview[J], 1965, 140(4A): A1133-A1138.

[21] Perdew J P, Burke K, Ernzerhof M.PhysicalReviewLetters[J], 1996, 77(18): 3865-3868.

[22] Wei Xin, Dong Chaofang, Chen Zhanghuaetal.TransactionsofNonferrousMetalsSocietyofChina[J], 2015, 25: 4102-4109.

[23] Guo F Y, Long C G, Zhang J,etal.AppliedSurfaceScience[J], 2015, 324: 584-589.

[24] Carrasco J, Hodgson A, Michaelides A.NatureMaterials[J], 2012, 11(8): 667-674.

[26] Tonigold K, Gross A.JournalofComputationalChemistry[J], 2012, 33(6): 695-701.

[27] Wei X, Dong C F, Chen Z H,etal.RSCAdvances[J], 2016, 6:79836-79843.

[28] Wei X, Dong C F, Chen Z H,etal.RSCAdvances[J], 2016, 6: 56303-56312.

[29] Digne M, Sautet P, Raybaud P,etal.JournalofCatalysis[J], 2012, 226(1): 54-68.

[30] Wei X, Dong C, Chen Z,etal.AppliedSurfaceScience[J], 2015, 355(15): 429-435.

[31] Fernandez N, Ferro Y, Kato D.ActaMaterialia[J], 2015, 94(1): 307-318.

[32] Kristinsdóttir L, Skúlason E.SurfaceScience[J], 2012, 606:1400-1404.

[33] Ferrin P, Kandoi S, Nilekar A U,etal.SurfaceScience, 2012, 606(7-8): 679-689.

[34] Mao W, Chikada T, Suzuki A,etal.JournalofNuclearMaterials[J], 2014, 455(1): 360-365.

[35] Simonovski I, Cizelj L.EngineeringFractureMechanics[J], 2013, 110: 364-377.

[36] Zhu L K, Yan Y, Qiao L J,etal.CorrosionScience[J], 2013, 77(77): 360-368.

[37] Pouillier E, Gourgues A F, Tanguy D,etal.InternationalJournalofPlasticity[J], 2012, 34(34): 139-153.

[38] Lu Xinzheng(陸新征), Lin Xuchuan(林旭川), Ye Lieping(葉列平).JournalofHuazhongUniversityofScienceandTechnology(華中科技大學(xué)學(xué)報)[J], 2008, 25(4): 76-79.

[39] Guo Yafang(郭雅芳), Wang Chongyu(王崇愚).MaterialsReview(材料導(dǎo)報)[J]，2001，15(7)：9-11.

[40] Wei X, Dong C F, Chen Z H,etal.RSCAdvances[J], 2016, 6:27282-27292.

[41] Xu A N, Dong C F, Wei X,etal.ElectrochemistryCommunications[J], 2016, 68: 62-66.

[42] Kong D C, Xu A N, Dong C F,etal.CorrosionScience[J], 2017, 116: 34-43.