Thermodynamics and Related Kinetics of Staging in Intercalation

Subscriber access provided by Nottingham Trent University

C: Physical Processes in Nanomaterials and Nanostructures

Thermodynamics and Related Kinetics of Staging in Intercalation Compounds Marion Chandesris, Damien Caliste, Didier Jamet, and Pascal Pochet J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.9b05298 • Publication Date (Web): 05 Sep 2019 Downloaded from pubs.acs.org on September 5, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Thermodynamics And Related Kinetics Of Staging In Intercalation Compounds

Marion Chandesrisa,b,*, Damien Calistec, Didier Jameta,b, Pascal Pochetc a

Univ. Grenoble Alpes, F-38000 Grenoble, France

b

Laboratoire de modélisation et suivi des performances, CEA Liten, DEHT, 38054 Grenoble, France

c

Department of Physics, IRIG, Univ. Grenoble Alpes and CEA, F-38000 Grenoble, France

*[email protected]

Abstract: Staging phenomena occurs in a wide range of materials. In this work, we investigate the thermodynamic properties of multi-layer free-energy models, which includes both intra-layer and inter-layer interactions. We show that a continuous screening effect for the repulsive interaction between layers is necessary to predict the occurrence of staged phases, without the non-physical nonpure ones. Included in a phase-field formulation, this high-order free energy model allows to simulate for the first time the complex staging kinetics between stage 3 and stage 2, without occurrence of the non-pure 3/2 phase, including island migration dynamics of the Daumas & Hérold model.

1 ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 19

1. Introduction Staging phenomenon corresponds to the appearance of an ordered sequence of intercalant layers separated by planes of the host-layered material in intercalation compounds. The stages are defined according to the number of host layers that periodically separates two successive intercalant layers. This phenomenon appears in a wide range of materials including clays1, transition metal dichalcogenides2, layered metal oxydes3,4, as well as graphite intercalation compounds (GICs)5,6,7,8. Intercalation is expected to allow property engineering in 2D materials9 requesting a deep understanding of the detailed intercalation sequence. In graphite, staging occurs for many intercalants, including lithium, to form LixC6, which is actually used as anode material in more than 97% of commercial lithium ion batteries10. Intercalation into such multi-layered materials proceeds through a sequence of stages with decreasing successive numbers, as suggested as soon as 1938 by Rüdorff and Hoffman while studying intercalation of alkaline metals in graphite11. For example, lithium intercalation into graphite proceeds from a dilute stage, where no internal structure can be observed, to stage 4, 3, 2L, 2, and finally stage 112. This staging phenomenon is directly visible on the equilibrium voltage of graphite, which presents successive staircases, characteristics of the successive two-phase co-existences of stages. It should be noticed that this equilibrium voltage is not symmetric around x = 1/2, x being the Li concentration in LixC6. This implies that the phase diagram of LixC6, is also not symmetric, with more regions of co-existing phases appearing below x = 1/2, as clearly evidenced by ex-situ13 or operando DRX studies14. This representation in stages, with empty and full layers, is very appealing from a thermodynamic perspective, but remains insufficient from a kinetic point of view to explain the staging transitions, especially between even and odd stages, since the diffusion of intercalants through the host layers is not possible. To account for staging kinetics, Daumas & Hérold5 introduced a model where lithium intercalants occupy all the layers as islands, locally maintaining the staging, as depicted Figure 1. Such formation of intra-layer islands has been confirmed in theoretical studies for intercalation of lithium into graphite15,16.

Figure 1 : Graphite staging structures. Schematic representation with (top) the full and empty layers of Rüdorff and Hoffman’s model, (bottom) islands domain of the Daumas and Hérold model.

Numerical simulations can help better predict the link between the crystallographic structure of these multilayered materials and their thermodynamic and kinetic properties. However, tracking of phase formation and kinetic evolution is very challenging for mathematical models. Ab initio and Monte-Carlo simulations can give access to the thermodynamic stability of stages in different configurations16,17,18, 2 ACS Paragon Plus Environment

Page 3 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


