COMSOL 6.4 - SEI Formation in a Lithium-Ion Battery

Side reactions and degradation processes may lead to a number of different undesirable effects causing capacity loss in lithium-ion batteries. Typically, aging occurs due to multiple complex phenomena and reactions that are occurring simultaneously at different places in the battery. The degradation rate varies during a load cycle, depending on potential, local concentration, temperature, and the direction of the current. Different cell materials age differently, and the combination of different materials may result in further acceleration of aging, for instance due to “cross-talk” between electrode materials.

This model demonstrates how to model aging in the negative graphite electrode in a lithium ion battery, where a parasitic solid-electrolyte-interphase (SEI) forming reaction results in irreversible loss of cyclable lithium. The model also includes the effect of increasing potential losses due to the resistance of the growing SEI film on the electrode particles, as well as the effect of reduced electrolyte volume fraction on the electrolyte charge transport.

Model Definition

Battery Chemistry and Aging reaction

The battery cell model is created using the Lithium-Ion Battery interface. This model uses the template model 1D Lithium-Ion Battery Model for the SEI Formation Tutorial, which contains the physics, geometry, and mesh of a lithium-ion battery. A more detailed description of how to set up this type of model can be found in the model example 1D Isothermal Lithium-Ion Battery. The template model 1D Lithium-Ion Battery Model for the SEI Formation Tutorial does not contain any SEI forming reactions or mechanisms. They are included in this model as described below.

In addition to the main graphite-lithium intercalation reaction on the negative electrode, the following parasitic lithium/solvent reduction reaction is also included in the model:

where S is the solvent (ethylene carbonate, EC) and PSEI is the product formed in the reaction. The production of PSEI results in loss of cyclable lithium in the battery, increase in the resistance of the SEI layer (Ref. 1 and Ref. 2), and in a reduction of the electrolyte volume fraction in the negative electrode.

The kinetic expression for the SEI forming reaction is based on a paper by Ekström and Lindbergh (Ref. 3). In this paper the SEI formation is assumed to be limited by kinetics in combination with a diffusion process through the formed SEI film, with the result that the aging slows upon thickening of the film. In addition, when the graphite electrode particles expand during Li intercalation into the negative electrode, the aging is accelerated due to “cracking” of the SEI film. The graphite expansion rate depends on both the state of charge and the intercalation current. The SEI forming reaction is assumed to be a reduction reaction, resulting in higher reaction rates for lower potentials (that is, lower battery state-of charge). The values of the model parameters were fitted, using a lumped zero-dimensional model, to experimental data obtained during cycling and calendar aging, for different state of charges, for a graphite/LFP cell, at 45°C. Only graphite aging effects were considered in the paper.

In this model example, the 0D model of the paper is expanded and applied to a 1D lithium-ion battery model using graphite/NCA electrodes. The kinetics of the parasitic reaction are described by the following kinetics expression for the local current density on the particle surface, iloc, SEI(SI unit: A/m2), in the negative graphite electrode:

(1)

Here

•

iloc,1C,ref (A/m2) is the local current density corresponding to a 1C discharge rate.

•

HK (1) is a dimensionless graphite expansion factor function (depends on the graphite state of charge). HK is zero during de-intercalation.

•

J (1) is the dimensionless exchange current density for the parasitic reaction.

•

α (1) is the transfer coefficient of the electrochemical reduction reaction.

•

ηSEI (V) is the overpotential, assuming an equilibrium potential of 0 V versus lithium.

•

qSEI (C/m2) is the local accumulated charge due to SEI formation.

•

f (1/s) is a lumped nondimensional parameter based on the properties of the SEI film.

The section of the Porous Electrode node is used to solve for an additional degree of freedom to keep track of the formed SEI concentration, cSEI (mol/m3), in the porous electrode according to

where νSEI is the stoichiometric coefficient of the SEI species in the reaction.

qSEI above is directly proportional to cSEI according to

(2)

where Av (1/m) is the electrode surface area.

Film resistance calculation

The thickness of the SEI layer, δfilm, is then calculated from the SEI concentration as

where MP (0.16 kg/mol) is the molar weight and ρP (1600 kg/m3) is the density of the product formed by the side reaction. The initial film thickness at t = 0, δfilm,0, is assumed to be 1 nm.

