基于模擬退火算法的真實多孔電極中熱-質傳遞的研究

時間：TIME\@"yyyy'年'M'月'd'日'"2022年3月29日學海無涯頁碼：第1-頁共1頁基于模擬退火算法的真實多孔電極中熱-質傳遞的研究電場作用下電極中電化學反應-熱質傳遞現(xiàn)象是典型的多尺度問題[1-3]，如圖1所示。微觀尺度(電子-離子遷移、晶格穩(wěn)定性)、介觀尺度(界面熱/動力學、熱-質傳遞、電流電壓分布)和宏觀尺度(散熱性能、充放電管理)的傳遞和反應特性，直接決定系統(tǒng)中的濃度分布和反應速率，最終影響儲能轉化效率。而系統(tǒng)中發(fā)生傳遞和反應過程的場所主要在多孔電極中，且多孔電極的不規(guī)則表面引起的限域效應顯著影響著其內的傳質和傳熱過程，所以深入研究介觀尺度下多孔電極中的熱-質傳遞現(xiàn)象及其耦合機制，對高性能電化學儲能器件的設計有重要意義[3]。但是，與平板電極不同，對多孔電極中熱-質傳遞現(xiàn)象的研究存在兩個關鍵問題：一方面，多孔電極豐富的孔道體系、多樣化的表面形態(tài)和分布復雜的催化活性位點，難以用傳統(tǒng)數學工具去準確描述其結構[4-6]；另一方面，多孔電極中發(fā)生的過程涉及在受限空間中的離子遷移、熱量產生和傳遞、雙電層形成、電子傳遞的耦合，極大地影響了離子和熱量在孔道中的傳遞過程，并難以直接用實驗進行原位表征[7-8]。因此，多孔電極中的熱-質傳遞的研究一直以來都難以獲得突破性進展。

圖1

圖1儲能和轉換過程中的電極尺度問題及與其他尺度問題的聯(lián)系

Fig.1Theelectrodescaleissuesinenergystorageandconversionprocessesandtheirrelationtootherscaleissues

目前，研究者們可以采用密度泛函理論(DFT)[9-10]和分子動力學(MD)模擬[11-12]等方法研究納米級多孔結構中的演化過程和充電機理。例如納米多孔電極與離子液體的結合是提高超級電容器能量密度的重要手段。然而，這總是伴隨著功率密度的降低，特別是考慮到高黏度和大空間位阻的離子液體。Gan等[11]利用MD模擬發(fā)現(xiàn)在閾電位作用下，具有疏離子孔的電極內的離子呈現(xiàn)出一種新的充電機制，即離子吸附。進一步得到了充電時間/電容與電壓/耐離子性能之間的定量關系，以評價協(xié)同提高超級電容器能量密度和功率密度的臨界條件。雖然DFT和MD模擬可以考慮微觀粒子間的作用力，但是受限于計算能力的瓶頸，只能模擬離子濃度和溫度在含周期性邊界條件的高度有序的納米級多孔結構(如分子篩和金屬有機框架)中納秒級的演化過程，并不能全面描述真實多孔電極中，尤其是常見的無定形碳電極中的離子和熱量的傳遞過程[9]。連續(xù)介質模型通過偏微分方程組描述離子濃度和電勢的狀態(tài)演化及其本構關系，可以建立納米級及以上空間尺度的復雜的多物理場耦合關系[13]。例如，Tao等[14]采用基于Poisson-Nernst-Planck(PNP)方程和Navier-Stokes方程的組合方法，通過控制氧化石墨烯膜的偏置角度，研究了不同層間距分布的氧化石墨烯膜中離子的輸運現(xiàn)象。d′Entremont等[15]從第一性原理推導出控制能量方程，并與修正PNP模型相結合，從而得到由于離子擴散、空間效應和混合熵的變化而產生的不可逆焦耳熱和可逆熱產生速率。但是，目前連續(xù)介質模型在多孔電極中的大規(guī)模運用仍存在兩個難點：準確描述多孔介質的復雜結構和快速求解含有復雜邊界條件的偏微分方程組。因此，多種簡化的多孔電極表示方法和經驗模型被廣泛運用。天然形成的和部分人工制造的多孔介質的微結構具有分形特征，可用分形維數表示。Sakaguchi等[16]利用PNP方程研究了分形多孔電極中的充電過程，并發(fā)現(xiàn)分形多孔介質充注過程的時間演化遵循冪律，指數與分形維數有關。目前大量研究都集中在單個孔隙中的充放電過程，因此可以通過等效介質近似(effectivemediumapproximation,EMA)方法[17]，用假設的、具有相同導電性的單個孔隙組成的均質網絡來替代無序多孔介質的非均質網絡中的孔隙，同時保持孔隙網絡的連通性。Lian等[18]利用該方法研究了平均孔徑、孔徑分布和孔連通性對多孔炭電極離子輸運特性(包括電導和電導率)的影響。更進一步地，孔隙網絡模型(porenetworkmodeling,PNM)采用不同大小的球體和圓柱代表多孔結構中的空腔和吼道[19]。因此，PNM可以考慮多孔電極中的孔徑分布和孔道連通性等性質，且計算成本非常低，可以在合理的時間內通過PNP方程計算儲能設備的電化學性能和尋找最佳電極結構[20-21]。此外，Lian等[8]報道了一種簡單且基于物理的堆疊電極模型來表示多孔電極，成功用PNP方程和等效電路模型解釋了超級電容器的緩慢充電動力學，并發(fā)現(xiàn)充電過程可以分為兩個不同弛豫時間尺度的階段。雖然上述簡化的多孔電極模型被成功運用在解釋儲能過程中的各種現(xiàn)象，但這些方法只能在一定程度上接近真實的多孔電極的孔徑分布，很難表示多樣化的表面形態(tài)和分布復雜的催化活性位點，因此限制了這些簡化模型在多孔電極中的應用。