but kinetic simulations remain limited to few micro seconds given their high computational costs. In the field of lithium-ion batteries, renewed attention has been given to Cahn-Hilliard phase field models for their ability to predict and simulate phase separation in two-phase materials such as LFP19,20,21,22. In the early 80’s Safran15,23, and later Hawrylak and Subbaswamy24 introduced thermodynamic lattice gas models to predict phase diagrams and later kinetic25 of multi-layered materials such as graphite. More recently, Bazant et al.22,26,27 extended this approach to analyze the kinetic of staging in graphite until stage 2, using a bilayer model and focusing on thermodynamically consistent reaction boundary conditions. In the present work, we focus on higher order free-energy models to capture the complex kinetic of staging between stage 2 and stage 3, and the prediction of the corresponding non-symmetric equilibrium voltage which occurs only for stages strictly greater than 2. First, we present the mathematical model, including the phase field formulation and the general form of the free-energy model that includes two contributions, an intra-layer and an inter-layer one. Key to predict staging is the form of the free-energy model and especially the inter-layer interactions. Therefore, models with increasing complexity of inter-layer interactions are presented together with their corresponding phase diagrams. This allows us to develop an original next-nearest neighbor term to capture stage 2 / stage 3 staging kinetics for both intercalation and de-intercalation without the occurrence of non-pure stages, as presented in the simulations of the last section. 2. Mathematical model: phase-field formulation We follow the conventional Cahn-Hilliard phase-field model applied28 to multilayer intercalation materials25,22,27. We consider an intercalation host material composed of N stacked layers along the zaxis. Intercalant can move between two layers of the host material, but cannot cross host material’s layers. This assumption is reasonable for many intercalant materials, such as graphite15 or layered metal oxydes3,29 given the high barriers for migration through the host material ( ≈ 10 𝑒𝑉 for graphite16). The total Gibbs free energy of the system, 𝐺, can be written as: 𝑁

𝐺=

∭𝑔 𝑑𝑥𝑑𝑦𝑑𝑧 = ∑𝐿 (∬𝑔 𝑑𝑥𝑑𝑦) 𝑗

𝑗

(1)

𝑗=1

where 𝑔 is the local Gibbs free-energy density and 𝑔𝑗is the Gibbs free-energy density of the 𝑗th layer of height 𝐿𝑗. The model chosen for this free-energy density determines most of the physical properties of the system. a. Thermodynamic model Following previous works23,25,27,30, we assume that the Gibbs free-energy density of the 𝑗th layer is the sum of an intra-layer energy 𝑔𝑙𝑎𝑦𝑒𝑟 which depends only on the configuration of the 𝑗th layer and an 𝑗 𝑖𝑛𝑡𝑒𝑟 inter-layer repulsion energy 𝑔𝑗 which accounts for the inter-planar interactions of intercalant clusters: 𝑔𝑗 = 𝑔𝑙𝑎𝑦𝑒𝑟 + 𝑔𝑖𝑛𝑡𝑒𝑟 𝑗 𝑗

(2)

Let 𝑐𝑗 be the dimensionless normalized concentration of intercalant on the 𝑗th layer. For the intra-layer energy, the simplest free energy model that allows phase-separation within the layer is the regular solution model:



Page 4 of 19

1 (𝑐𝑗,∇𝑐𝑗) = 𝑔0𝑗 (𝑐𝑗) + 𝜅(∇𝑐𝑗)2 𝑔𝑙𝑎𝑦𝑒𝑟 𝑗 2

[

𝑔0𝑗 (𝑐𝑗) = 𝑁𝑉 𝑘𝐵𝑇(𝑐𝑗 ln(𝑐𝑗) + (1 ― 𝑐𝑗) ln(1 ― 𝑐𝑗) ) + Ω𝑎𝑐𝑗(1 ― 𝑐𝑗) + 𝑐𝑗

(3) 𝜇𝑒𝑞 𝒩𝐴

]

(4)

where 𝑁𝑉 is the number of intercalation sites per unit volume for one layer, 𝑘𝐵 the Boltzmann’s constant, 𝑇 the temperature, Ω𝑎 is the average energy density of the interaction between intercalant particles and vacancies, 𝜇𝑒𝑞 is the equilibrium chemical potential of the regular solution model. For a given intercalation compound, 𝑁𝑉 can be estimated by 𝑁𝑉 = 𝒩𝐴𝑐max, where 𝒩𝐴 is the Avogadro number and 𝑐max is the maximum intercalant concentration of the host material. The first term in Eq. (3), 𝑔0𝑗 , is the classical part of the free energy. The second term is a non-local term introduced by Cahn & Hilliard28, following Van der Walls theory. It allows both to regularize the interface description between phases and to associate an energy to this interface via the energy penalty coefficient relative to the concentration gradient, 𝜅. In the regular solution model, the classical part of the free energy, 𝑔0𝑗 , is composed first, of an entropic contribution which favors the mixing of the system, second, by an enthalpic term which promotes the separation of the system in two phases of high and low concentration. For the inter-layer interactions, different model exists in the literature. To study the kinetic intercalation in graphite, Hawrylak and Subbaswamy25 considers an interlayer repulsion energy that is a sum of interactions between all the layers, the repulsion energy decreasing with distance between the layers by an exponent -4: 𝑔𝑖𝑛𝑡𝑒𝑟 = 𝑗

