上個(gè)月給一個(gè)儲能項(xiàng)目做18650圓柱鋰電池的熱仿真要求在COMSOL 5.6里把電化學(xué)過程和產(chǎn)熱耦合起來對比0.5C、1C、2C三種放電倍率下的溫升表現(xiàn)。模型參數(shù)前后配置了一個(gè)多星期參數(shù)化掃描倒是跑得很順利最后得到的規(guī)律很直觀倍率從1C翻到2C最高溫度增幅遠(yuǎn)不是線性翻倍而是接近三四倍。這篇就把整個(gè)建模思路、參數(shù)表配置、參數(shù)化掃描設(shè)置和踩過的坑完整梳理一遍給準(zhǔn)備入手電化學(xué)-熱耦合仿真的同學(xué)一份可以直接參考的路徑。1. 18650不是普通圓柱體卷繞結(jié)構(gòu)給仿真挖了哪些坑1.1 18650的幾何特征與卷繞電極18650這個(gè)名字本身就暴露了規(guī)格直徑18mm、高度65mm0代表圓柱形。它最典型的特征是內(nèi)芯不是實(shí)心塊而是一層層“正極片-隔膜-負(fù)極片”疊好后卷繞成螺旋卷筒業(yè)內(nèi)叫jelly roll。這個(gè)卷繞結(jié)構(gòu)直接決定了熱仿真里最重要的一件事情——導(dǎo)熱系數(shù)是各向異性的。徑向方向和軸向方向的傳熱路徑完全不一樣。徑向的熱量要穿過幾十層薄電極和隔膜每一層界面接觸熱阻都在做貢獻(xiàn)等效徑向?qū)嵯禂?shù)通常只有0.31W/(m·K)左右軸向的熱量主要沿著銅箔、鋁箔流動等效軸向?qū)嵯禂?shù)能到2030W/(m·K)。如果仿真時(shí)圖省事把18650直接當(dāng)成各向同性的圓柱體來算溫度場結(jié)果會差很多尤其是2C、3C這種大倍率下中心熱堆積會被明顯低估。從建模角度說完整還原每層卷繞幾何在COMSOL里非常痛苦螺旋結(jié)構(gòu)網(wǎng)格數(shù)量大、收斂也難。常見做法有兩種第一種是二維軸對稱模型把卷繞體在橫截面上等效成若干層同心圓弧正極、隔膜、負(fù)極交替排列必要時(shí)通過折合厚度來匹配真實(shí)體積第二種是三維模型但把jelly roll整體簡化成等效圓柱徑向和軸向分別賦予不同的等效導(dǎo)熱系數(shù)。我個(gè)人做倍率對比研究時(shí)首選二維軸對稱因?yàn)閰?shù)化掃描要跑很多組二維軸對稱模型計(jì)算量小一大截而且電極間的濃度梯度和溫度梯度都能體現(xiàn)出來。1.2 電化學(xué)-熱耦合的三熱源模型電池產(chǎn)熱不是單一來源嚴(yán)格來說分三部分不可逆的極化熱、歐姆熱以及可逆的熵變熱。極化熱來自電極反應(yīng)的過電位放電時(shí)實(shí)際電位低于平衡電位這部分損失的能量直接變成熱數(shù)值上等于局部電流密度乘以活化過電位。歐姆熱來自離子在電解液中的傳導(dǎo)和電子在集流體、活性物質(zhì)顆粒中的傳導(dǎo)電阻它與電流平方成正比所以倍率升高時(shí)這一項(xiàng)增長非常快?？赡鏌醽碜噪姌O材料嵌鋰和脫鋰過程中的熵變化放電時(shí)一般是放熱充電時(shí)一般是吸熱它與電流一次方成正比低倍率時(shí)占比很高。COMSOL中如果用“鋰離子電池”接口配合“固體傳熱”接口電化學(xué)計(jì)算出的局部產(chǎn)熱會被當(dāng)作熱源傳遞給溫度場同時(shí)溫度又反饋給電解液電導(dǎo)率、反應(yīng)速率常數(shù)、平衡電位等參數(shù)形成完整的電化學(xué)-熱雙向耦合。這個(gè)耦合在機(jī)理上不復(fù)雜但實(shí)際建模時(shí)溫度反饋的開關(guān)經(jīng)常被漏掉做出來的溫度曲線一路不變那就是沒耦合上需要回多物理場節(jié)點(diǎn)里檢查。2. COMSOL 5.6模型搭建幾何簡化與參數(shù)配置清單2.1 幾何簡化方案二維軸對稱 vs 三維標(biāo)題明確是圓柱鋰電池那就要先把坐標(biāo)系統(tǒng)想清楚。標(biāo)準(zhǔn)18650的旋轉(zhuǎn)軸就是圓柱中心軸橫截面是一個(gè)圓環(huán)建立二維軸對稱模型利用旋轉(zhuǎn)對稱性既能把電極厚度方向的變化真實(shí)還原計(jì)算量又比三維模型小一個(gè)量級。二維軸對稱幾何一般把jelly roll簡化成有限的幾層同心圓域正極、隔膜、負(fù)極、正極……再加上銅箔和鋁箔集流體。簡化時(shí)要注意總厚度與真實(shí)卷繞厚度的匹配尤其隔膜厚度不能拍腦袋因?yàn)楦裟︿囯x子傳輸阻力影響巨大。如果目標(biāo)是研究電池殼表面溫度分布以及端部極耳散熱那就要考慮三維模型通常把jelly roll當(dāng)成一個(gè)各向異性等效塊體正極、負(fù)極、隔膜不再逐層展開而是通過等效參數(shù)把電化學(xué)產(chǎn)熱映射到宏觀塊體上。它的好處是能加外殼、極耳、絕緣墊片壞處是電化學(xué)場的局部信息丟了一部分。以三種放電倍率對比這種研究需求來說我強(qiáng)烈建議先用二維軸對稱模型跑通確認(rèn)規(guī)律后再升級三維。二維軸對稱模型參數(shù)調(diào)整方便參數(shù)化掃描速度快踩完坑再上三維在工程上是最省時(shí)間的路徑。2.2 物理場接口與耦合方式COMSOL 5.6里需要用到三個(gè)接口鋰離子電池li-ion、固體傳熱ht以及多物理場耦合節(jié)點(diǎn)。鋰離子電池接口負(fù)責(zé)求解電極固相鋰濃度、電解液鋰濃度、固相電勢、電解液電勢和局部電流密度。需要分別定義正極和負(fù)極的活性材料包括初始鋰濃度、最大嵌入濃度、反應(yīng)速率常數(shù)、傳遞系數(shù)以及OCV隨嵌鋰分?jǐn)?shù)的變化曲線。這些參數(shù)是電化學(xué)仿真的核心也是最花時(shí)間的地方。固體傳熱接口負(fù)責(zé)求解溫度場需要給定每個(gè)域的密度、比熱容、導(dǎo)熱系數(shù)。這里特別提醒二維軸對稱模型里jelly roll域的徑向和軸向?qū)嵯禂?shù)要分開放不能填同一個(gè)值。耦合方式在COMSOL 5.6里比較直接多物理場節(jié)點(diǎn)選擇“電化學(xué)-熱耦合”它會自動把鋰離子電池計(jì)算出的熱源傳遞給固體傳熱同時(shí)讓溫度影響電化學(xué)參數(shù)。需要注意的是耦合節(jié)點(diǎn)創(chuàng)建后電解液電導(dǎo)率的溫度依賴通常需要手動補(bǔ)充表達(dá)式COMSOL自帶的材料庫未必覆蓋。我一般用Arrhenius形式補(bǔ)上電導(dǎo)率隨溫度變化的表達(dá)式這部分在后面第五章會具體講。2.3 參數(shù)表最容易出錯(cuò)的單位換算參數(shù)配置是這個(gè)模型里最容易被軟坑絆倒的地方。我的習(xí)慣是把所有參數(shù)放進(jìn)全局參數(shù)表用符號管理不要散落在各物理場里。下表是當(dāng)時(shí)用的核心參數(shù)例子具體數(shù)值可以根據(jù)你手里的電芯調(diào)整參數(shù)值說明C_cell2.5 Ah電池額定容量V_nom3.6 V標(biāo)稱電壓C_rate1放電倍率參數(shù)掃描變量A_geo0.0185 m2電極投影面積單位換算用sigma_e_ref1.1 S/m電解液參考電導(dǎo)率Ea_sigma15000 J/mol電導(dǎo)率活化能x_init_neg0.85負(fù)極初始嵌鋰分?jǐn)?shù)y_init_pos0.25正極初始脫鋰分?jǐn)?shù)最容易出錯(cuò)的是電流密度和容量單位。1Ah等于3600庫侖如果你在邊界條件里施加電流單位是A/m2就必須先算清楚電極的投影面積電流密度表達(dá)式是C_cell乘以C_rate再除以A_geo。算錯(cuò)一個(gè)小數(shù)點(diǎn)溫度曲線整條都會偏。另一個(gè)坑是鋰濃度與SOC的關(guān)系。COMSOL里初始濃度由嵌鋰分?jǐn)?shù)和最大鋰濃度決定負(fù)極初始嵌鋰分?jǐn)?shù)、正極初始脫鋰分?jǐn)?shù)一定要對應(yīng)同一個(gè)SOC否則初始開路電壓從第一步就漂移后面所有倍率的結(jié)果都沒意義。搞反初值是很多新手一上來發(fā)現(xiàn)電壓曲線不對的頭號原因。3. 三種放電倍率的參數(shù)化掃描從電流表達(dá)式到求解設(shè)置3.1 倍率參數(shù)與電流負(fù)載表達(dá)參數(shù)化掃描的前提是讓模型里的關(guān)鍵量用參數(shù)表達(dá)。我定義了全局參數(shù)c_rate默認(rèn)值1然后在電流負(fù)載邊界上寫表達(dá)式I_disch c_rate * C_cell / A_geo這里的A_geo需要根據(jù)幾何模型來定。二維軸對稱模型里電流施加在電極端面邊界上A_geo就是該邊界的面積三維模型如果多個(gè)極耳要把總電流拆到各個(gè)極耳面上。還要注意符號約定。放電電流方向和充電電流方向不同COMSOL里習(xí)慣把放電定義為正電流還是負(fù)電流取決于你建的是正極域還是完整電池模型。建議先跑一個(gè)0.5C的短時(shí)間測試看電壓曲線是下降還是上升電壓下降且SOC減小說明方向?qū)α巳绻妷荷仙央娏鞣柗催^來就行。3.2 參數(shù)化掃描與瞬態(tài)求解的配合COMSOL 5.6中在“研究1”下面添加“參數(shù)化掃描”節(jié)點(diǎn)選擇全局參數(shù)c_rate值列表寫0.5、1、2。這樣同一個(gè)瞬態(tài)研究會被依次執(zhí)行三次每次都使用不同的c_rate值結(jié)果會自動按參數(shù)值區(qū)分。這里有一個(gè)值得注意的點(diǎn)參數(shù)化掃描是對整個(gè)研究循環(huán)執(zhí)行所以掃描列表里的每個(gè)值都會獨(dú)立完成從初始條件開始的瞬態(tài)計(jì)算。如果模型網(wǎng)格密三次非線性瞬態(tài)會非常耗時(shí)。我的經(jīng)驗(yàn)是如果初始SOC相同可以在掃描設(shè)置里開啟“重用前一個(gè)參數(shù)值的解”作為初值能節(jié)省不少冷啟動時(shí)間。但要注意這并不總是合適如果參數(shù)變化導(dǎo)致物理過程差異很大復(fù)用舊解可能帶入誤差。求解器方面瞬態(tài)用BDF方法如果發(fā)現(xiàn)1C能收斂、2C不收斂通常要限制最大時(shí)間步長比如設(shè)置成10s以內(nèi)并且打開“初始步長控制”。理論上參數(shù)化掃描不會額外增加求解器難度難點(diǎn)都在單次瞬態(tài)的收斂上。3.3 不同倍率的時(shí)間窗口與截止條件三種倍率放電到同一截止電壓的時(shí)間完全不同。0.5C大概是120分鐘1C大概是60分鐘2C大概30分鐘具體數(shù)值取決于電芯容量和化學(xué)體系但量級就是這樣。如果參數(shù)化掃描里所有倍率都用同一個(gè)時(shí)間范圍[0, 3600s]0.5C還在放電中期2C可能已經(jīng)跑完并且電壓跌破截止值。更好的做法是把時(shí)間也參數(shù)化定義t_max 3600 / abs(c_rate)這里3600對應(yīng)1C放電一小時(shí)然后求解器的時(shí)間范圍寫成range(0, 50, t_max)。這樣每個(gè)倍率有各自的結(jié)束時(shí)間積分終點(diǎn)大致對齊到滿放條件。固定時(shí)間范圍無法精確停在截止電壓。如果需要嚴(yán)格按2.5V截止來判斷就要用事件接口或者停止條件比如檢測全局電壓低于閾值時(shí)終止。但事件會讓參數(shù)化掃描和求解配置復(fù)雜不少。對于多數(shù)倍率對比研究先算到固定SOC終點(diǎn)再對比溫度曲線已經(jīng)足夠事件屬于錦上添花的高級選項(xiàng)。4. 結(jié)果解讀0.5C、1C、2C下生熱行為的差異4.1 溫升曲線對比倍率與產(chǎn)熱的非線性參數(shù)化掃描跑完后在結(jié)果里可以繪制全局溫度最大值隨時(shí)間的曲線。我習(xí)慣用點(diǎn)探針監(jiān)測jelly roll中心點(diǎn)溫度這比只看邊界溫度更有代表性因?yàn)殡姵刂行牡臏囟韧ǔＷ罡?。模型?shí)測下來0.5C放電結(jié)束時(shí)中心溫升大約58K1C大約1218K2C能到3545K。這些數(shù)字會隨電芯容量和化學(xué)體系變化不需要照搬但規(guī)律性很明顯2C溫升不是1C的兩倍而是接近三倍。原因就是歐姆熱和極化熱都強(qiáng)烈依賴電流電流從1C翻到2C時(shí)各項(xiàng)熱源疊加放大而對流散熱只有一個(gè)與溫差線性相關(guān)的通道熱量累積速度超過了散熱的線性增長。這也是為什么電池倍率性能的瓶頸往往不在電化學(xué)而在熱大電流下熱量堆積嚴(yán)重如果散熱設(shè)計(jì)跟不上電芯內(nèi)部溫度會迅速逼近材料的安全上限。做電化學(xué)生熱研究第一眼要看的圖就是這個(gè)溫升對比圖。4.2 熱源分解高倍率下誰在主導(dǎo)發(fā)熱COMSOL 5.6里可以分別提取極化熱、歐姆熱、可逆熵變熱三項(xiàng)熱源。通過對整個(gè)jelly roll域做體積積分可以定量比較每一項(xiàng)的貢獻(xiàn)。我觀察到的典型趨勢是0.5C時(shí)可逆熱占比很高有的區(qū)間能到40%以上因?yàn)殡娏餍?、過電位和歐姆壓降都小可逆熱作為一次方項(xiàng)相對突出到了2C極化熱和歐姆熱占比加起來超過80%可逆熱的絕對數(shù)值雖然也變大但相對占比被稀釋。這個(gè)拆解很有用。判斷電池發(fā)熱瓶頸時(shí)不能只看總產(chǎn)熱還要看熱源位置和隨倍率的變化。比如負(fù)極歐姆熱在高倍率下集中在近隔膜一側(cè)電解液電導(dǎo)率不足時(shí)局部產(chǎn)熱會畸形升高。拆項(xiàng)還有一個(gè)額外好處如果某個(gè)倍率下可逆熱出現(xiàn)負(fù)值且數(shù)值很大說明熵變參數(shù)符號填錯(cuò)了或者OCV隨SOC求導(dǎo)方向有問題這類問題單看總產(chǎn)熱曲線很難察覺。4.3 溫度場分布中心熱堆積與邊界散熱隨著放電倍率升高徑向溫度梯度會明顯增大。中心區(qū)域的熱量要向外傳遞必須穿過卷繞層溫度最高外表面通過對流散熱溫度最低。軸向方向兩個(gè)端面暴露在環(huán)境中散熱相對快極耳區(qū)域是金屬集流體局部導(dǎo)熱強(qiáng)會形成低溫帶。這類結(jié)果最直接的用途就是判斷電芯內(nèi)部是否存在熱點(diǎn)。如果2C時(shí)中心溫度比表面高810K說明電芯內(nèi)部散熱結(jié)構(gòu)設(shè)計(jì)不夠BMS里表面溫度傳感器的讀數(shù)不能直接代表內(nèi)部最大溫度。這也是18650做電化學(xué)生熱研究最實(shí)際的價(jià)值之一。畫溫度場時(shí)我會同時(shí)輸出徑向溫度分布線圖和各時(shí)刻等溫線圖。等溫線如果沿著中軸方向被拉長說明軸向?qū)嵴純?yōu)如果看到等溫線大量聚集在近表面區(qū)域說明表面散熱對徑向梯度影響大。參數(shù)化掃描后這些圖會自動按倍率分組對比非常直觀。5. 實(shí)操避坑與收斂性調(diào)優(yōu)5.1 網(wǎng)格劃分策略映射網(wǎng)格和邊界層對于二維軸對稱模型網(wǎng)格優(yōu)先用映射網(wǎng)格。正極、隔膜、負(fù)極和集流體都是規(guī)則的環(huán)狀域映射網(wǎng)格規(guī)整可控。沿電極厚度方向一定要多分幾層我一般至少安排5到10個(gè)單元因?yàn)榇蟊堵史烹姇r(shí)鋰離子濃度梯度集中在電極厚度方向網(wǎng)格太粗會出現(xiàn)濃度振蕩。隔膜域的網(wǎng)格可以比電極稍粗但要保證與兩側(cè)電極的界面節(jié)點(diǎn)對齊。界面不對齊會造成插值誤差局部電流密度分布出現(xiàn)鋸齒狀進(jìn)而影響產(chǎn)熱分布的平滑性。大倍率下電解液濃度在電極和隔膜界面的邊界層效應(yīng)明顯可以加邊界層網(wǎng)格。COMSOL的邊界層屬性直接加在域邊界上厚度和層數(shù)可以控制。如果只做0.5C和1C對比邊界層網(wǎng)格不是必須的但到2C、3C就值得加。5.2 大倍率不收斂的常見原因與對策遇到2C不收斂先不要懷疑網(wǎng)格。最常見的根本原因是參數(shù)空間本身已經(jīng)讓局部電流逼近極限電流密度電解液鹽濃度在局部趨近于零數(shù)值上表現(xiàn)為濃度負(fù)值、NaN或者電壓劇烈抖動。我的排查順序是這樣的先看端電壓曲線是否在某個(gè)時(shí)間點(diǎn)突變再看是否只有高倍率有問題。如果是局部濃度問題調(diào)大電解液體積分?jǐn)?shù)、提高初始鹽濃度只能算壓參數(shù)不能解決物理模型本身的限制。更合理的做法是檢查活性材料的有效擴(kuò)散系數(shù)是否低估或者給電解液補(bǔ)充濃度相關(guān)的電導(dǎo)率模型而不是常數(shù)。時(shí)間步長也常背鍋。BDF求解器的最大時(shí)間步長如果全自動瞬態(tài)過程中電流變化劇烈時(shí)容易跳過物理過程。手動限制最大步長到25s有時(shí)就解決問題但代價(jià)是求解時(shí)間明顯變長所以先用1C跑通過再開2C。最后還有初始條件一致性初始化的問題。鋰離子電池接口的初始濃度場如果不做初始兼容處理第一個(gè)時(shí)間步里就會產(chǎn)生大量虛假產(chǎn)熱有時(shí)表現(xiàn)為起始瞬間溫度跳變。打開一致性初始化或者對初始值先做穩(wěn)態(tài)計(jì)算能消除這類偽源。這個(gè)坑很隱蔽排查順序建議放在電壓突變之后。5.3 材料物性溫度依賴與結(jié)果校核電化學(xué)-熱耦合仿真里最關(guān)鍵的參數(shù)其實(shí)是電解液電導(dǎo)率隨溫度的變化。如果沒有溫度依賴溫度場再高也不會反饋給電化學(xué)所謂耦合只是單向的電化學(xué)產(chǎn)熱加傳熱2C下的歐姆熱會被低估因?yàn)閷?shí)際溫度升高后電解液電導(dǎo)率上升歐姆電阻反而下降。建立溫度依賴表達(dá)式后要注意一個(gè)反直覺現(xiàn)象熱耦合導(dǎo)致電導(dǎo)率上升高倍率后期的歐姆熱可能不再單調(diào)上升溫升曲線出現(xiàn)先急后緩的轉(zhuǎn)折。這符合物理但首次看到會以為是數(shù)值誤差。結(jié)果校核方面最有效的是對照電壓曲線。模擬得到的端電壓-時(shí)間曲線與實(shí)測放電曲線的偏差如果控制在5%以內(nèi)溫度趨勢基本可信。如果沒有實(shí)驗(yàn)數(shù)據(jù)可以用能量守恒做粗校核總產(chǎn)熱積分減去表面散熱量等于電池系統(tǒng)內(nèi)能增量近似等于m乘以Cp乘以ΔT。如果這個(gè)平衡差得很遠(yuǎn)多半是熱源項(xiàng)或者邊界條件標(biāo)錯(cuò)了。6. 模型擴(kuò)展與實(shí)際應(yīng)用的個(gè)人經(jīng)驗(yàn)三種放電倍率只是第一層參數(shù)化。模型跑通后把SOC或者初始充電狀態(tài)也作為第二個(gè)掃描參數(shù)就能得到倍率-SOC二維矩陣。這對研究低SOC大倍率下的溫升特別有意義因?yàn)榈蚐OC時(shí)電池內(nèi)阻增大同樣倍率產(chǎn)熱更高。COMSOL的參數(shù)化掃描支持多參數(shù)嵌套先掃SOC再掃倍率結(jié)果組織也很方便。如果后續(xù)要支撐實(shí)時(shí)BMS熱管理策略可以把高保真仿真結(jié)果降階成集總參數(shù)模型。具體做法是先跑幾組典型倍率和典型環(huán)境溫度下的仿真提取時(shí)間常數(shù)和溫升系數(shù)再擬合為狀態(tài)空間方程。這個(gè)過程中COMSOL的探針表批量導(dǎo)出數(shù)據(jù)特別好用。我的個(gè)人體會是仿真項(xiàng)目最花時(shí)間的不是模型搭建而是參數(shù)校準(zhǔn)和對實(shí)驗(yàn)數(shù)據(jù)的擬合。COMSOL里參數(shù)配置再漂亮如果電芯材料參數(shù)來自不同文獻(xiàn)、不同批次結(jié)果偏差很容易超過20%。所以每個(gè)參數(shù)都得能說清出處尤其是OCV-SOC曲線、熵變系數(shù)、電解液電導(dǎo)率這三項(xiàng)。最后分享一個(gè)小技巧參數(shù)化掃描結(jié)果別急著合并曲線先單獨(dú)檢查各倍率的電壓曲線起始段是否重疊。同一初始SOC下三種倍率初始開路電壓應(yīng)當(dāng)幾乎相同如果有較大偏差說明初始鋰濃度設(shè)置或SOC換算有問題這個(gè)誤差會在0.5C、1C、2C三組結(jié)果里同步放大。先把這個(gè)校核點(diǎn)過了再去看溫度結(jié)果能省下很多反復(fù)排查的時(shí)間。