為了對真實的多孔電極(尤其是非均質、各向異性材料)中的熱-質傳遞過程進行深入的理解，需要對孔隙和骨架結構進行高分辨率三維表示。目前可以通過微尺度X射線計算機斷層掃描[22-23](computedtomography,CT)或掃描電子顯微鏡[24](scanningelectronmicroscope,SEM)以非破壞性的方式獲取多孔介質的孔隙結構。但是，CT不能分辨微米以下的結構，SEM的缺點是它只提供二維信息。獲得多孔電極的三維結構表示的一種流行的替代方法是隨機重構?？偟貋碚f，隨機重構包括處理可用的信息，并生成符合所導出的介質特性的隨機結構兩個步驟?；旧嫌袃煞N重構方法。第一種依賴于基于流程的建模。這種方法試圖模擬自然材料在其原始環(huán)境中的形成過程[25]，或創(chuàng)造工程材料的過程[26]。然而，現(xiàn)有的基于過程的方法只考慮了某些方面，難以準確描述材料結構。一種試圖解決這些問題的方法是相場方法[27]。第二種類型的重構方法采用了一種更加數據驅動的方法，并使用了統(tǒng)計數據。這些統(tǒng)計數據主要來自通過圖像分割識別出孔隙相和骨架相的各類圖片，使用兩點統(tǒng)計信息，如自相關或概率函數[28]、線性路徑函數[29]和聚類函數[26]。最常用的兩點統(tǒng)計量重構方法是模擬退火(simulatedannealing,SA)算法[30-32]。與可能非常復雜，并且可能需要大量超參數調優(yōu)的基于多點統(tǒng)計[33]的方法相反，模擬退火算法的簡單性和準確性在實際運用中很有吸引力。Wu等[34]采用SA算法重構商業(yè)鋰離子電池陰極LiCoO2的三維微觀結構，包括活性材料相、孔隙相和添加劑相。重構陰極的表征提供了重要的結構和輸運性質，包括兩點相關函數、體積比表面積、曲折度和單個相的幾何連通性。Habte等[35-36]研究了正極材料的微觀結構形態(tài)對鋰離子電池性能的影響，首先通過模擬退火算法生成了球形電極的結構，然后計算出結構參數并代入阻抗譜經驗公式中，最后得到的阻抗值與實驗結果符合較好。目前研究大多都停留在多孔電極的重構、孔結構參數[34,37]和等效傳遞系數[38-39]的計算上，真實多孔結構的形貌和孔隙結構對離子傳遞和電極導熱的影響尚未被研究過。

