A MODEL FOR ELECTROCHEMICAL INSERTION LIMITED BY A PHASE TRANSITION PROCESS ( EILPT )

This paper deals with electrochemical insertion into a cathodic material. New results on modelling of the influence of a solid phase transformation on the shape of voltamograms are presented. The original experiments concern the insertion of sodium into carbon during the cathodic reduction of molten NaF at 1020 °C, but in the present manuscript emphasis on the theoretical aspects of the work is put. Phase transformations during electrochemical insertion are taken into account, with various values for parameters such as the thermodynamic biphase equilibrium potential, the compared diffusion and phase transformation kinetics, and the electrode thickness. The voltamograms calculated present very specific features; some of them have already been observed experimentally in literature.


Introduction
Many electrochemical processes involve phase transition phenomena.The electrodeposition of metals and the creation of passive films are well known examples for which theoretical models are well established [1,2].
During electrochemical insertion, the incorporation of "guest" species into the lattice of the host material, accompanied by the injection (or removal) of charge compensating electrons into (or from) the electronic band, leads to the creation of a new homogeneous phase characterized by a new crystallographic structure and electronic properties.If the new phase is not too different from the original material, the insertion process will take place as if it were in a perfectly non-stoichiometric system.Various theoretical studies concerning such a phenomenon have been published.However, if structural and electronic changes are sufficiently important to create limiting compositions and induce big changes in the transport properties of the guest species, biphase domains may be created [3].In that case, the kinetics of insertion process will depend on the mobility of the interface boundary between the two phases, which is controlled by a phase transition phenomenon.If the latter process is slow compared to ionic diffusion within both phases, the kinetics of insertion will not obey classical diffusion laws.
When studying the kinetics of hydrogen insertion into palladium by chronoamperometric technique Chen et al. [4] obtained particular chronoamperograms showing a quasi-plateau region.They attributed this behavior to the occurrence of a two-phase equilibrium between α (poor) and β (rich) phases of palladium hydride.Using a different technique (electrochemical potential spectroscopy) to study the insertion of H into LaNi 4.5 Al 0.5 , Li into γ-Fe 2 O 3 [5] and Li into C 60 [6], Chabre et al. obtained voltamograms presenting narrow and asymmetric non-diffusional peaks, characterized by their stripping peak shape.These results were also ascribed to the control of the insertion process by a phase transition phenomenon.Although many experimental studies on the electrochemical insertion process involving a phase transition phenomenon have been published, the interpretations remained very qualitative and theoretical aspects have been paid little attention so far.R. Ash and R.M. Barer [7] studied the diffusion problem in a biphase medium, considering diffusion in each phase as a quasi-stationary process with a linear concentration profile.They assumed that the concentration of the intercalated species is constant at the interface boundary, which leads to the phase boundary displacement rate proportional to t 1/2 .M. W. Verbrugge and B. J. Koch [8] proposed a model for the diffusion of Li in carbon microfiber taking into account the variation of activity of Li with the concentration in the host material.
E. Wicke and H. Brodowski [9] studied the classical diffusion of atomic hydrogen in Pd and proposed various experimental means for studying and characterizing these phenomena.
The aim of the present work is to develop a mathematical model which can be used to account for the experimental results of electrochemical insertion processes limited by a phase transition phenomenon (EILPT) and allow an accurate quantitative analysis.It will be mainly applied to modeling linear sweep voltamograms.