Ω0𝑐𝑗𝑐𝑖

∑|𝑖 ― 𝑗|

4

(5)

𝑖>𝑗

Hess31 uses a similar expression but considering an exponent of -5, while Millman et al.32,33 analyze the impact of different exponent to represent c-axis spacing typical of either graphite or large intercalant species such as transition metal chlorides. They also introduced the notion of “strongly screened” potentials, where the repulsive interaction between two occupied layers is set to zero if there is some occupied intercalant layers between them and study the impact of this hypothesis on the material phase diagram. More recently, Fergusson & Bazant22 and Smith et al.27 introduced a first neighbor repulsive interaction in their bilayer model incorporating 2- and 4- body interactions that allowed them to simulate the intercalation kinetics between stage dilute, stage 2 and stage 1. The form of the interlayer interaction directly affects the phase diagram of the studied system and several expressions will be studied in this work. b. Bulk transport In such multi-layer phase-field models, the approach is semi-continuous, in the sense that inside a given layer ((x,y) plane), we follow the coarse-grained dimensionless intercalant density, 𝑐𝑗, but in the z-direction, each layer of the material is considered. The kinetic evolution of this system is governed by the mass balance equation in each layer 𝑗: 𝒩𝐴

∂𝑐𝑗 ∂𝑡

= ― ∇𝑥,𝑦.𝐽𝑗

(6)


Page 5 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


where the intercalant flux, 𝐽𝑗, in the 𝑗th layer is linearly related to the gradient of the local chemical potential 𝜇𝑗, via the mobility 𝑀: (7)

𝐽𝑗 = ― 𝑀∇𝑥,𝑦𝜇𝑗 Here, the chemical potential is computed as the variational derivative of the free energy: 𝜇𝑗 =

1 𝛿𝑔 𝑐max𝛿𝑐𝑗

(8)

3. Inter-layer interactions In order to develop an inter-layer interaction model able to predict the complex non-symmetric transitions to stages greater than 2, we follow the methodology presented by Schon et al 2 and analyze the thermodynamic properties of the total Gibbs free-energy considering various inter-layer interaction models with increasing complexity. The models analyzed in the present work are not aimed to provide quantitative values for physical quantities of interest of given multi-phase layered materials. The aim is to provide a conceptual framework to analyze staging and phase diagrams depending on inter-layer interactions. Therefore, the values of the different physical parameters of the models are chosen arbitrarily, for illustration purpose. a. No inter-layer interactions If we assume no-interlayer interactions, the system reduces to the classical Cahn-Hilliard model with a Gibbs free energy based on the regular solution of the lattice gas model. Spontaneous phase separation in two phases occurs when the system is unstable, i.e., the second derivative of the Gibbs free energy versus the filling fraction is negative. The zero values delimit the spinodal region. For the regular solution model (Eq. (4)), the spinodal curves exist only when Ω𝑎 > 2𝑘𝐵𝑇 , and the corresponding filling fractions 𝑐𝑠𝑝 ± are defined by: 𝑐𝑠𝑝 ± =

(

2𝑘𝐵𝑇 1 1± 1― 2 𝛺𝑎

)

(9)

The spinodal curves are presented in dashed lines (Figure 2(b)), together with the corresponding concentrations, 𝑐𝑠𝑝 ± , obtained at room temperature in red dots. In this figure, the values of 𝛺𝑎 and 𝜇𝑒𝑞 have been set to arbitrary values. Once phase separation occurs, the composition of the two phases is not arbitrary and can be determined by studying the thermodynamic equilibrium condition of the system under the constraint of mass conservation. This methodology is sometimes referred as the common tangent construction, 𝑒𝑞 as presented (Figure 2(a)). The equilibrium concentrations in the two phases, 𝑐𝑒𝑞 + and 𝑐 ― , are such 0( 𝑒𝑞 ) 0( 𝑒𝑞 ) that the classical part of their chemical potentials are equals, 𝜇 𝑐 + = 𝜇 𝑐 ― , and they verify the 0 𝑒𝑞 𝑒𝑞 𝑒𝑞 following equation, 𝑔0(𝑐𝑒𝑞 + ) ― 𝜇𝑒𝑞𝑐 + = 𝑔 (𝑐 ― ) ― 𝜇𝑒𝑞𝑐 ― . These values can be computed and are presented in the phase diagram of Figure 2(b) in plain lines, the values obtained at room temperature 𝑒𝑞 being in green dots. The system separates into two phases of concentration 𝑐𝑒𝑞 + and 𝑐 ― , and the ratio of the different phases is deduced from the knowledge of the total filling fraction of the system.