本文提出一套研究介觀尺度下真實多孔電極中的熱-質傳遞的研究框架：首先采用一種基于改進的狀態(tài)更新的隨機重建方法和動態(tài)退火系數相結合的模擬退火算法，將圖像分割后的二維SEM圖重構為真實介觀尺度的三維多孔電極，重構生成的多孔電極的結構和真實多孔電極截面上的結構在統(tǒng)計學意義上是一致的，并在此基礎上建立真實多孔電極中的離子傳遞和電極導熱模型。

1研究思路

本文的研究思路如圖2所示。通過改進的模擬退火算法和有限元計算，研究了真實多孔電板中的離子傳遞和電極導熱現(xiàn)象。

圖2

圖2基于模擬退火算法重構的真實多孔電極中的離子傳遞和電極導熱

Fig.2Theiontransportandelectrodeheatconductioninrealporouselectrodebasedonsimulatedannealingalgorithm

2真實多孔電極的三維重構及其結構參數計算

如圖3所示，多孔電極三維重構過程如下：首先利用SEM掃描多孔電極，得到參考圖像；然后對其進行預處理，即通過圖像分割識別孔隙相和骨架相，得到二值化圖像序列，并將二值化圖像序列作為參考模型；然后生成與參考模型的孔隙率相同的隨機初始點云，利用改進的模擬退火算法進行隨機重構，將兩點相關函數作為優(yōu)化目標，采用一種改進的體素交換方法更新狀態(tài)，得到重構的多孔電極的點云；最后將點云實體化為立方體構成的空隙相和骨架相的結構，并計算了其結構參數。

圖3

圖3多孔電極三維重構流程圖

Fig.3Theflowchartofthree-dimensionalreconstructionofporousmedia

2.1圖像預處理

目前研究和應用較多的電極材料主要集中在多孔炭材料、過渡金屬氧化物和導電高分子等方面。其中多孔炭材料具有比表面積大、導電性好、結構穩(wěn)定、資源豐富、價格低廉等優(yōu)勢，既可直接作為電極，又可與其他材料復合，起到傳輸電子、舒緩體積膨脹以及優(yōu)化界面反應等作用，被廣泛用為電極材料[40]。多孔炭材料的SEM圖來自于Kunanusont等[24]的工作，如圖3(a)所示。他們采用超臨界二氧化碳干燥法，在10.0~20.0MPa和40℃條件下制備了乙炔炭黑和聚偏氟乙烯基多孔電極，并研究了壓力對電極形態(tài)、孔隙率和電化學性能的影響。

通過SEM掃描得到的灰度圖像序列中，孔隙相和骨架相的灰度值差異明顯。因此，首先通過閾值分割算法將灰度圖像序列轉化為只有孔隙相(黑色區(qū)域，灰度值為0)和骨架相(白色區(qū)域，灰度值為255)的二值化圖像序列。采用SEM圖中所有像素的灰度值的平均值(本文為144)作為閾值。當某處的灰度值高于閾值時，則將此處的灰度值設為255；當灰度值小于等于平均值時，則設灰度值為0。處理后的圖像中含有斑點噪聲和大面積的孔洞，采用形態(tài)學運算填充大面積孔洞，同時去除圖像中的雜質和亮斑等噪聲干擾[41]，從而得到預處理之后的二值化圖像序列，如圖3(b)所示。

2.2圖像結構描述符計算

多孔電極的結構與其結構描述符密切相關，在對多孔電極重構時，需要考慮其結構描述符等信息，使重構模型與參考模型更加符合[42]。采用特征函數Irj定義不同的結構描述符。在多相系統(tǒng)中，假設第j相所占的區(qū)域為Vj，其在系統(tǒng)中的濃度為φj，則第j相的特征函數為

Irj=1r∈Vj0r?Vj(1)