Background of phase transition kinetics
Thermodynamically, every phase transition can be considered as a changeover from configuration characterized by its initial Gibbs function (G i ) to a new structural arrangement having Gibbs function (G f ).The difference (G f -G i ) is the driving force (free energy) for the phase transition.It depends, for a fixed temperature and pressure, on the difference of composition between the actual state and the biphase equilibrium state.To describe the evolution of the system from the initial to the final state, a kinetic model can be established using considerations on transition processes on the atomic scale.A phas e transition is, strictly speaking, made of two steps.The first one involves the germination process, which consists in the formation of nuclei as discrete centers.This first step, which can be instantaneous or progressive, is commonly considered as thermally activated.The second step concerns the phase growth, which may be a complex process, mainly when the transition is accompanied by composition changes.In that case, the phase growth involves both diffusional transport in each phase and mass transfer process across the interface boundary.Phase growth kinetics can be, in a simplified manner, divided into two categories depending on whether it is limited by diffusion phenomena within the phases or by mass transfer across the phase boundary.In the first case, diffusion within the phases is slower than mass transfer across the interface.According to the Fick's laws, it induces a concentration gradient in both phases.Local phase equilibrium can be assumed at the interface boundary.The flux conservation at the interface can be expressed by equation ( 1): (1 where ξ is the interface abscissa and D α , D β are diffusion coefficients in α and β phases (see Fig. 1).
In the second case, when diffusion of the species involved is faster than their transfer across the interface, the concentration gradients in both phases may be neglected.The growth process is similar to a thermally activated germination phenomenon [10,11], and the interface displacement rate is linked to the phase transition activation energy (∆G * ) and to the phase transition driving energy (∆G α β ) by the following equation: where ν β and ν α are jump frequencies from α to β and in the opposite direction, respectively.
The δ parameter is the induced interface displacement following the transfer of an atom from α to β.
The commonly used assumption of jump frequencies equality (ν β =ν α ) allows equation (2a) to be simplified to: The k ph coefficient is called phase transition rate constant.
The two categories of phase growth kinetics presented above correspond to borderline cases.In general, phase transition is partially controlled both by diffusion and transfer across the interface.Therefore, in a general model, it is necessary to combine diffusion laws with equation (2b) corresponding to the transfer process across the interface.

Model formulation
Let us consider a simple reduction process, consisting in the electrochemical insertion of guest species (M) into a host lattice structure <H> (with thickness e) and suppose the occurrence of a phase transition process (α-β) (Figure 1).In order to

Model for the electrochemical insertion process accompanied by a phase transition phenomenon.
simplify the mathematical model and keep the computational work at a reasonable level, we have made the following assumptions: A1-the phase transition front is considered as an ideal planar surface, A2-phase (α) is considered as the pure initial structure and the concentration of guest species is neglected.This assumption can never be satisfied experimentally but it is necessary for enabling the solution of the mathematical system, A3-phase (β) is considered as an ideal solid solution in which electrostatic interactions between guest species are neglected.Expansion of the host structure is also taken as negligible, A4-ionic diffusion in the electrolyte (D=10 -5 cm 2 s -1 ) is high compared with the ionic diffusion in the host material (D=10 -8 to 10 -12 cm 2 s -1 ).Diffusion in the electrolyte will not limit the electroinsertion process and the concentration of (M + ) ions at the electrolyte-electrode interface can therefore be considered as constant (C M + ), A5-electronic conductivity of the host material is sufficiently high to neglect the contribution of transport by migration within the host structure.
No assumption concerning the phase transition kinetics is made.The general model will be solved combining diffusion in the phases and transfer kinetics laws across the two-phase boundary.