Page 6 of 19

Such Gibbs free-energy models without inter-layer interactions are very well suited to describe phase separation for non-layered materials and have been extensively used to study two-phase materials such as LiFePO4. Implemented in phase-field simulations they allow to analyze the intercalation dynamics of lithium focusing on the kinetics of phase separation34,19,20,21,22. Nevertheless, to study staging in intercalation compounds, inter-layer interactions must be added.

Figure 2 : (a) Gibbs free energy curve of the general model with no interlayer interactions, for 𝛺𝑎 = 64,3𝑚𝑒𝑉 and 𝜇𝑒𝑞/𝒩𝐴 = 1,8𝑚𝑒𝑉 at room temperature. (b) Corresponding phase diagram. Plain line: equilibrium curve with corresponding equilibrium concentrations at room temperature in green; dashed line: spinodal curve with corresponding spinodal concentrations at room temperature in red.

b. First neighbors interaction Our approach to analyze staging is to consider first only inter-layer interaction to the first neighbors. This configuration allows studying the occurrence of only three stable phases: stage 1, stage 2 and dilute stage 1. To conduct mathematical developments, we choose the simplest form for the first neighbor interaction, with a term that energetically penalizes high concentration of intercalants in adjacent layers: 1 𝑔𝑖𝑛𝑡𝑒𝑟 = 𝑁𝑉 Ω (𝑐𝑗 + 1𝑐𝑗 + 𝑐𝑗 ― 1𝑐𝑗) 𝑗 2 𝑏

(10)

With a first neighbor interaction, one can consider only two layers to analyze the model’s thermodynamic properties. The Gibbs free energy density of the model depends only on the filling fractions 𝑐1and 𝑐2 in these two layers and is plotted Figure 3(a), without the contribution from the non6 ACS Paragon Plus Environment

Page 7 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


local term (concentration gradient). For this system, the spinodal curves, which delimits the filling fraction in each layer for which spontaneous phase separation occurs, are determined from the thermodynamic stability criteria: ∂2𝑔 𝑑𝑒𝑡 1. The obtained filling fractions are displayed in blue in Figure 3(a). It should be noticed that the “full” and “empty” filling fractions of stage 2 are lower when stage 2 is in equilibrium with dilute stage 1 (domain 𝑐1 + 𝑐2 < 1), than when it is in equilibrium with stage 1 (domain 𝑐1 + 𝑐2 > 1). From a thermodynamic point of view, this is important, otherwise, we would observe the same equilibrium potential for the two two-phase regions. It implies that stage 2 undergoes a small solid solution behavior, with increasing filling fraction in both layers, between the two occurrences of the two-phase regions. The phase diagram of this model is obtained by solving the thermodynamic equilibrium conditions for different temperatures and is displayed in Figure 3 (b). In this figure, the filling fraction is the mean filling fraction of the whole system, i.e. equals to (𝑐1 + 𝑐2)/2. The phase diagram is solved first in the 𝑐1 + 𝑐2 < 1 domain, which corresponds to inter-layer phase separation in dilute stage 1/stage 2, then in the 𝑐1 + 𝑐2 > 1 domain, which corresponds to inter-layer phase separation in stage 2/stage 1. The 7 ACS Paragon Plus Environment


Page 8 of 19

proposed thermodynamic model allows to predict the successive occurrence, with increasing total filling fraction, of dilute stage 1, co-existence of dilute stage 1 and stage 2, stage 2 only, co-existence of stage 2 and stage 1, and finally only stage 1. In the upper part of the diagram, the solution corresponding to inter-layer phase separation does not exist. In this phase diagram, full and empty filling fractions play the same role in the thermodynamic model, which leads to a symmetric diagram around 𝑐 = 0.5. Kinetics of staging between these three phases can be studied by implementing this Gibbs free-energy model in the general phase-field formulation presented in section 2, but since we are interested by the complex stage 2/stage 3 staging kinetics, we go further and study the inter layer interactions to the second neighbors.

Figure 3 : (a) Room temperature contour plot of the Gibbs free energy model in [J/mol] with interlayer interactions to the first neighbor (𝛺𝑎 = 64,3 𝑚𝑒𝑉, 𝛺𝑏 = 5.14 𝑚𝑒𝑉 and 𝜇𝑒𝑞 = 0) as function of the filling fractions c1 and c2 of the two adjacent layers. Red lines are the spinodal curves solution of equation (12), green (blue) dots are the equilibrium filling fractions in the dilute stage 1/stage 2 configuration (stage 2/stage 1 configuration) respectively. (b) Corresponding phase diagram as function of the total filling fraction.