式中，r是系統(tǒng)中任意一點。因為本文中的多孔電極僅由孔隙相和骨架相兩相構成，所以令j=1表示骨架相，j=2為孔隙相，則多孔電極的孔隙率為φ=1-φ1=φ2=I2rˉ。

本文采用的第一個結構描述符是兩點概率函數S2jr1,r2，簡寫為S2，它給出了位置r1和r2的兩個像素屬于同一相位的概率[30]。S2jr1,r2定義為

S2jr1,r2=Ir1jIr2j(2)

在平衡態(tài)系統(tǒng)中，粒子之間由于相互作用而存在一定的相關性。兩點概率函數作為一種常用的結構描述符，可以很好地描述系統(tǒng)的微觀結構特征。第二個描述符是兩點線性路徑函數L2jr1,r2，簡寫為L2。這個函數計算從r1到r2的線段完全位于第j相的概率[29]。值得注意的是，S2和L2僅在執(zhí)行重構時使用。而為了評估重構質量，本文使用兩點聚類函數C2jr1,r2，簡寫為C2，它計算在位置r1和r2上的兩個像素屬于同一簇的概率。簇被定義為屬于同一相并相互連接的像素的集合。

2.3模擬退火算法模擬退火算法是一種通用的概率優(yōu)化算法，起源于金屬的退火過程，用于在一個很大的搜尋空間中尋找出最優(yōu)解，可以求解復雜的非線性優(yōu)化問題。模擬退火算法其實是一種貪心算法，每次都選擇一個當前最優(yōu)解，因此只能搜索到局部的最優(yōu)值。但是它的搜索過程引入了隨機因素。模擬退火算法以一定的概率來接受一個比當前解要差的解，因此有可能會跳出局部的最優(yōu)解，達到全局的最優(yōu)解。

2.3.1模擬退火算法基本框架圖4是由Yeong等[31]首次提出的模擬退火算法運用在多孔介質重構中的基本框架。首先獲得參考模型w0[圖3(b)]的能量函數Ew0,設置初始溫度為T=T0。然后隨機產生與w0的孔隙率相同的初始結構w，并計算其能量函數Ew。通過像素交換的方式產生新的結構w，計算其能量函數Ew以及能量增量ΔE=Ew-Ew。若ΔE0，則接受新結構，即w=w,Ew=Ew；否則以概率P=e-ΔE/T接受新結構。在當前溫度下新結構被多次拒絕后，通過退火過程降低溫度T，然后重新搜索最優(yōu)結構，直到滿足算法的終止條件。

圖4

圖4模擬退火算法流程圖

Fig.4Theflowchartofsimulatedannealingalgorithm

退出模擬退火算法的條件如下：①系統(tǒng)的能量函數小于閾值Eth；②相鄰兩步能量下降的絕對值小于ΔEth；③連續(xù)Ncon次結構無更新，即連續(xù)拒絕隨機產生的新結構Ncon次；④達到最大的迭代次數Niter。

2.3.2能量函數原則上，任何結構描述符或描述符的組合都可以用來計算能量函數E。一般來說，定義E時考慮的描述符越多，計算量就越大，且過多描述符的加入并不能顯著提高重建質量。因此，本文僅使用S2和L2的組合S2-L2來計算E[31-32]。設能量函數E為重構模型與參考模型在各個方向上S2-L2的差值平方和的范數，計算方法為

E=∑d∑rS2,djr-L2,djr-S2ref,djr-L2ref,djr22(3)

式中，d指模型維度。對于二維圖像，d包含4個維度，即沿x軸和y軸，以及兩個主對角線方向dxy和dyx；對于三維圖像，d包含9個維度，除了二維的4個方向，還有沿z軸，以及dxz、dzx、dyz和dzy方向。

2.3.3基于DPN值的像素交換模擬退火算法的最重要的過程是在保持多孔電極孔隙率?不變的前提下，通過一定規(guī)則產生新結構。將孔隙相和骨架相的像素進行位置交換，等效于產生新結構的過程。