Basic equations
*Thermodynamic relation: on the basis of assumption (A3), guest species will be randomly distributed over equivalent sites and then obey Fermi-Dirac statistics [3].The chemical potential of (M) inserted in phase β is given by: is the relative occupancy of insertion sites in phase β and is the standard chemical potential of the (M + ) species in the β phase.
*Transport laws: assuming pure linear diffusional transport and taking into account the insertion isotherm, the motion of (M + ) in the host structure can be treated according to the classical Fick's law: The boundary conditions at the electrode-electrolyte interface (x=0) and at the twophases boundary (x=ξ) are: (5) (6) where i is current density.Equation (6a) is derived from equation ( 1) by setting C α =0 and it expresses the absence of atoms accumulation at the two-phase separation, whereas (6b) expresses the absence of species diffusion across the electrode bottom (x=e) when the phase transition is completely achieved.
*Phase transition kinetics: phase growth kinetics is assumed to obey Perturbating the biphase equilibrium by increasing the concentration of one of its components (M) induces growth of the phase where (M) species have the lowest chemical potential, and the system reaches its initial equilibrium state.Hence, there is a force, which tends to bring the system back to its equilibrium state.The driving free energy for the phase transition can be taken as the deviation of Gibbs function value in the phase equilibrium state: for =e (6b) thermally activated law given by equation (2b): where ∆G β α is driving energy of the phase transition.
where E eq is the phase equilibrium potential.
In the particular case when diffusion is very fast compared to phase transition kinetics, the concentration of the guest species is constant (θ β (x=ξ)=θ β (x=0)) and equation ( 8) simplifies to: In the general case we must use equation (8).*Charge transfer kinetics: given the thermodynamic model, the charge transfer flux at the electrode-electrolyte interface can be written: (10) where E is the electrode potential, k 0 is the kinetic constant of charge transfer, C (M + ) is the guest ion concentration in the electrolyte at the electrode interface, and is the standard bulk concentration of (M + ) in the electrolyte.Now we have to find numerical solution of the problem represented by equations (4, 5, 6, 8 and 10) for a given triangular variation of the electrode potential with time.There is a major difficulty in the way of treating this problem by finite difference method.In fact, the position of the two-phase separation surface is not known in advance and it has to be determined in the course of computation.To overcome this difficulty, a number of authors have utilized the similarity transformation to solve similar problems [13][14][15].The singularity can be handled by immobilizing the moving boundary owing to the introduction of a new variable given by: η = x/ξ (11) In addition, during numerical simulation it is more convenient and advantageous to transform equations into dimensionless form before solving them.By combining the variable defined by Eq.( 11) and dimensionless parameters introduced in Table 1, Eqs. (4, 5, 6, 8 and 10) become expression as given in equation 12.This set of differential equations and associated boundary conditions can be solved numerically using finite differential methods [15,16].The implicit finite differential formulae based on the Crank-Nicolson discretization scheme [17] has been used to approximate Eqs.(4, 5, 6 and 10).For the discretization of Eq.( 8), we preferred the use of an explicit method, which permits to avoid hard complexity of computation happening if Crank-Nicolson approximation is used.In general, this process may lead to the loss of the stability insurance provided by the use of an implicit approximation.However, in the present case this is found to be without a measurable incidence provided discretization steps are not too small (>50).1)exp( E ) The successive steps of computation are as follows: 1-initialize all concentration and phase separation abscissa by setting them to zero at t * =0, 2-estimate ξ at t * =dt * from the explicit discretization of Eq.( 8), 3-if ξ ∈]0, e * [, then use the Choleski [16] algorithm to the tridiagonal system obtained by discretizing Eqs.(4, 5, 6a and 10).This leads to a concentration profile at t * =∆t * , 4-if ξ=e * then the phase transition is achieved and the electrode is constituted of a single homogeneous phase.So, (dε/dτ = 0) and the system ( 12) become similar to that of the classical diffusion phenomenon in a finite media.The concentration profile is also calculated by the Choleski method applied to the finite difference transformation of Eqs.(4, 5, 6b and 10).5-if ξ<0, this corresponds physically to the absence of phase β and the presence of only phase α, which is assumed to be identical to the initial structure.Then, there is no diffusion in the bulk of the electrode and the concentration profile remains flat and equal to zero.Similarly, if the solution at t * = n∆t * is known, then we can reach its value at t * =(n+1)∆t * using the successive computation steps described above.

Resulting calculated voltamograms
The characteristics of the theoretical voltamograms, obtained from the numerical solution of the system of Eqs.(12), depend upon several factors, namely phase transition rate, difference between the biphase equilibrium potential and the standard potential, charge transfer rate, and the electrode thickness.To obtain useful information and characterize this dependence, all these parameters were varied in the calculations.

Effect of phase transition rate
To characterize, solely, the effect of phase transition rate we carried out digital simulations by taking into account reversible charge transfer, fixing E eq * and e * and assuming different values for k ph * that correspond to a dimensionless ratio between the phase transition rate and diffusion coefficient.Typical simulated curves are shown in Fig. 2. On the one hand, if k ph * <2.10 -2 the diffusion is markedly faster than species transfer across the interface boundary and the process kinetics is limited by the latter phenomenon.This leads to a small charge transfer flux at the electrode-electrolyte interface.On the other hand, when k ph * increases, the rate of species transfer across the two-phase separation becomes progressively faster until it overruns the diffusion rate.Diffusion into the β phase becomes the limiting phenomenon for k ph * >1.In that case, voltamograms have the classical shape and characteristics obtained for intercalation in a perfectly non-stoichiometric and semi-infinite system without phase transition.Since several theoretical studies have been devoted to the latter case, we will focus on the process governed rather by species transfer across the two-phase separation.
In this case, the height of the cathodic peak is smaller than that of the anodic one (Fig. 3).This is directly due to the exponential form of the phase transition law, which means that to create a phase is more difficult than to make it disappear.Moreover, dimensionless current is shown to be proportional to k ph * .This implies the invariability of voltammograms when the scan rate varies.Analysis of a large number of calculated curves allowed us to relate the magnitude of the cathodic current peak to k ph * : However, the height and position of the anodic peak depend on the switching potential.Therefore, it is not possible to present a quantitative correlation for the anodic portion.Nevertheless, the anodic current peak is found to depend linearly on the phase transition rate k ph * .

Effect of phase equilibrium potential
Phase equilibrium potential is the upper limit of potential, below which the phase transition occurs.Its value expresses thermodynamic difficulty to create the new phase β.In fact, typical voltamograms obtained by varying E eq * and fixing all other parameters (Figure 4), show that the beginning of current flow may be delayed to lower potentials when E eq * is decreased.The total absence of current flow when E* is larger than E eq * is directly linked to the stability of phase α only above E eq * and to the assumption (A2) where we supposed the imperviousness of phase α to guest species.In reality, this latter assumption is not true and the experimental current flow is expected to start growing smoothly from the initial potential and it has a sharp slope change whenis reached.For E eq * >2, calculated curves show a cathodic peak attributed to the filling of the interfacial insertion sites.In contrast to diffusion phenomena, the cathodic current reaches a plateau after the peak.This is related to the occurrence of a steady insertion rate governed by the interface displacement.For E eq * <2, the cathodic peak disappears and the steady insertion rate is rapidly reached.Furthermore, the anodic portion is highly affected by the E eq * value.For E eq * <-2 the anodic portion presents a "stripping" peak.The abrupt decrease of the anodic current corresponds to the complete vanishing of inserted species, which can be oxidized, and to the time when the interface boundary reaches the electrodeelectrolyte interface.

Effect of electrode thickness
Under semi-infinite conditions, the electrode thickness is very large compared to the average distance traveled by the phase transition front, and the process takes place normally without any interruption.On the other hand, during insertion into a thin electrode, the phase transition front may reach the other face of the electrode.Thus, the electrode becomes constituted of a single homogenous phase and the transport process shifts from a phase transition-limited phenomenon to the classical diffusion-governed one.This leads to discontinuity of the current intensity observed in the calculated voltamograms and it corresponds to the change in the transport limiting process (Figure 5).The effect of electrode thickness can be shown more clearly using modeling of current response to a cathodic potential pulse below E eq * (Fig. 6).At first, the current response to such a pulse reaches rapidly a plateau that expresses the limitation of the process by the phase growth rate.After that, when phase transition is achieved, the current flow decreases progressively in agreement with the Cottrel's law.This leads to typical curves, which were experimentally observed during hydrogen insertion into palladium [4].

Effect of charge transfer kinetics
On the one hand, the height of the cathodic peak does not depend on the kinetic parameter (k 0* ) but it does depend on the value of the factor of symmetry α (Figure 7).On the other hand, the peak potential depends slightly on k 0* .This behavior is qualitatively similar to that encountered with quasi-reversible insertion in perfectly non-stoichiometric system.The anodic peak is more affected by the change of k 0* and α values.Its potential increases and its magnitude decreases when k 0* is decreased.Moreover, a decrease of the factor of symmetry leads to an increase of the height of the anodic peak and to a reduction of its half-height width.

Conclusions
The theoretical model proposed should be useful in any electrochemical process where a phase transformation plays a kinetic role.However, in the present state of the model, the application range is limited because, for mathematical reasons, the authors had to assume that one of the phases was pure.The generalization of the model to the case when both phases can have any composition is currently studied.
The calculations have in particular been done for the case of linear sweep voltametry, for thin or bulky electrodes.According to the values of parameters such as thermodynamic equilibrium potentials, diffusion and phase transformation kinetics and electrode dimensions, the calculated voltamograms present some very specific features, some of which have already been observed experimentally and reported in literature.
In a subsequent paper, the authors will present an experimental work in the view of validating this theoretical model for the particular case of sodium insertion into pyrocarbon at high temperatures.