Modelo de Potencial Químico Aprimorado para Deformação Rápida em Eletrodos de Nanopartículas de Silício de Baterias Li-Ion

Resumo

Apresentamos uma expressão de potencial químico refinada que remove a suposição de baixa taxa de deformação, permitindo modelagem precisa de tensão induzida por difusão e distribuição de lítio em eletrodos de nanopartículas de silício. Comparando a nova expressão com a tradicional sob ciclagem potenciostática e galavanostática, demonstramos divergência pronunciada durante operação potenciostática, enquanto resultados galavanostáticos permanecem amplamente inalterados. Um raio crítico — onde fluxo de Li e tensão hidrostática são essencialmente idênticos para ambos modelos — é identificado pela primeira vez.

Introdução

Eletrodos de alta capacidade como silício, que sofrem mudanças volumétricas de até 400% durante litiação, são fundamentais para eletrônicos portáteis, veículos elétricos, e armazenamento em escala de rede [1-3]. No entanto, a tensão induzida por difusão gerada por expansão volumétrica homogênea pode desencadear fratura frágil, encurtando vida útil e reduzindo capacidade [4]. Eletrodos compostos exibem morfologias complexas que desafiam descrições puramente teóricas; a maioria dos modelos trata propriedades do material como espacialmente variáveis enquanto negligencia efeitos interfaciais. Geometrias típicas de eletrodo consideradas na teoria são esférica, cilíndrica, e planar, com modelos unidimensionais ou bidimensionais de acordo.

Investigações recentes em eletrodos de nanopartículas de silício destacaram a interação entre transformação de fase induzida por litiação, evolução morfológica, e geração de tensão [5-7]. Estes estudos baseiam-se no equilíbrio termodinâmico de sólidos, onde difusão atômica é influenciada por tensão, e vice-versa. Estrutura de Larche e Cahn [8] introduziu potencial químico dependente de tensão assumindo deformações pequenas e isotrópicas, enquanto formulação de Wu [9] incorporou tensor momento Eshelby. Cui et al. [10] estenderam para deformações finitas, ainda sua derivação permanece precisa apenas para baixas taxas de deformação relativas à difusão. Dada a grande expansão volumétrica do silício durante carregamento rápido, esta suposição pode levar a erros significativos.

Neste trabalho, derivamos expressão de potencial químico válida para deformação rápida sem suposição de baixa taxa, independente de morfologia do eletrodo. Aplicando o novo modelo a eletrodos de nanopartículas de silício sob controle potenciostático e galavanostático, quantificamos impacto na tensão e distribuição de lítio. Também identificamos raio crítico onde os dois modelos convergem, fornecendo insight sobre limites espaciais da influência de tensão.

Métodos

Equações Mecânicas

Inserção de lítio induz mudança volumétrica. Descrevemos deformação usando coordenadas Lagrangeanas e Eulerianas:
$$\mathbf{U}=\mathbf{x}-\mathbf{X}\quad (1)$$ onde U é campo deslocamento. Tensor gradiente de deformação é
$$\mathbf{F}=\frac{\partial \mathbf{x}}{\partial \mathbf{X}}=\mathbf{I}+\mathrm{Grad}\,\mathbf{U}\quad (2)$$

Para partícula esférica com coordenadas (R,Θ,Φ) (Lagrangeana) e (r,θ,φ) (Euleriana), tensor gradiente simplifica para
$$\mathbf{F}=\begin{bmatrix}{F}_R&0&0\0&{F}_{\Theta}&0\0&0&{F}_{\Phi}\end{bmatrix}=\begin{bmatrix}1+\partial u/\partial R&0&0\0&1+u/R&0\0&0&1+u/R\end{bmatrix}\quad (3)$$

Deformação é decomposta em mudança de forma induzida por litiação e deformação elástica reversível:
$$\mathbf{F}=\mathbf{F}^e\mathbf{F}^c\quad (4)$$ com F^c isotrópico:
$$\mathbf{F}^c=(1+\Omega C)^{1/3}\mathbf{I}\quad (5)$$ onde Ω é volume molar parcial. Consequentemente,
$$\mathbf{F}^e=(1+\Omega(R)C)^{-1/3}\begin{bmatrix}1+\partial u/\partial R&0&0\0&1+u/R&0\0&0&1+u/R\end{bmatrix}\quad (6)$$