利用傳統(tǒng)的隨機像素交換方法產生新結構時，圖像中所有像素的交換概率相同，這種方法容易破壞已經重構好的局部結構，從而增加了尋找最優(yōu)結構的時間。同時，圖像中孤立的像素點可以被視為噪點，使用傳統(tǒng)方法也難以去除。本文引入DPN(differentphaseneighbors)[43]的概念來改進傳統(tǒng)的隨機像素交換方法。像素的DPN值指在該像素的鄰域內與該像素不同相的像素個數。像素的DPN值反映了該像素在圖像中的孤立程度，DPN=0表明該像素鄰域內所有的像素都與其同相。對于二維圖像，選取8鄰域；對于三維圖像序列，選取26鄰域。圖5列舉了二維圖像中幾種不同的3×3大小的結構，并指出了中心像素的DPN值。

圖5

圖5中心像素在不同結構中的DPN值

Fig.5DPNvaluesofcenterpixelindifferentstructures

2.3.4動態(tài)退火系數退火過程就是體系溫度下降的過程。溫度越低，體系跳出當前解的概率就越低。溫度降低太快，容易陷入局部最優(yōu)；溫度降低太慢，算法運行效率太低。因此，為了平衡精度和速度，需要選取合適的退火方法。常用的退火方法為Tk=λTk-1，其中，k為迭代的次數，λ為退火系數。傳統(tǒng)模擬退火算法的退火系數是常數，無法考慮系統(tǒng)在不同重構階段中能量函數的變化，容易陷入非最優(yōu)結構或者增加尋優(yōu)時間。本文采用基于動態(tài)退火系數的退火過程，考慮模擬退火過程中系統(tǒng)能量的變化情況，退火系數的計算方法[44]為

λ=λminEmin/EavλminλmaxEmin/EavλmaxEmin/Eavotherwise(4)

式中，λmin和λmax分別是退火系數的最小值和最大值，用于對退火系數限幅；Emin和Eav是迭代過程中能量函數的最小值和平均值。

2.3.5參數設置、點云實體化和運行環(huán)境根據前文模擬退火算法的退出條件，設置具體參數如表1所示。

表1模擬退火算法退出迭代的參數設置

Table1Parametersettingofsimulatedannealingalgorithmtoexittheiteration

EthΔEthNconNiter1×10-61×10-85001×106新窗口打開|下載CSV

通過模擬退火算法產生的點云無法通過有限元方法模擬其內的離子傳遞和電極導熱過程，需要被實體化為由面構成的三維結構。本文采用Blender自帶的API(applicationprograminterface)，將點云實體化為由小立方組成的多孔電極，最后提取表面結構。圖3(c)、(d)分別展示了200×200×200的多孔電極的骨架相和孔隙相，其中一個小立方體的邊長為8.2nm。

本文的計算平臺為4xIntel(R)Xeon(R)Platinum9242CPUat2.30GHz和384G內存的服務器，模擬退火算法和實體化算法通過Python3.8編程實現(xiàn)。

2.4多孔電極的結構參數孔隙率φ是保證模擬退火算法隨機構造的多孔電極結構符合實際的前提之一。φ的計算方法是孔隙相的立方體數量Np與多孔結構所有立方體數量Nt的比值，即φ=Np/Nt。因為本文構建的多孔結構為實際結構的一小部分，所以構建的多孔結構的六個方向的面為截面而非表面。因此，本文比表面積SV的計算方法為除去截面面積的骨架相表面積與骨架相的體積之比。

多孔電極是分形物質，它的孔隙結構復雜度與其孔隙分形維數等特征參數密切相關。目前常用的計算多孔結構分形維數的方法有盒計數法[45]、豪斯道夫方法[46]和隨機游走法[47]。盒計數法由于其精度更高、運算速度更快而被廣泛用于計算分形維數。對于多孔電極的孔隙體積Vp，有冪律關系Vpε∝εd，其中，Vpε是在對應長度尺度ε下測得的孔隙體積，d為多孔介質孔隙相的分形維數。對于三維多孔電極,用尺度為ε的立方體測量孔隙體積，則孔隙分形維數計算方法為d=lnNpε/lnε，其中，Npε為在長度尺度ε下識別出孔隙的立方體數量。