c. Second neighbors interaction In order to study staging phenomena to higher stages, it is necessary to introduce inter-layer interactions in a longer range. An inter-layer repulsive free energy model is thus introduced which energetically penalizes high concentration of intercalants in the next neighbor layers in addition to the adjacent ones: 1 𝑔𝑖𝑛𝑡𝑒𝑟 = 𝑁𝑉 [Ω𝑏(𝑐𝑗 + 1𝑐𝑗 + 𝑐𝑗 ― 1𝑐𝑗) + Ω𝑐(𝑐𝑗 + 2𝑐𝑗 + 𝑐𝑗 ― 2𝑐𝑗)] 𝑗 2

(13) 8

ACS Paragon Plus Environment

Page 9 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


With first and second neighbor’s interactions, six layers are necessary to analyze the thermodynamic properties of the system. While it is not possible to simply represent this free-energy model in six dimensions, one can still compute the corresponding phase diagram by studying the thermodynamic equilibrium condition of the system under the constraint of mass conservation. The problem can be solved numerically, using the same methodology as in the previous section and the obtained result is presented in Figure 4 with the following set of parameters: Ωa = 61,7 meV, Ωb = 23.1 meV, Ωc = 4.1 meV and µeq = 0. With these parameters, occurrences of dilute stage 1, stage 3, stage 2 and stage 1 are indeed predicted at room temperature. However, the symmetric form of the second neighbor repulsive energy also implies the existence of a non-pure stage, called stage “3/2” following the notation of Millman & Kirczenow32, to refer to a periodic structure in which the period consists of three host layers and two layers with high intercalant concentration. The occurrence of non-entire staging phases has been reported in other theoretical studies25,23,30,32,24, supporting the view that interlayer interaction is effective at long range. Experimentally, the fractional 3/2 stage is not observed for lithium intercalation in graphite, but has been observed for potassium and rubidium intercalation in graphite under high pressure35,36.

Figure 4 : Phase diagrams of the general model with interlayer interactions to the first and second neighbors for 𝛺𝑎 = 61,7 𝑚𝑒𝑉, 𝛺𝑏 = 23.1 𝑚𝑒𝑉, 𝛺𝑐 = 4.1 𝑚𝑒𝑉 and 𝜇𝑒𝑞 = 0. Colored lines indicate equilibrium between: black: dilute-stage 3 and



Page 10 of 19

stage3/2-stage1; red: stage 3-stage 2 and stage2-stage 3/2; blue: dilute-stage 2 and stage2-stage1. The dashed lines correspond to metastable solutions, energetically less favorable than the ones drawn with plain lines. (a) “standard” repulsive interaction at second neighbor, (b) dissymmetric repulsive interaction at second neighbor with continuous screening effect.

In the phase diagram of Figure 4(a), one can observe that above a critical temperature 𝑇𝑐, the existences of stage 3 and stage 3/2 do not occur any more. When increasing the filling fraction from dilute stage 1, the system displays an inter-layer phase separation directly to dilute stage 1/stage 2; and, after the single-phase stage 2 solution, the system phase separates to stage 2/stage 1. We have observed numerically that this critical temperature depends on the ratio Ω𝑐/Ω𝑏 the other parameters, Ω𝑎 and Ω𝑏, being fixed. For a given temperature, above a given ratio of Ω𝑐/Ω𝑏 = 𝑋, the repulsive energy to the second neighbor is high enough to lead to the formation of stage 3, but the drawback is the occurrence of a non-pure stage 3/2 (See Figure 5). Below this value the repulsive energy to the second neighbor is not sufficient and only inter-layer phase separation to stage 2 occurs, without formation of stage 3, but the advantage is the disappearance of the non-wanted non-pure stage 3/2, as schematically presented Figure 5. This observation leads us to propose another second neighbor inter-layer model, where high concentrations of intercalant in the next neighbor layers are energetically penalized, except if intercalant is present in between: 1 𝑔𝑖𝑛𝑡𝑒𝑟 = 𝑁𝑉 [Ω𝑏(𝑐𝑗 + 1𝑐𝑗 + 𝑐𝑗 ― 1𝑐𝑗) + Ω𝑐(𝑐𝑗 + 2(1 ― 𝑐𝑗 + 1)𝑐𝑗 + 𝑐𝑗 ― 2(1 ― 𝑐𝑗 ― 1)𝑐𝑗)] 𝑗 2

