Explaining 5 equations of P2D model of Li-ion battery

某天在极度无聊中写了这篇1万7千多字的笔记。

豆瓣居然不能显示图片，有图片版本的：

http://note.youdao.com/noteshare?id=2abb6068e5a3f8cc40af1875e0a65d4f

The original equations are from Newman's book. I got the screenshot from Gregory Plett's book. I hope this note can help to understand P2D model of Li-ion battery.

A P2D model is pseduo-2D model of battery. There is a x dimension in the thickness direction and a r dimension for an assumed particle at each x location, representing all solid diffusion behavior of ions in the specific x location. in the x direction, the potential phi is different (phix), the Li-ion concentration of electrolyte is also different (cex), and the Li-ion concentration at each x is also dependent on current density.

This model can be described at 2 scales. Level 1 is called mirco-scale. Level 2 is called cotinuum scale.

At micro-scale, the cell is not space-averaged. Solid domain and liquid domian has clear boudaries. The equations are modeled in mutiple space dimensions, x,y, and z. There are five basic equations discribing dynamics and kinetics of the ions and charge. If (x,y,z) falls in solid domain, charge and mass conversation equations in solid domain is used. If (x,yz) falls in electrolyte doman, charge and mass conservation equaition in liquid domain is used. The interaction between solid and liquid is discribed by Bulter-Volmer kinetic equation. A summary of the five equations are listed below:

We now look at equation 1.

This is the diffusion equation of charge in solid particle. In solid prticle there is only electrons moving (ions cannot move), so the charge here means electrons.

Well, it is not entirely true that positive charge cannot move in solid. They can move inside the lettice of solid particles (such as LPFO), only from particle surface to particle center or in the oppositive direction.That cannot generate current flow (not in x direction)

IS here is current vector generated inside solid particles and between particles. It is equal to he gradient of potential (a scale) mutiply by electri conductivity of solid. The direction is from peak to valley, so there is a negative sign here. The divergence of IS is 0, meaning that all current vectors are parrellel in particle- no source, no sink.

Is there really no source? When a electron is combined with a ion inside a parricle,is that a sink? Now the eletron is not free electron. It is not moving. It seems to be a sink. Similarly, if a electron is seperated with a ion and moves towards the electric wire, it seems to be a source.

To understand this, we need to look at the equation 1a that derives equation 1:

As we mentioned earlier, ions can moved from liquid to particle or reversely (not between particles), so if there a ion in, and a electron in, that means IS is 0. That does not mean it is a sink. It just means current vector is 0.

Therefore, I need to understand that in the process of current flows from solid to collector to tab to wire to tab to colleor to solid, this equaion only describes electrons passing through a particle. In this context, we should notice that there a assumption that makes this valid, which is the accurate form of equation 1 (we call this equation 1b):

1a is the integral form. 1b is differential form of 1a. 1 is approximation of 1b. How does this approximation come out? We need to think it this way. Electrons moves extremely fast, and ions or atoms move slowly. In a fixed regon dV, we would assume rhov is constant, as the steady state is easy to reach.The assumption of rhov=constant is valid as long as time scale is larger than ms. We usually only look at the scale of second, so this will be a valid assumption.

We can understand this this way, the electrons are moving so fast, that ions and atoms cannot move or change within dt. In the scale of second, electrons move in and out really fast, and rhov is not impacted by passing of electrons, which means rhov= constant, drhov/dt=0. In this way, equation 1 is well explainined. (if you look at time scale smaller than ms, you see charge change due to electrons, you can keep the drhov/dt term. This is totally fine)

Equation 2 is very straightfoward: the rate of change of Li in solid domain is the divergence of flux vector of Li flow mutiple by a diffusion coefficient, and the Li flux vector here is the gradient of Li concentration. In a word, equation 2 means the rate of change of concentration is dependent on the spatial distribution of the concentration in solid.

The derivation of equation 2 is very similar to derivation of equation 1: first get the integral form and then use definition of divergence to derive the differential form.

Stating from equation 3, we need a lot of knowledge in electrocheimstry.

For equation 3, I am trying to explaining it in this way.

First I need to understand the structure of the equation. The left term is the rate of increase of concentraion of Li-ion in electrolyte. On the RHS, there are 3 terms. The first term is the diffusion term. This is very close to a diffusion equation like solid domain. But it is not the solid domain. The movement of Li-ion can not only be driven by concentration gradient. In liquid Li-ion can form ion current flow, so there is a second term called migration term. Of course, i-ion current flow is due to existence of electric field. In addition, liquid is moving due to density gradient, so there is a 3rd term called convetion term. This means ions are moving together due to the movement of solvent.