在實際的多孔電極中，孔道的形態(tài)一般是彎曲的，離子的傳輸路徑總長度為l，一般遠大于電極的厚度H。曲折因子τ表示結構對擴散過程的影響程度，定義為τ=l/H。本文采用皮爾遜隨機行走模擬[48]來計算曲折因子。

所有的結構參數均采用Python3.8編程計算。表2展示了尺寸分別為50×50×50、100×100×100和200×200×200個立方體的多孔電極的結構參數。結果表明多孔電極的尺寸對結構參數影響不大，驗證了重構的多孔電極結構是穩(wěn)定且有效的。

表2通過模擬退火算法重構的多孔電極的結構參數

Table2Structuralparametersofporouselectrodereconstructedbysimulatedannealingalgorithm

多孔電極尺寸/

(個立方體)

邊長/μm結構參數孔隙率比表面積/m-1分形維數曲折因子50×50×500.410.5451.493×1073.282.36100×100×1000.820.5451.223×1073.192.01200×200×2001.640.5451.306×1073.112.13新窗口打開|下載CSV

3三維多孔電極中的傳質與傳熱原理

為了研究多孔電極中的傳質和傳熱現(xiàn)象，參考實際超級電容器的工作原理[8]，構建了如圖6所示的兩個模型。受限于計算機的計算資源，離子傳遞模型采用兩個100×100×100大小的立方體多孔電極(邊長約為0.82μm)的孔隙相與200×100×100大小的長方體體相相接。值得注意的是，兩個孔隙相互為鏡像關系，即右邊的多孔電極是左邊多孔電極經過x、y、z三個方向的鏡像處理后生成的。而電極導熱模型采用100×100×100個格子大小的立方體多孔電極的骨架相。

圖6

圖6真實多孔電極中的離子傳遞和電極導熱模型

Fig.6Themodelofiontransferandelectrodeheatconductioninrealporouselectrodes

3.1充電過程

剛開始時整個體系的離子密度分布是均勻的，電勢為零。當左右兩個多孔電極施加相反的電勢后，體相和孔隙相中的陰陽離子受到電場力的作用，分別朝著相反的方向定向移動到多孔電極表面，從而形成電流并最終形成雙電層(electricaldoublelayer,EDL)儲存能量。雙電層的厚度通過德拜長度來表示。

λD=ε0εrkBT/2e2ρbulk(5)

式中，ε0是真空介電常數，F(xiàn)·m-1;εr是相對介電常數；kB是Boltzmann常數，J·K-1；T是溫度，K;e是基本電荷，C;ρbulk是體相離子數密度,m-3。由于在這個過程中離子電遷移導致電解液中的離子濃度不均勻，體系存在濃度梯度，所以充電過程中離子的運動還受到擴散的影響。

本文考慮到生成的雙電層不會發(fā)生重疊以及求解的復雜度兩個因素，于是采用Poisson-Nernst-Planck(PNP)方程來描述上述復雜孔隙結構中的離子充電動力學。電勢?V和離子的數密度ρm-3是PNP方程要求解的變量，其控制方程形式為：

?ρ±?t+??J±=0(6)J±=-D±(?ρ±+z±eρ±kBT??)(7)ε0εr??(-??)=eq(8)

式中，±分別代表陰陽離子；t是時間,s;?=??xi+??yj+??zk，是哈密頓算子;J是離子通量，m-3·s-1;D是離子擴散系數，m2·s-1;z是離子的價態(tài);q=z+ρ++z-ρ-，是電荷密度。式(7)左邊的項為體系的離子通量，它由右邊第一項的擴散通量和第二項的電場中電遷移通量共同組成。

設由兩個多孔電極的孔隙相和體相組成的離子運動空間為Ω。剛開始時整個空間的離子數密度ρ是均勻分布的，電勢?為零，即

I.C.ρ±Ω,t0=NAc0

?Ω,t0=0?x?Ω,t0=0

式中，NA是阿伏伽德羅常數。