Tensor deformação Green-Lagrange E é
$$\mathbf{E}=\frac{1}{2}\left(\mathbf{F}^T\mathbf{F}-\mathbf{I}\right)\quad (7)$$ resultando componentes radial e tangencial
$$E_R^e=\tfrac{1}{2}\Big[\tfrac{(1+\partial u/\partial R)^2}{(1+\Omega(R)C)^{2/3}}-1\Big]\quad (9)$$
$$E_{\Theta}^e=E_{\Phi}^e=\tfrac{1}{2}\Big[\tfrac{(1+u/R)^2}{(1+\Omega(R)C)^{2/3}}-1\Big]\quad (10)$$

Assumindo elasticidade linear, densidade energia de deformação é
$$W=\tfrac{1}{2}\mathbf{E}^e:\mathbf{C}:\mathbf{E}^e=\det(\mathbf{F}^c)\frac{E_h}{2(1+\nu)}\Big[\tfrac{\nu}{1-2\nu}\big(tr\,\mathbf{E}^e\big)^2+tr\,\mathbf{E}^e\mathbf{E}^e\Big]\quad (12)$$ onde E_h e ν são módulo Young e razão Poisson. Tensão primeira Piola-Kirchhoff segue
$$\mathbf{P}=\det(\mathbf{F}^c)\frac{E_h}{1+\nu}\Big[\tfrac{\nu}{1-2\nu}tr\,\mathbf{E}^e+\mathbf{E}^e\Big]\mathbf{F}^e(\mathbf{F}^c)^{-1}\quad (13)$$

Componentes:
$$P_R=\left(1+\Omega C\right)^{1/3}\frac{E_h}{2(1+\nu)(1-2\nu)}\left(1+\tfrac{\partial u}{\partial R}\right)\Big[(1-\nu)E_R^e+2\nu E_{\Theta}^e\Big]\quad (14)$$
$$P_{\Theta}=P_{\Phi}=\left(1+\Omega C\right)^{1/3}\frac{E_h}{2(1+\nu)(1-2\nu)}\left(1+\tfrac{u}{R}\right)\Big[\nu E_R^e+E_{\Theta}^e\Big]\quad (15)$$

Equilíbrio sem forças de corpo requer
$$\frac{\partial P_R}{\partial R}+\frac{2(P_R-P_{\Phi})}{R}=0\quad (16)$$ com condições contorno
$$u(0,t)=0,\;P_R(R_0,t)=0\quad (17)$$

Equação Transporte de Massa

Concentração lítio C(X,t) e fluxo J(X,t) satisfazem
$$\frac{\partial C}{\partial t}+\mathrm{Div}\,\mathbf{J}=0\quad (18)$$. Sob simetria esférica, reduz a
$$\frac{\partial C(R,t)}{\partial t}+\frac{1}{R^2}\frac{\partial\big(R^2J(R,t)\big)}{\partial R}=0\quad (19)$$. Fluxo segue lei Fick com força motriz potencial químico:
$$J=-\frac{CD}{R_gTF_{11}F_{11}}\frac{\partial \mu}{\partial R}\quad (20)$$

Potencial químico deriva de densidade energia interna total
$$\Pi=\varphi(C)+W(C,\mathbf{E}^e)\quad (22)$$, dando
$$\mu=\frac{\partial\Pi}{\partial C}=\mu_0(C)+\tau(\mathbf{E}^e,C)\quad (23)$$ onde τ é parte dependente de tensão. Para solução diluída,
$$\mu_0(C)=\mu_0+R_gT\ln(\gamma C)\quad (24)$$ com γ=1.

Tradicionalmente, abordagem energia livre Helmholtz resulta em
$$\tau_H=\frac{\partial W}{\partial C}= -\det(\mathbf{F}^e)\sigma_m\Omega\quad (25)$$ levando a
$$\mu={\mu}_0+R_gT\ln C-\det(\mathbf{F}^e)\Omega\sigma_m\quad (26)$$. Aqui, σ_m é tensão hidrostática Cauchy
$$\sigma_m=\tfrac{1}{3}tr\,(\det^{-1}\mathbf{F}\mathbf{PF}^T)\quad (27)$$. Esta expressão (26) é o modelo tradicional de potencial químico.