(14)

Such an approach allows mimicking two values for the second neighbor parameter:  Ω𝑐 (1 ― 𝑐𝑗 + 1) ~ Ω𝑐 for low filling fraction value  Ω𝑐 (1 ― 𝑐𝑗 + 1) ~ 0 for high filling fraction value

Figure 5 : Schematic representations of staging structures encountered for the general model with interlayer interactions to the first and second neighbors with (top) a repulsive energy to the second neighbor high enough to obtain stage 3, but leading also to non-pure stage 3/2 formation; (bottom) repulsive energy to the second neighbor low enough to get rid of non-pure stage 3/2 at high filling fraction; (middle) screened repulsive energy to second neighbor.


Page 11 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


The phase diagram corresponding to the dissymmetric expression of the inter-layer interaction (14) can be computed numerically and the obtained result is presented Figure 4(b). The phase diagram is now dissymmetric. For mean filling fraction above 0.5, the equilibrium filling fractions corresponding to the phase separation in stage 2-stage 3/2 and stage 3/2-stage 1 have been determined numerically using the same procedure as in the previous section. They are presented in dashed lines, since they are energetically less favorable. Furthermore, one can note that there is no thermodynamic solution for the occurrence of the non-pure stage 3/2. The more energetically favorable case is the interlayer phase separation in stage 2-stage 1. On the contrary, for filling fraction below 0.5, the system displays successive inter-layer phase separations in dilute stage 1/stage3 and then stage 3/stage 2. The proposed dissymmetric expression for the second neighbor repulsive interaction (14) allows obtaining a phase diagram, such that, at room temperature, only phase separations in pure stages (1, 2 and 3) are obtained. The proposed expression for the second neighbor interaction (14) corresponds to the introduction of a continuous screening of the repulsion energy of the form (1 ― 𝑐) between a pair of intercalant layers, where 𝑐 is the dimensionless concentration of intercalant in the layer in between. The introduction of a screening effect has been suggested by Millman et al.32,33, but in a discontinuous way, considering a zero repulsion energy if the layer in between is occupied. The purpose of this screening effect was also to predict the occurrence of staged phases, without the non-physical non-pure ones. The continuous expression proposed in this work corresponds to a screening effect that is proportional to the intercalant concentration. This is physically more coherent even though the physical nature of this effect (i.e. elastic, electrostatic, …) has still yet to be determined. Furthermore, it is also numerically more convenient when studying not only phase diagram but also the kinetics of staging with phase field simulations. d. Equilibrium potentials In this framework, the equilibrium potential evolution with filling fraction is not a model’s input, but an emerging property of the free energy model. At a given temperature, the equilibrium voltage, 𝑉𝑒𝑞, can be determined from the following relation18,27: 𝑉𝑒𝑞 = 𝐸0 ―

𝜇eff 𝐹

(15)

where 𝐸0 is the standard redox potential, 𝐹 is Faraday’s constant, 𝜇eff is equal to the classical part of the chemical potential, 𝜇0, in the single-phase regions and the equilibrium potential solution of the equilibrium thermodynamic condition, 𝜇𝑒𝑞, in the two-phase ones. The results obtained for the four different energy models studied in this work are displayed Figure 6. In solid solution regions (single phase), the equilibrium voltage varies with filling fraction, as the Gibbs free-energy density from which it is derived. In two-phase regions, plateaus are observed. They are characteristics of the two-phase transformations. They correspond to the fact that, at thermodynamic equilibrium, the chemical potential is the same in the two phases, and does not evolve as long as the two phases are co-existing.



Page 12 of 19

Figure 6 : Equilibrium potential curves at room temperatures as function of the mean filling fraction corresponding to the four studied free-energy models, assuming 𝐸0 = 0: no interlayer interaction with 𝛺𝑎 = 64,3𝑚𝑒𝑉 and 𝜇𝑒𝑞/𝒩𝐴 = 1,8𝑚𝑒𝑉 (a); first neighbor interaction with 𝛺𝑎 = 64,3 𝑚𝑒𝑉, 𝛺𝑏 = 5.14 𝑚𝑒𝑉 and 𝜇𝑒𝑞 = 0 (b); first and second-neighbor interaction with 𝛺𝑎 = 61,7 𝑚𝑒𝑉, 𝛺𝑏 = 23.1 𝑚𝑒𝑉, 𝛺𝑐 = 4.1 𝑚𝑒𝑉 and 𝜇𝑒𝑞 = 0 for (c) “standard” repulsive interaction at second neighbor, (d) dissymmetric repulsive interaction at second neighbor with continuous screening effect.