The resistance of the SEI layer, Rfilm (Ω·m2), used in the negative electrode, is then calculated from:

Load cycle, events, and Applied Current

The switching between the different stages of the load cycle is modeled by the event-based Charge-Discharge Cycling boundary feature.

The battery is cycled in the following sequence of operation modes:

Constant current charge at 1C until the cell voltage exceeds 4.1 V

Constant voltage charge at 4.1 V until the charge current drops below C/20

Constant current discharge at 1C until the cell voltage drops below 2.5 V

Time accelerating factor

Typically, a battery will need to be cycled many times to show any capacity loss. The incremental cycle-to-cycle differences in cycling behavior can therefore usually be assumed to be very small.

By assuming that every simulated charge-discharge-cycle in the model represents an average aging behavior for a larger number of cycles τ, and by further assuming that, over one complete charge-discharge cycle, all lithium captured in the SEI layer can be seen as stemming from the negative electrode, the capacity loss can be accelerated by rewriting the stoichiometry of the SEI forming reaction by adding the reaction formula

resulting in

τ can here be seen as a time accelerating factor, representing how many real cycles each simulated battery cycle should represent. In the model, τ is set to 100.

Postprocessing (plotting)

This model defines filter variables in the definitions, which are used to determine the onset of the discharge cycle (dch_start_filter), charge cycle (ch_start_filter), first cycle (first_cycle_filter), and last cycle (last_cycle_filter) for the load cycle applied during the charge-discharge of the battery. These variables are used while plotting the results after solving the model. See the postprocessing in the Modeling Instructions below.

Results and Discussion

Figure 1 shows the cell voltage during discharge for different cycle numbers.

Figure 1: Cell voltage during discharge.

Figure 2 and Figure 3 show the relative capacity versus time and cycle number, respectively. Both the capacity based on the total amount of cyclable lithium and the nominal 1C discharge capacity (based on the time spent during the 1C discharge part of the load cycle) decrease continuously, but with a higher capacity fade rate during the first cycles. Both capacities decay similarly, about 20% during the 2000 cycles of the study, indicating that the main contributor to the 1C discharge capacity fade is the loss of lithium, and not a significantly increased internal resistance due to film formation or worsened ion transport in the negative electrode.

Figure 2: Capacity versus total accumulated cycling time.

Figure 3: Capacity versus cycle number.

Figure 4 depicts the change in electrolyte volume fraction. Due to a, on average, higher SEI forming rate close to the separator, the reduction in the electrolyte volume fraction is more pronounced there. A similar trend is seen in Figure 5, which depicts the change in the potential drop over the SEI film. Due to a higher amount of formed SEI close to the separator, the increase in potential drop is higher there.

Figure 4: Electrolyte volume fraction versus cycle number.

Figure 5: SEI film potential drop versus cycle number.

Finally, the local states of charge in the electrodes at the boundaries facing the separator are shown in Figure 6.

Figure 6: Local state-of-charge on the separator-electrode boundaries.

References

1. P. Ramadass, B. Haran, P. Gomadam, R. White, and B. Popov, “Development of first principles capacity fade model for li-ion cells,” J. Electrochem. Soc., vol. 151, no. 2, pp. A196–A203, 2004.

2. G. Ning, R. White, and B. Popov, “A generalized cycle life model of rechargeable Li-ion batteries,” Electrochim. Acta, vol. 51, pp. 2012–2022, 2006.

3. H. Ekström and G. Lindbergh “A model for predicting capacity fade due to SEI formation in a commercial Graphite/LiFePO4 cell,” J. Electrochem. Soc., vol. 162, pp. A1003–A1007, 2015.

Battery_Design_Module/Lithium-Ion_Batteries,_Aging_and_Abuse/sei_formation

Modeling Instructions

Start this tutorial by opening a seed file that contains a 1D battery model, without any SEI formation reactions or mechanisms added.

Application Libraries

From the menu, choose .

In the window, select > > in the tree.

Click

Lithium-Ion Battery (liion)

Run this model in Constant Current-Constant Voltage (CCCV) mode using the Charge-Discharge cycling boundary node.

Charge-Discharge Cycling 1

In the window, expand the node.