因為多孔炭電極是良好的導體，所以可以假設整個多孔炭電極是一個等勢體，即多孔電極和電解液的界面處的電勢相等。圖6中紅色和藍色區(qū)域代表電解液和多孔電極的界面以及電解液和集電極的界面，分別帶正電勢Vm和負電勢-Vm，其中Vm=kBT/e為熱電勢。由于重構出的多孔電極并不是周期性結構，不能采用周期性邊界條件，所以本文采用不可滲透且電勢為零的邊界條件，即灰色區(qū)域。因此。設帶正電的邊界為S+，帶負電的邊界為S-，其他邊界為S0，則

B.C.?S+,t=+Vm

?S-,t=-Vm?x?S+,t=0?x?S-,t=0JiS0,t=0

本文的研究重點是在T=298.15K時NaCl水溶液中的電解質。鈉離子(+)和氯離子(-)的價態(tài)分別為z+=1和z-=1。本文使用以下參數設置[49-50]：D±=1.6×10-9m2?s-1、εr=71、c0=10mol?m-3。因此，Vm=0.02526V，λD=2.87nm，τRC=2.94×10-6s。模擬時間為(10-3~10-1)τRC，其中τRC=H+LλD/D為表征平板電極充電完成的弛豫時間[10]，H=0.82μm和L=0.82μm分別代表電極的厚度和體相的半長。

本文基于COMSOLMultiphysics5.4軟件，采用有限元法求解上述含有復雜邊界條件的PNP方程，從而得到離子在多孔電極中的遷移行為。計算平臺為4xIntel(R)Xeon(R)Platinum9242CPUat2.30GHz和384G內存的服務器。

3.2導熱過程

在實際的充放電過程中，由于電流的生成和離子的重排等因素[15]，整個體系會產生熱量。體系溫度主要受到熱生成、電解液內部的熱傳導和熱對流，以及多孔電極內的熱傳導影響。而溫度又會影響離子的輸運參數和電極的導電性能，從而影響整個充電過程。充電、產熱和傳熱三種現(xiàn)象相互耦合，極大地增加了計算量。此外，本文研究體系的電極導熱的弛豫時間遠小于充電的弛豫時間(或產熱的弛豫時間)，且模擬時間較短，可以將充電現(xiàn)象和電極導熱現(xiàn)象解耦，從而分開研究。因此，為了簡化電極導熱模型，只考慮電極中的熱傳導。本文采用傅里葉定律描述熱傳導過程。

?T?t=kρfcp?2T=a?2T(9)

式中，k是熱導率，W·m-1·K-1;ρf是骨架質量密度，kg·m-3;cp是骨架相的比定壓熱容，J·kg-1·K-1;a是熱擴散系數，m2·s-1。

設骨架相空間為γ，骨架相的初始溫度為T0，即

I.C.Tγ,t0=T0(10)

熱量通過電解液和電極的界面從電解液傳入電極中，設電解液和電極的界面為ΛT，并假設ΛT上產生不隨時間和位置發(fā)生變化的溫升ΔT0。同2.1節(jié)，因為不能設置周期性邊界條件，所以假設其他面Λ0是熱絕緣邊界。因此邊界條件為

B.C.TΛT,t=T0+ΔT0(11)-n?-k?TΛ0,t=0(12)

式中，n是法線方向的單位向量。

Kunanusont等[24]只測量了多孔炭電極的電導率。因此，根據文獻[51]，設置多孔炭電極的k=0.81W?m-1?K-1。假設ρf=1100kg?m-3，cp=1100J?kg-1?K-1，則a=6.69×10-7m2?s-1，τa=1.01×10-6s，ΔT0=0.1K。模擬時間為(10-3~10-1)τa，其中τa=H2/a為表征厚度為H的平板電極在表面加溫后整個電極內部的溫度達到均一時的時間[51]。

電極導熱過程同樣基于COMSOLMultiphysics5.4軟件，采用有限元法求解。計算平臺為4xIntel(R)Xeon(R)Platinum9242CPUat2.30GHz和384G內存的服務器。

4三維多孔電極中的傳質與傳熱原理

本節(jié)首先驗證了改進后的模擬退火算法的高效性，然后模擬了如圖6所示的多孔電極模型中的離子傳遞和電極導熱現(xiàn)象。

4.1改進的模擬退火算法結果