When multiple two-phase regions are predicted, (Figure 6 (b), (c) and (d)), single-phase regions exists in between. In these transition regions, only one phase is present as depicted schematically in Figure 6. These single phases display a solid solution behavior as the mean filling fraction is increased, with increasing concentration in each layer and decreasing equilibrium potential. The concentrations in the “full” and “empty” layers of these single phases are varying in these transition regions, with different equilibrium values when equilibrium is reached with the two adjacent phases. This concentration variation in the solid solution regions can also be intuited given their non-zero filling fraction span. The potential difference between two successive plateaus, as well as the stoichiometric ranges of the single-phase transition regions are directly linked to the values of the three parameters of the energy model 𝛺𝑎, 𝛺𝑏, and 𝛺𝑐. This dependence provides a direct link between microscopic (intra-layer and inter-layer repulsive energies) and the macroscopic properties of the studied material. Symmetric equilibrium potential curves around the filling fraction 0.5 are obtained for symmetric first and second neighbor inter-layers interactions (Figure 6 (b), (c)), while the continuous screening at second neighbors gives a non-symmetric equilibrium potential with two plateaus in the region below 𝑐 = 0.5 and only one plateau in the region above. This later case is more representative of what is observed for Li intercalation into graphite12, even though the present model predicts occurrence of stage 3 around a mean filling fraction of 0.33, while it is reported to occur around 0.212,14. A comprehensive model for Li intercalation into graphite would require accounting also for the formation of stage 2L. Nevertheless, this simple modeling approach allows predicting the main thermodynamic properties of the system, given the form of the inter-layer interactions and can easily be extended to higher order staging, including a screening effect at the desired neighbors depending 12 ACS Paragon Plus Environment

Page 13 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


on the symmetry properties of the measured equilibrium potential and observed phases. For example, occurrence of fractional stage 3/2 at high pressure for potassium and rubidium in GICs35,36 suggests that screening effect by intercalant becomes ineffective when increasing the pressure leading to an effective repulsion to second neighbor for these intercalants. Given the physical grounds of the thermodynamic model, ab initio computations could be used to test this hypothesis or to quantify the repulsive energies to first and second neighbors of different intercalants. 4. Staging simulations In order to study the kinetics of staging between stage 2 and stage 3, the inter-layer interaction model with continuous screening effect (Equation (14)) is introduced in the phase-field formulation of section 2. A system consisting of six layers is considered with intra-layer filling fraction independent of the 𝑦 direction. Mass conservation and transport equations (6) and (7) are solved in the six layers, considering only the 𝑥 direction and a planar configuration of length 25µm. The equations in the six layers are coupled via the dependence of the intra-layer chemical potential 𝜇𝑗, to the filling fractions in the other layers. Periodic conditions are used in the 𝑧 direction while a zero flux boundary condition is applied at 𝑥 = 25µ𝑚. At 𝑥 = 0 , where intercalation occurs, the incoming flux in layer j, 𝐽𝑗 ,has been assumed to be of the following form: 𝑖𝑛 (16) 𝐽𝑗 = 𝐹 where the intra-layer intercalation current, 𝑖𝑛, is modeled without accounting for the dependence to the local potential, while keeping a dependence to the local concentration, 𝑖𝑛 = 𝐼𝑖𝑚𝑝𝑐𝑗(1 ― 𝑐𝑗), using a constant pre-factor called 𝐼𝑖𝑚𝑝. The idea is to introduce as simply as possible the dependence of the intra-layer intercalation reaction to the surface local concentration recommended in references 26,27, without focusing on the definition of thermodynamically consistent reaction boundary conditions, which is out of the scope of this work. Finite volume difference is used to integrate the equations in the 𝑥 direction, ensuring that the grid spacing is less than the scale of the interfacial width, which can

(

be estimated using the following expression: 𝜆~

1/2 𝜅 . (Ω𝑎𝒩𝐴𝑐max/6)

)

Time integration is carried out using