Right-click > and choose > .

Select Boundary 4 only.

In the window for , locate the section.

From the list, choose .

In the Crate,dch text field, type -C_rate.

In the Vmin text field, type E_min.

Locate the section. From the list, choose .

In the Crate,ch text field, type C_rate.

In the Vmax text field, type E_max.

Select the checkbox.

In the Iupper text field, type I_upper.

Locate the section. From the list, choose .

Run the model at this stage, before adding any aging reactions, in order to check that the battery model cycles correctly.

Study 1

Step 2: Time Dependent

In the window, under click .

In the window for , locate the section.

In the text field, type range(0,180,10000).

In the window, click .

In the window for , locate the section.

Clear the checkbox. For this model we are not interested in the default plots.

In the toolbar, click

. The model should solve in a few seconds.

Results

Using the variables E_cell and I_cell, defined at the > > node, you can create a plot of the first load cycle as follows:

Load cycle

In the toolbar, click

In the window for , type Load cycle in the text field.

Click to expand the section. From the list, choose .

Locate the section. Select the checkbox.

Global 1

Right-click and choose .

In the window for , locate the section.

In the table, enter the following settings:


Expression	Unit	Description
I_cell	A	Cell current

Filter 1

Right-click and choose .

The seed model contained a number of definitions for filter variables, defined under the > > node. These variables have either the value 1 or 0, depending on the charging state of the battery. In the postprocessing part of this tutorial we will use these filter variables extensively to filter out different time ranges of the solution data.

In the window for , locate the section.

In the text field, type first_cycle_filter.

Global 2

In the window, under > right-click and choose .

In the window for , locate the section.

Select the checkbox.

Locate the section. In the table, enter the following settings:


Expression	Unit	Description
E_cell	V	Cell potential

In the toolbar, click

Note the time for the 1C discharge, ~3570s. This value will be used as reference in the analysis of the 1C discharge capacity loss later when the SEI formation is accounted for.

Lithium-Ion Battery (liion)

Now return to the battery model and add the needed features for defining the SEI formation reaction. Use the Dissolving-Depositing species functionality to define the SEI layer thickness on the negative electrode.

Porous Electrode 1

In the window, under > click .

In the window for , click to expand the section.

Click

In the table, enter the following settings:


Species	Density (kg/m^3)	Molar mass (kg/mol)
sei	rho_sei	M_sei

You can control if the volume change induced by a dissolving of depositing species should affect the electrolyte and/or electrode volume fractions. In this case we will assume the formed SEI reduces the electrolyte volume fraction only.

Clear the checkbox.

Also add a film resistance that depends on the deposited film thickness.

Click to expand the section. From the list, choose .

In the s0 text field, type dfilm_0.

From the Δs list, choose .

In the σfilm text field, type kappa_film.

Porous Electrode Reaction 2

Add a second porous electrode reaction on the negative electrode to account for the parasitic lithium/solvent SEI forming reduction reaction, and set the stoichiometric coefficient for the degradation reaction.

In the toolbar, click

and choose .

In the window for , locate the section.

From the Eeq list, choose . Locate the section. From the iloc,expr list, choose . In the associated text field, type I_SEI. The I_SEI variable was already defined in the seed file. You find the definition on the > > node.

Locate the section. In the νLiθ text field, type -(t_factor-1). The t_factor parameter is used to speed up the capacity fade per simulated cycle. You can read more about how the parameter is defined in the model documentation above.

In the table, enter the following settings:


Species	Stoichiometric coefficient (1)
sei	t_factor

Click to expand the section. From the list, choose . This is just a cosmetic setting to avoid the Materials node reporting missing material properties. The Heat of Reaction settings are not used in the model.

Definitions

Also add a domain integration operator with the default name intop1. It will be used during postprocessing to integrate the amount of cyclable lithium in the battery in order to calculate the remaining capacity.

Integration 1 (intop1)

In the toolbar, click

and choose .

In the window for , locate the section.

From the list, choose .

Study 1

The physics part of the model is now complete. Tweak the solver settings before computing.

Step 1: Current Distribution Initialization

Disable the aging reaction in the initialization study step.

In the window, under click .

In the window for , locate the section.

Select the checkbox.

In the tree, select > > > .