I don't need to remember all the detials of derivation. But I should remebmer the framework of the derivations and important assumptions.

The equation is derived from a theory called "concentrated solution theory". It is derived in 4 steps:

1. Examine net force of collisions on a species (Maxwell–Stefan);

2. Connect this force to the electrochemical potential of the species;

3. Then, show how this produces a flux of ions;

4. Finally, solve for mass balance via divergence of flux.

Now let's elaborate a litte bit on the 4 steps.

1. Examine net force of collisions on a species (Maxwell–Stefan);

In this step, derivation is for ideal gas as it is easier to understand. But the process is the same for ideal gas, liquid and polymer. The difference being that for ideal gas ideal gas law is used during the derivation, but for liquid we use the equation of electrochemical potentials to get the same results.

In fact, nothing important is derived except that we defined a K cofficient that desbrites the collision force between each pair 2 particle species in the electrolyte by assuming a relationship called Maxwell–Stefan relationsip:

12 is for species 1 and 2. Usually we will have 0 for neutral particles and + for positive particles and - for negative particles. F1,V means volume V-averaged force.

Here the force is like a drag force as it is velocity difference dependent. K is like a friction coefficient. K can be equivalently written as a Maxwell–Stefan diffusivity:

Here 12 means between particle 1 and 2. x is concentration ratio ci/cT, and p is just pressure being substituted by pV=nRT for ideal gas.

2. Connect this force to the electrochemical potential of the species;

We assume in another angle that

Before I continue, I will spend a little time to explain Gibbs free energy G. We know internal energy U describing random movement of molecules of a object. We know enthalphy H describing the energy of working fluid in a cycling system. Well H=U+pV literatlly means energy needed to create a system and and work needed to make space from the environment. G=U+pV-Tds means the energy needed to create a system and make a room for it subtracting energy that will spontaneously ben contributed from the ambient environment Tds via heat transfer. So it means all energy that a system in a specific environment can contribute.

The concept of G is not here, we just need to know how we use the equation above.

By making them equivalent. We come to the relationship that

12.1 will be used to writing gradient of mue into K and v in part 1 of step 3. That is all step 1 and 2 is doing.

Here electrochecmial potential mubar shows up. The next equation we will have chemical potential mu. The definition of mubar is dG/dn, so the unit is J/mol/kg. It is energy per unitl of matter. It includes all energy the matter has (internal energy, electric, gravitational, magneto...). Chemical potential is only the internal energy energy. So if we ignore gravitational and magnetic energy, we have electrochemcial potential=chemical potential+electric potential:

Here zi is charge number of species. F is sA/mol, and phi is V. In P39, we also need to use

The equatiopn above works because salt is neutral and ziFphi term will be canceled.

3. Then, show how this produces a flux of ions;

This means we need ot derive the equation 12.8. This is the most important step to do.

There are 3 terms that need to be studied seperately: the chemical potential mue, the electrolyte average diffusivity D, and transferance number t+0.

The first term mue is written as drag coefficient and velocity of +,- and 0 particles using 12.1 from step 2. This actually means gradient of chemical potential is equal to drag force of diffusion.

The second term is simply a definition constructed and eventuall expressed as drag coefficient and concentration:

Here nu is stoichiometric coefficient of positive and negative charge. c0 is concentration of neutral particles, and K is drag coefficient, reciprocal of diffusivity.

Transference number is an important concept as it showsup in the final form of the equation and is a model input parameter of P2D.

A transference number states the fraction of ionic current i carried by a certain ion when there is no gradient in chemical potential. I should be contributed by two parts, cation i+ and anion i-. There should be i=i++i-.

The transference number of the cation is proportional to the drag experienced by the anion:

This equation means that if the drag is larger for anion, then more current is contributed by cation, so t+0 is larger. This quite makes sense. t++t-=1. Typcially we assume t+0=0.4. This means in a current flow, the drag force experienced by anion is relatively smaller than cation, and 40% of the current is cation current.

We put all terms on the left together, and it gives RHS=c+v+=N+, and the equation 12.8 is proved.

4. Finally, solve for mass balance via divergence of flux.

There are 2 important things we do here. First, electrochemcial potential is an inconvinient term. We will use the definition of abosulute activity to replace it with concentration. With the equation below:,

we can tell that mu will be eventually replaced by concentration c, and that is exactly what we want.

During the derivation of the equation above, we need to use the 2 equations describing relationship between mubar and mu+ mu- and mue .