a non-linear GMRES method. Initial filling fractions are set to 0.24 (resp. 0.76) for intercalation (resp. deintercalation) simulations. Very small perturbations are introduced between the initial filling fractions of the six layers to break the perfect numerical symmetry of the simulation. Table 1 summarized the physical parameters used in the simulations. The value of 𝐼𝑖𝑚𝑝 has been chosen to simulate regimes of roughly 2C for both intercalation and de-intercalation Results for the two simulations, intercalation and de-intercalation, are presented Figure 7 (two movies of the simulations are included in the supplement). For intercalation, given the initial filling fraction, the system becomes unstable at very early times at the particle surface and transition to stage 2 can be observed Figure 7 (a) looking at the pattern with the islands formation predicted by the Daumas & Hérold model. Nevertheless, since this configuration in stage 2 is not the most stable one, transition to stage 3 occurs rapidly and coexists with dilute stage 1 (Figure 7 (b)). As predicted by the thermodynamic analysis, the filling fraction of dilute stage 1 is higher than the low filling fraction of stage 3 (and 2). At later times, stage 2 forms at the surface of the particle and co-exist with stage 3 (Figure 7 (c)), propagates to the center of the particle via formation and migration of intercalant islands (Figure 7 (d)) and occupy all the particle until transition to stage 1 occurs at the particle surface. Stage 1, then co-exist with stage 2 (Figure 7 (e)), until complete filling of the particle (Figure 7 (f)). Migration and fusion followed by formation of islands is only observed during the stage 3/stage 2 transition. 13 ACS Paragon Plus Environment


Parameter 𝛀𝒂 𝛀𝒃 𝛀𝒄 𝑴 𝑫 𝜿 𝑻 𝑰𝒊𝒎𝒑 𝒄max

Page 14 of 19

Value 2.4 𝑘𝐵𝑇 = 61,7𝑚𝑒𝑉 0.9 𝑘𝐵𝑇 = 23.1𝑚𝑒𝑉 0.32 𝑘𝐵𝑇 = 8.2𝑚𝑒𝑉 𝐷 𝑘𝐵𝑇 5e-13 m2/s 3e-6 J/m 298 K ± 6 mA/cm² 30 000 mol/m3

Table 1 : Physical parameters used in the staging simulations

For de-interacalation, the system also becomes unstable at very early times given the initial filling fraction, leading to formation of stage 2 at the particle surface (Figure 7 (g)). Stage 2 co-exist with stage 1. The lower concentration of the “full” state of stage 1 compared to the “full” state of stage 2 predicted by the thermodynamic study is visible on this figure. De-intercalation then proceeds by the extension of the stage 2 phase by growing size of the islands and formation of new ones (Figure 7 (h)). Once stage 2 occupies the entire particle (Figure 7 (i)), we can observe a decrease of the concentration of the “full” state, before transition to phase 3 occurs at the particle surface by migration, formation and disappearance of some islands (Figure 7 (j)). Once stage 3 fills all the particle, transition to dilute stage 1 proceeds by progressive disappearance of the islands (Figure 7 (k-l)).

5. Conclusion In this paper, we have presented several high order free energy models focusing on staging phenomena in multi-layer intercalation materials. The general form of the free-energy models includes two contributions: an intra-layer and an inter-layer one, the latter being key to predict staging. Models with increasing complexity of inter-layer interactions have been presented, including their corresponding phase diagram and equilibrium potential at room temperature. We have shown that a continuous screening effect in the next-nearest neighbor interaction term is crucial for the formation of the simple stages observed experimentally for lithium intercalation in graphite without the occurrence of the non-pure ones. This continuous screening effect leads to non-symmetric phase diagrams, as reported experimentally for several multi-layers materials. With few parameters, the proposed model allows to simulate the complex stage 2 / stage 3 staging dynamic for both intercalation and de-intercalation, including the intercalant islands formation predicted by the Daumas & Hérold model. Furthermore, since the approach is based on physical considerations, it could be further parameterized using ab initio calculations to better decipher between the various properties of layered intercalation compounds. Indeed, this mesoscopic approach provides a clear connection between the structural properties of the atomistic scale and the measured macroscopic properties of these compounds.


Page 15 of 19 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(a)

(g)

(b)

(h)

(c)

(i)

(d)

(j)

(e)

(k)

(f)

(l)

Figure 7 : Simulation of the intercalation (left), respectively de-intercalation (right) processes for a six-layer geometry of a 25µm particle. Intercalant concentration at six successive times during intercalation. Between white planes of the host layers, intercalant content is referred with a white (low)/black (high) scale. Intercalation times (from top to bottom): 20s, 60s, 320s, 620s, 1080s and 1500s, de-intercalation times (from top to bottom): 20s, 310s, 660s, 1030s, 1290s 1480s (movies in supplement). The colors on the top are assigned based on the local mean concentration over the six-layer to give indications of stages: yellow for stage 1 (c> 0.75), red for stage 2 (0.43

Thermodynamics and Related Kinetics of Staging in Intercalation

Recommend Documents