Right-click and choose .

Step 2: Time Dependent

In the window, click .

In the window for , locate the section.

In the text field, type range(0,180,t_cycling).

Solution 1 (sol1)

In the window, expand the > node.

Add a stop condition that makes use of a cycle number variable comp1.cycle_no defined at the > > node. This variable in turn is based on an internal variable defined by the Charge-Discharge Cycling node.

In the window, expand the node.

Right-click and choose .

In the window for , locate the section.

Click

In the table, enter the following settings:


Stop expression	Stop if	Active	Description
comp1.cycle_no>(no_cycles-0.5)	True (>=1)	√	Stop expression 1

Locate the section. From the list, choose .

By enabling storing the solution before and after events, the solution data for the first and last time step of the constant charge current, constant voltage and constant discharge current stages will be stored.

Clear the checkbox.

In the window, click .

In the window for , click to expand the section.

Change the nonlinear method to use a constant damping factor and a minimal jacobian update scheme. This makes the solver less robust, but faster.

From the list, choose .

In the text field, type 0.8.

In the toolbar, click

. The model should solve in about a minute or so.

Results

The following creates a comparison plot of the constant current discharge curves (Figure 1).

Discharge curve comparison

In the toolbar, click

In the window for , type Discharge curve comparison in the text field.

Locate the section. From the list, choose .

In the text area, type Color Legend: Number of cycles.

Locate the section.

Select the checkbox. In the associated text field, type Cell potential (V).

Global 1

Right-click and choose .

In the window for , locate the section.

In the table, enter the following settings:


Expression	Unit	Description
E_cell	V	Discharge voltage

Locate the section. From the list, choose .

In the text field, type (t-liion.cdc1.t_dch_start).

The Charge-Discharge Cycling node stores the latest time for switching to discharge mode in the liion.cdc1.t_dch_start variable. The above expression hence defines the time elapsed since the discharge started within each cycle.

Filter 1

Right-click and choose .

In the window for , locate the section.

In the text field, type dch_filter.

Locate the section. Clear the checkbox.

Color Expression 1

In the window, right-click and choose .

In the window for , locate the section.

In the text field, type cycle_no.

In the toolbar, click

SEI layer potential drop at 1C

The following creates a plot of the potential drop over the SEI layer, using the liion.deltaphi variable defined by the Film Resistance section of the Porous Electrode 1 node (Figure 5).

In the toolbar, click

In the window for , type SEI layer potential drop at 1C in the text field.

Locate the section. From the list, choose .

Locate the section.

Select the checkbox. In the associated text field, type Potential drop over SEI layer (V).

Point Graph 1

Right-click and choose .

Select Boundary 1 only.

In the window for , locate the section.

In the text field, type -liion.deltaphi.

Locate the section. From the list, choose .

In the text field, type cycle_no.

Select the checkbox. In the associated text field, type Cycle number.

Click to expand the section. Select the checkbox.

From the list, choose .

In the table, enter the following settings:


Legends
At negative electrode-current collector

Filter 1

Right-click and choose .

In the window for , locate the section.

In the text field, type dch_start_filter.

Point Graph 2

In the window, under > right-click and choose .

In the window for , locate the section.

Click

Select Boundary 2 only.

Locate the section. In the table, enter the following settings:


Legends
At negative electrode-separator

In the toolbar, click

Electrolyte volume fraction

In the window, right-click and choose .

The electrolyte volume fraction will decrease as a result of the formed SEI. Plot the volume fraction (Figure 4) as follows:

In the window for , type Electrolyte volume fraction in the text field.

Locate the section. In the text field, type Electrolyte volume fraction (1).

Point Graph 1

In the window, expand the node, then click .

In the window for , locate the section.

In the text field, type liion.epsl.

Point Graph 2

In the window, click .

In the window for , locate the section.

In the text field, type down(liion.epsl).

The down() operator indicates in this case that the value on the electrode side of the point should be used, not the default mean of the values on each side.

In the toolbar, click

Capacity vs. time

Now plot the capacity versus time (Figure 2). We will compute the remaining capacity in two ways: 1) based on the cyclable charge between 0% and 100% SOC and 2) based on the time elapsed during the 1C discharge.

In the toolbar, click