為了驗證改進后的模擬退火算法重構的準確性，圖7展示了重構結構w的能量函數E(w)隨迭代次數的變化曲線。結果表明相較于隨機的像素交換策略和靜態(tài)的退火系數的經典模擬退火算法，采用DPN值的像素交換策略和采用動態(tài)的退火系數的改進模擬退火算法在相同的迭代次數下的重構精度更高(即E(w)更小)。但由于每次迭代需要計算每個像素點的DPN值，所以花費的總時間也更長。但是，由于改進的模擬退火算法在每次迭代時E(w)下降的值更大(即斜率更大)，所以單位時間內E(w)的減小量可能也更大（即重構效率更高）。

圖7

圖7重構結構的能量函數隨迭代次數的變化

Fig.7Theenergyfunctionofthereconstructedstructurevarieswiththenumberofiterations

4.2離子傳遞

陰陽離子在電勢的驅動下，分別朝著正極和負極移動。10-1τRC時整個超級電容器的電勢分布?/Vm如圖8(a)所示，正極部分的孔隙相為高電勢，負極部分的孔隙相為低電勢，體相的電勢則整體線性變化。圖8(b)、(c)分別展示了正極部分的孔隙相與體相的截面處(x=-L)以及此截面前一個小立方體邊長距離的截面處(x=-L+0.01H)的?/Vm。如邊界條件設置一樣，多孔電極的骨架相和體相的接觸面的電勢?始終為Vm，即圖8(b)中紅色區(qū)域?/Vm=1。此外，電勢在yz平面并不是均勻的，骨架相附近的電勢較高，且電勢隨著x的變化而存在較大的變化。

圖8

圖810-1τRC時不同位置的電勢分布?/Vm

Fig.8Thepotentialdistribution?/Vmatdifferentpositionat10-1τRC

圖9展示了陰陽離子濃度分布的歸一化平均值(c++c-)/(2c0)。因為模擬時間只有10-1τRC，所以陰陽離子主要朝與體相接觸的骨架相截面遷移，較少地往孔中遷移。此外，結果還表明，陰陽離子更傾向于從截面的邊緣往中心遷移，如圖9(b)所示，邊緣的(c++c-)/(2c0)大于中心的。而遷移的離子來自于體相對應的位置，如圖9(c)所示，邊緣的(c++c-)/(2c0)小于中心的。

圖9

圖910-1τRC時不同位置的濃度分布(c++c-)/(2c0)

Fig.9Theconcentrationdistribution(c++c-)/(2c0)atdifferentpositionat10-1τRC

4.3電極導熱

圖10為多孔電極中的歸一化溫度分布[(T-T0)]/ΔT0隨時間變化的等值面圖，包括在10-3τa、10-2τa和10-1τa時刻下。圖中黃色區(qū)域，即[(T-T0)]/ΔT0=0，代表低溫;白色區(qū)域，即(T-T0)/ΔT0=1，代表高溫。當(T-T0)/ΔT0=0時，表明熱量還沒傳導到這個位置；當(T-T0)/ΔT0=1時，說明此時溫度已經上升到能達到的最高溫度T0+ΔT0。結果表明：在10-1τa時，多孔電極的整個骨架相幾乎都是白色，表明從導熱邊界傳遞給骨架相的熱量已經完全傳遞到各個角落，骨架相的溫度已經均勻分布。此時的時間10-1τa遠小于平板電極的溫度弛豫時間τa，原因是熱量不僅僅是從骨架相與體相接觸的界面?zhèn)鲗牍羌芟嗟?，還包括孔道內部的壁面，從而減小了實際的導熱距離，最終使在多孔電極中的熱弛豫時間遠小于平板的弛豫時間τa。

圖10

圖10不同時刻下的溫度分布(T-T0)/ΔT0的等值面圖

Fig.10Theisosurfacediagramoftemperaturedistribution(T-T0)/ΔT0atdifferenttime

5結論

（1）一種改進的模擬退火算法重構出的多孔電極結構是穩(wěn)定且有效的。

（2）多孔電極內的電勢分布不均勻，體相中的電勢沿x方向整體線性變化。

（3