Alternativamente, tratando energia livre Gibbs resulta em
$$\tau_G=\frac{\partial W}{\partial C}\quad (28)$$ e consequentemente
$$\mu={\mu}_0+R_gT\ln C-\frac{\partial W}{\partial C}\quad (29)$$, nosso modelo desenvolvido. Sistema transporte massa acopla (19), (20), (26), e (29). Condições contorno para operação potenciostática (C(R₀,t)=C_max) e galavanostática (J(R₀,t)=j₀(1+u/R)²) são
$$C(R_0,t)=C_{\max},\;J(R_0,t)=j_0\big(1+u/R\big)^2\quad (30-31)$$ com inicial
$$C(R,0)=0\quad (32)$$.

Implementação Numérica

Soluções analíticas são intratáveis; resolvemos EDPs acopladas usando COMSOL Multiphysics, adotando parâmetros materiais listados na Tabela 1. Variáveis adimensionais usadas para clareza. Estado de carga (SOC) calculado como
$$SOC=\frac{\int_0^{R_0}C(R,t)R^2dR}{\int_0^{R_0}C_{\max}R^2dR}\quad (33)$$. Fluxos derivados dos dois modelos potencial químico são
$$J_H=\frac{\partial\tau_H}{\partial R},\;J_G=\frac{\partial\tau_G}{\partial R}\quad (34)$$

Resultados e Discussão

Figura 1 ilustra concentração lítio, tensão radial, e tensão circunferencial para nanopartícula Si esférica sob ciclagem galavanostática. Modelos desenvolvido (linhas sólidas) e tradicional (triângulos) produzem distribuições concentração quase idênticas across todos SOCs. Perfis tensão diferem apenas ligeiramente, com elevação central modesta a 46,7% e 65,5% SOC.

Figura 2 apresenta dados análogos para ciclagem potenciostática. Aqui, divergência entre dois modelos é mais pronunciada, particularmente em tensões centrais em SOCs intermediários. Como operação potenciostática impõe concentração superficial constante, a taxa de carga — e consequentemente taxa de deformação — é maior que no modo galavanostático. Deformação total para um dado SOC permanece similar, ressaltando que fator chave é taxa de deformação em vez de deformação absoluta.

Figura 3 (de Ref. [17]) demonstra que nanopartícula Si 620-nm pode expandir totalmente dentro de um minuto a 2 V versus Li, destacando relevância de modelagem deformação rápida. Em tais regimes, os dois modelos potencial químico predizem estados de tensão marcadamente diferentes, embora dados experimentais tensão permaneçam indisponíveis para validação quantitativa.

Figura 4 examina razão fluxo J_H/J_G sob ciclagem galavanostática em várias correntes carga. Razão permanece próxima à unidade across partícula, com leve aumento perto do centro e diminuição perto da superfície. Dispersão alarga com corrente mais alta, indicando sensibilidade mais forte à taxa deformação. Todas curvas intersectam perto de raio crítico — denominado região independente de potencial químico (CIR) — onde fluxos convergem independentemente do modelo escolhido. CIR desloca mais próximo da superfície conforme corrente aumenta.

Figuras 5 e 6 comparam distribuições tensão hidrostática Cauchy para os dois modelos. No CIR, σ_m é quase zero durante maior parte do período litiação, explicando por que modelos concordam lá. Fora CIR, σ_m diverge, amplificando diferenças em potencial químico e logo fluxo.

Conclusões

Introduzimos expressão potencial químico válida para deformação rápida, eliminando suposição de baixa taxa inerente a modelos tradicionais. Comparação com formulação tradicional revela diferenças desprezíveis sob ciclagem galavanostática mas desvios significativos sob controle potenciostático, impulsionados por taxa deformação. Raio crítico — o CIR — emerge onde fluxos e tensões são virtualmente idênticos, movendo-se mais próximo da superfície em correntes mais altas. Estes achados destacam necessidade de usar modelo desenvolvido para predição precisa de evolução tensão em eletrodos nanopartículas silício durante carregamento rápido.