We will briefly explain the aboslute activity here. A Cation might attracts some anion surrouding it and form a aggregated particle group due to Coulomb forces. So is anion. Due to these groups, the movement of ions become more difficult. This means the ionic conductivity is lower. It feels like the concentration is lower. This is called "Debye–Hückel–Onsager theory" The definition is

Here 0<fi<1, ai=ci*fi is a effective concentraion or activity. ai- is the coefficient between absolute activity and activity. We don't need to know too much, because gamma will be hidden from the final form.

The way we eliminate gamma is by substitue the equation above into 12.8, and making a definiion:

Now this D is true diffusivity that is used in final form. greek D is weight averaged diffusivity made from Maxwell–Stefan diffusivity. In reality we measure D directly (remeber this is still micro-scale, in continum-scale we will do further average).

Another important thing we do here is we will use continuity equation

In fact we will use this equation twice, The first time we just use it to replace LHS of 12.8 so that we get

Note that here a divergence operator is added to N+, so each term of 12.8 has delta operator added.

For the second time we use it again

and we get

This helps to re-write the migration term

Other than that, we need to know that this equation 11.14 implies that charge can neither be stored nor created nor destroyed in a electrolyte solution. We also need to know that 11.14 is derived from mass continuity equation in solution (Newman 11.3)

After that we finally come from 12.8 to the final form and simplification:

This is the final form

This is simplification 1, by assuming dlnc0/dlinc=0. This assumption comes from Doyle.

This is another simplification, but realizing that most electrolyte (here we mean salt) is LiPF6. That gives Li+ and PF6-, which means v+=1, z+=1, even if it is not LiPF6, this can still be true. Also, the c above means Li+ in electrolyte, so it is ce, and here we come to the final form:

We can make a final explaination here. Rate of concentration change of Li+ (and so is LPF-) in electrolyte depends on spatial divergence of 3 terms: diffusion term - migration term - convection term. Diffusion term adds Li+ as concentraion difference drives Li+ in. Ionic current drives Li+ away, and the maginitude is dependent on vector product of the current flow and the gradient of transferance number of Li+. Convection of neutrual particles (solvent) aslo drives Li+ away through grandient of the solvent partilce velocity. The 3 terms here drives N+ flow, and the divergence of this flow is the increase of ce with respect to time.

An interesting question here is that if t+0=0.4, the migration term will be 0. We will discuss this later.

Equation 3 talks about mass conservation. There is a migration term involving ionic curernt ie, which will drive Li+ away to reduce ce. But we did not show how ie is governed. Equation 4 does this. We derive eq. 4 in 2 steps. The first step is to prove the following 12.27:

The second step is to replace mue with concentration. Also we need the divergence of i=0 equation, and then we can have Eq. 4 proved.

We first work to get 12.27. We need to understand what is s+. It is signed stoichimetry number. It is different from v because first it is signed, and second, it is the amount of specicies participating in a electrode reaction:

In contrast, nu is just for salt in electrolyte.

The derivationof 12.27 starts fron equation of energy balance:

And that simplifies to

Note that here is the introduction of the electrolyte potential, the driving force of the current flow i.

In the equation above there is mu+,mu- and mu0. The goal here to get 12.27 is to replace mu0 with mue, replac mu++mu- with mue and mu- and eventally replace mu- with mue and i.

To replace mu0 with mue, we need to use a "Gibbs–Duhem relationship" we have not mentioned before. That equation is simplied derived from dG=dU+d(pV)-Tds.

When T=const, p=const, only the right term is left, and we can replace mu0 with mue.

mu- is replaced with mue and i vector by using the follwing equation (P40) and 12.8.

This equation is from 12.1. 12.8 states the realationship between v0-v- , v+-v- and mue and i. This finally leads to 12.27. During the process a definition called eletrolyte conductivity is also defined as kappa:

This kappa has the unit of molar conductivity. It is important and is a function of concentration. It will be directly used as model input parameter in P2D modeling.

After we have 12.27, we replace mue with absolute activity, and we get the final form:

We can remove s by making the equation specific for Li-ion cell:

Negative electrode

Positive electrode

For both electrodes we have s+=-1, s-=0, s0=0 (no solvent involved in reaction), n=1.

Adding that divergence of ie=0, we have:

The f here is activity coefficient showing effective concentration a=f*c, and 0<f<1. Usually we lnf+-/lnce=0, and we can define an additonal conductivity kD which is related to kappa:

We then have final form of equation 4, which is really simpified:

Equation 5 and equation 2 are most important, because they are the only 2 equations needed in SPM.

Equation 1 and 2 are charge and mass PDEs inside solid (1 and 2 seem not to be connected, strange, right?). Equation 3 and 4 are mass and charge PDEs inside electrolyte. 4 describes how electric field phi creates the i in the migration term of equation 4. Equation 5 connects 12 and 34 by describing the movement of lithium between the interface.

The derivation is in 2 steps.

The first step is to get total electrode current density equation for 1 electrode. The 2nd step is to find the exchange current denstiy i0. The idea of the BV equation is that both reductoin and oxidization happen at the same time and i is the combination of the t terms. i0 here is equilibrium current density which means when iox=i0, ired=-i0, i=0. Specific format of i0 will be found in step 2.

The way to calculate iox and ired ix is based on reaction rate. See below, Reduction generates electrons, and oxidation generate positive charge.

Therefore,

The definition of reaction rate r is that the rate of overall product concentration change. Here we use one product case as an example:

Here k is reaction rate constant and is calculated from reference value usign Arrhenius equation.

By now we have connected i with r with activation energy. Next we need to connect activation energy with potential.

We do this based on a theory called "Activiated complex theory". It is easy to understand using the figure below in combination of the electrode reaction formula.

The idea is that for an electrode, both red and ox can happen, but the level of difficulty will change depending on if there is external potential added. If yes, activation energy will change, and tendency of reaction direction will change.

Pick the green curve for instance, there is now no external electric field. Both reaction can happen but to do that, the molecues need to overcome certain level of energy barrier (through thermaldynamic movement ). The probability is different. From left to right (reduction) it is more difficult because Ea,r0 is larger, and for oxidation it is easier to over come the energy barrier because Eao,0 is smaller phi0 is the reference potential.

Now we add external electric field here. The potential phi makes deltaG to be different, which changes Ea,r and Ea,o (see the blue curve). Now red is easier (left to right).

Having understand that, we need to connect activation energy with potential. Now we can calculate iox and ired:

We need to explain the charge transfer-coefficient here since it appears in the final form. It is defined as

From the plot above it is easy to see deltaEao, deltaEa,r is not shown but we can imagine that. Remember 0<alpha<1, so is 1-alpha. The concept of charge-transfer coefficient is that how easier it become for reduction (Li+ enters particle and become Li) to happen (r in Ear is reduction), or how difficult for oxidization to happen when you add external potential.

We play a few more tricks here. First we can define phi0=0. Note phi here means potential difference between solid and liquid, i.e. phis-phie. Second we think about a potential that makes the equilibrium to happen, i.e. phi=phirest. We also define a new potential called overpotential eta=phi-phirest. At equilibrium eta=0. (Sound familiar? In fact, phirest is OCP of the electrode!) This will lead to i0.

We can then have the intermidiate formart:

Next we look for i0 at equibrium. We do this by making iox=-ired, and eventually we get a format as:

For k0, we just need to know it is effective reaction rate constant. It is dependent on alpha and temperature (and material type, of course). So it can sitll be a constant of tempeature is not changed.

The above equation is for single product. We now look at equations for mutiple products:

We will make this specific for Li-ion cell. Now don't get confused here. cox means concentration of oxidant. It is the reactants for reduction reaction (refer to the reaction formulus above). The reactant for reduction includes ce, and available "holes" in the solid partilce for lithium, i.e. cs,max-cs,e. The reductant concentration is c,s,e. So We have the final i0 to be:

Eventually we get the BV equation for current denstiy at interface for any electrode.

People like to remove F by making j=i/F, which changes current flow into Li particle flow (unit becomes mol/m^2/s)

We will finally comment about Eq.5 in P2D and SPM. It is the current term in electrolyte charge transfer equation at the solid-electrolyte boundaries. It is also the boundaries for Equation 1 (phis and current)

Lastly, remember

To apply those 5 equations however, we need to make it from micro-scale to contimuum scale by doing a volume average. We also need some additional notes for OCP, capacity and SOC.

We will list them here, and discuss about it next time.

Charge conservation in solid:

Mass conservation in the solid (we don't need any volume avergae for it)

Mass conservation in the electrolyte

Charge conservation in electrolyte

Also nothing will change for Bulter-Volmer equation:

Here there are major 2 changes. First this is term as appearing again and again for all flux term j. This is

Another change is that all coefficient is now effective coefficient after volume averaing.It is defined using a coefficient called Bruggeman's coefficient (typically 1.5). For example:

OCP is measured directly.

Capacity:

SOC: