The coupled thermal-electrochemical procedure is intended for the analysis of battery
electrochemistry applications that require solving simultaneously for temperature,
electric potentials in the solid electrodes, electric potential in the electrolyte,
concentration of ions in the electrolyte, and concentration in the solid particles used
in the electrodes.
A coupled thermal-electrochemical analysis:
is used for applications where the thermal, electrical, and ion concentration
fields affect each other strongly;
requires the use of coupled thermal-electrochemical elements;
allows for transient or steady-state solutions for temperature and ion
concentration and for steady-state solutions for solid and electrolyte
electric potentials;
can include the specification of a fraction of electrical and electrochemical
energy that is released as heat; and
The primary example of a battery electrochemistry application is the charging and
discharging of lithium-ion battery cells. During a charging cycle, the lithium ions are
extracted from the active particles of the positive electrode (cathode). The ions move
through the electrolyte by migration and diffusion from the positive electrode to the
negative electrode across the separator. At the negative electrode (anode), the ions
intercalate into the active particles. Heat is generated during the flow of electrons, ion
migration, and the intercalation process. During discharge, the cycle is reversed.
Coupled Thermal-Electrochemical Analysis
Rechargeable lithium-ion batteries are widely used in a variety of applications, including
portable electronic devices and electric vehicles. The performance of a battery
highly depends on the effects of repeated charging and discharging cycles, which can
cause the degradation of the battery capacity over time. The porous electrode theory
(Newman et al., 2004) is commonly
accepted as the leading method for modeling the charge-discharge behavior of
lithium-ion cells. The method is based on a homogenized Newman-type approach that
does not consider the details of the pore geometry. The porous electrode is based on
the concurrent solution of a highly coupled multiphysics-multiscale formulation:
At macroscale, the porous electrode is modeled as a homogenized medium
consisting of the superposition of an active solid electrode particle
phase and a liquid electrolyte solution phase, with known volume
fractions.
At microscale, a collection of spherical particles is assumed for which a (nonlinear)
lithium-ion diffusion model is solved. The solid particles are connected
by a conductive binder and together form the solid electrode phase.
You can use the coupled thermal-electrochemical procedure to model other battery
electrochemical processes that use the governing equations described in the sections
below. Although the discussion in the following sections is written in terms of
lithium cells, the theory is general and can be used for other cell types.
A full stack lithium-ion battery cell consists of an anode collector, porous anode, porous
separator, porous cathode, and cathode collector. The porous electrode can be
described as a material with a solid skeleton and a distribution of interconnected
pores that are filled with the electrolyte. The solid skeleton consists of a large
collection of active spherical particles that are in contact with the electrolyte.
The spherical particles are electrically conductive and are in electrical contact
among themselves, but they do not allow for interparticle diffusion. The spherical
particles are also chemically active and allow for intercalation and deintercalation
of lithium ions.
Figure 1 shows a typical jelly roll configuration of a cylindrical
battery. The porous parts of the battery are immersed in an electrolyte bath that
facilitates the movement of ions from the cathode to the anode during the charging
cycle (Figure 2). The porous electrodes are connected to nonporous
collectors that are connected to the external terminal of the battery. The external
terminals of the battery are connected to a power source during the charging cycle.
The electrons move in the solid phase, while the ions move in the liquid phase.
At the cathode solid-liquid interface, lithium ions can move out of the solid electrode
particles and enter the electrolyte through a deintercalation (or extraction)
process that can be written as:
Here, represents lithium in the solid electrode, which through the
deintercalation reaction, transforms into lithium ion in the electrolyte, , plus an electron, , and leaves behind a vacancy (or intercalation site) in the solid, . The reverse reaction corresponds to the intercalation of lithium
into the solid electrode.
During the charging cycle, deintercalated lithium ions at the cathode solid-liquid interface
move in the electrolyte through the separator and reach the anode. At the anode, the
lithium ions are intercalated (inserted) into the solid electrode particles. The
separator is porous and allows ions to pass through; however, it is an electrically
insulated material. Therefore, the electrons cannot pass through, which avoids the
possibility of a short circuit in the cell.
The governing equations used to model lithium-ion cells span both multiple scales and multiple
physics and involve solving a fully coupled thermal and electrochemistry problem.
The primary solution variables that are solved for at the macroscale (or continuum
scale) are the electric potential in the solid phase (), temperature (), electric potential in the electrolyte (), and the lithium ion concentration (). At microscale, the model solves for the concentration of lithium () in the particle.
Governing Equations
The porous electrode theory consists of five main parts:
Intercalation at the solid-liquid interface.
Conduction of electrons in the solid phase.
Diffusion and migration of lithium ions in the liquid phase.
Diffusion of lithium in the microscale particle.
Coupling of the continuum scale and microscale.
An optional sixth contribution is associated with Joule heating generation.
Intercalation at the Solid-Liquid Interface
The intercalation process at the solid-liquid interface is governed by electrochemical
kinetics and is defined by the Butler-Volmer equation. The Butler-Volmer equation
describes the relationship between the electric current density through the solid
electrode-electrolyte interface, , and the voltage difference between the electrode and electrolyte
according to:
where and are the exchange current density and overpotential, respectively; is the mean stress or hydrostatic pressure; and represents the partial molar volume of lithium in the electrode
material. The term captures the effects of mechanical stress on the electrical response of
the electrode material. The overpotential, , is defined as:
where is the solid electrolye interface resistance.
The exchange current density, , can be written as:
In the above equations,
is Faraday's constant;
is the gas constant;
is the charge number of the lithium ion
battery;
is temperature;
is the absolute zero temperature,
are the cathodic and anodic rate constants,
respectively;
are the cathodic and anodic transfer coefficients,
respectively;
is the theoretical maximum lithium capacity;
is the normalized surface lithium concentration in the
microscale particle;
is the reference value of lithium ion concentration in the
electrolyte; and
is the open circuit potential (OCP) as a function of and .
The normalized concentration, , is defined as , where is the lithium concentration in the particle. The normalized
surface concentration, , is defined as , where is the lithium concentration on the surface of the
particle.
Conduction of Electrons in the Solid Phase
The equilibrium equation for the current density, , in the solid phase is given by:
where
is the effective electrical conductivity of the solid phase
computed using the Bruggeman relationship for tortuosity (). Tortuosity represents an intrinsic property of a porous
material characterizing the ratio of the actual flow path length to the straight
distance between the ends.
For isotropic tortuosity:
For anisotropic tortuosity, the above expression can be written as:
where
is the electric conductivity of the bulk solid material (which can be a function
of normalized concentration);
is the tortuosity of the solid porous
electrode;
is the volume fraction of the solid phase;
is the unit tensor;
is the Bruggeman exponent for isotropic tortuosity;
are the Bruggeman exponents for anisotropic tortuosity; and
is the wetted particle surface area per unit volume, computed by default as , where is the outer radius of the spherical particle
associated with the microscale model.
Diffusion and Migration of Lithium Ions in the Liquid Phase
The current density in the electrolyte, , and the lithium ion concentration, , in the liquid phase are solved for using the following
governing equations:
where
In the above equations,
is the effective electrical conductivity in the electrolyte
phase,
is the effective ion diffusivity in the electrolyte phase,
is the transference number that defines the fraction of the
total electrical current that is carried by lithium ions in the
electrolyte, and
is the molar activity coefficient that accounts for
deviations from the ideal behavior in a mixture of chemical
substances.
and are computed using the Bruggeman relationship for
tortuosity.
For isotropic tortuosity:
For anisotropic tortuosity, the above expression can be written as:
where
are the electric conductivity and ion diffusivity of the
electrolyte, respectively;
is the unit tensor;
is the tortuosity of the liquid (electrolyte)
phase;
is the volume fraction of the liquid phase;
is the Bruggeman exponent for isotropic tortuosity; and
are the Bruggeman exponents for anisotropic tortuosity.
Diffusion of Lithium in the Microscale Particle
The active material at the electrodes is assumed to consist of
spherical microparticles. The particles are assumed to have electrical contact
between themselves but no interparticle diffusion. It is assumed that the
microparticle material is a good electronic conductor (transport number=1) and
has spherical symmetry. The conservation of lithium inside the particle is
governed by Fick’s second law of diffusion, written in spherical
coordinates:
where
is the radial coordinate of the spherical particle,
is the lithium concentration in the particle,
and
is the concentration-dependent diffusion constant of the
solid material.
The time evolution of lithium concentration in the particle is determined by diffusion and
the intercalation current density, , from the electrochemical model in the continuum scale. The
microscale particle is modeled as a sphere over which a one-dimensional finite
element mesh in spherical coordinates is created internally by Abaqus, and the diffusion problem is solved. You can specify the discretization for
the internal mesh and the type of meshing over the microscale particle (see
Particle Layer Mesh).
Coupling of the Continuum Scale and Microscale
The coupling of the microscale with the continuum scale happens through the Butler-Volmer
current density, , from the electrochemical model:
Heat Transfer and Joule Heating
The total heat generation in a battery is attributed to the contribution of five different
sources: flow of current at the solid-liquid interface, flow of current in the
solid phase, flow of current in the liquid phase, flow of ions, and entropy
generation. The amount of energy converted into heat from each of the source
terms can be scaled using the conversion factors, (=1–5):
where
is the ohmic loss at the solid-liquid
interface;
is the entropy generation, where is the derivative of the open circuit potential
with respect to temperature that you can specify in tabular
form;
is the ohmic loss in the electrode;
is the ohmic loss in the electrolyte; and
is the ohmic loss due to ion diffusion.
Fully Coupled Solution Scheme
The coupled thermal-electrochemical analysis in Abaqus uses an exact implementation of Newton’s method, leading to an unsymmetric
Jacobian matrix in the form:
Steady-State Analysis
Steady-state analysis provides the steady-state solution by neglecting the transient terms in
the continuum scale equations. It can be used to achieve a balanced initial state or
to assess conditions in the cell after a long storage period.
In the thermal equation, the internal energy term in the governing heat transfer
equation is omitted. Similarly, the transient term is omitted in the diffusion
equations for the lithium ion concentration in the electrolyte. Electrical transient
effects are not included in the equations because they are very rapid compared to
the characteristic times of thermal and diffusion effects. A steady-state analysis
has no effect on the microscale solution; the transient terms are always considered
in the solution of the lithium concentration in the solid particle.
Transient Analysis
In a transient analysis, the transient effects in the heat transfer and diffusion equations
are included in the solution. Electrical transient effects are always omitted
because they are very rapid compared to the characteristic times of thermal and mass
diffusion effects.
Initial Conditions
By default, the initial values of electric potential in the solid, temperature, electric
potential in the electrolyte, and ion concentration of all nodes are set to zero.
You can specify nonzero initial values for the primary solution variables (see Initial Conditions).
The typical set of initial conditions includes uniform but different values for the two
electrodes. The initial lithium concentration in the particles and solid and fluid
electric potentials in the electrodes are chosen such that the overpotential in the
electrodes is zero. The ion concentration in the electrolyte is typically assumed to
be uniform in the cell.
The initial condition for the microscale spherical particle is:
You can prescribe the following boundary conditions:
Electric potential in the solid, (degree of freedom 9).
Electric potential in the electrolyte, (degree of freedom 32).
Temperature, (degree of freedom 11).
Ion concentration in the electrolyte, (degree of freedom 33) at the nodes.
You can specify boundary conditions as functions of time by referring to amplitude curves.
A boundary without any prescribed boundary conditions corresponds to an insulated
(zero flux) surface.
The typical boundary condition consists of only grounding (setting to zero) the solid
electric potential at the anode. Thermal boundary conditions vary.
At the microscale, the boundary condition for the spherical particle is:
Abaqus applies the microscale boundary conditions automatically.
Loads
You can apply thermal, electrical, and electrochemical loads in a coupled
thermal-electrochemical analysis.
You can prescribe the following types of thermal loads (as described in Thermal Loads):
Concentrated heat flux.
Body flux and distributed surface flux.
Convective film and radiation conditions.
You can prescribe the following types of electrical loads on the solid (as described in Electromagnetic Loads):
Concentrated current.
Distributed surface current densities and body current densities.
You can prescribe the following types of electrical loads on the electrolyte (as
described in Electromagnetic Loads):
Concentrated current.
Distributed surface current densities and body current densities.
You can prescribe the following types of ion concentration loads (as described in Thermal Loads):
Concentrated flux.
Distributed body flux.
The typical loads include specification of a solid electric flux (current) at the
cathode. Thermal boundary conditions vary but typically include convective film on
the exterior surfaces. Customarily, no loads are applied on the concentrations and
electrolyte potential.
Predefined Fields
Predefined temperature fields are not allowed in coupled thermal-electrochemical analyses. You
can use boundary conditions to specify temperatures. You can specify other
predefined field variables in a coupled thermal-electrochemical analysis. These
values affect only field variable–dependent material properties.
Material Options
The thermal and electrical properties for both the solid and the electrolyte are active in
a coupled thermal-electrochemical analysis. All mechanical behavior material models (such as
elasticity and plasticity) are ignored in a coupled thermal-electrochemical analysis. The
electrochemistry framework requires that the material definition contain the complete
specification of properties required for the porous electrode theory, as described below and
in the sections that follow. In addition, the material name must begin with "ABQ_EChemPET_"
to enable the coupled micro-macro solution at the different electrodes. Special-purpose
parameter and property tables of type names starting with “ABQ_EChemPET_” are required in
these material definitions.
Thermal Material Properties
You must define thermal conductivity for the heat transfer portion of the analysis. In
addition, you must define the specific heat for transient problems. Thermal
expansion coefficients are not meaningful in a coupled thermal-electrochemical
analysis because the deformation of the structure is not considered. You can
specify internal heat generation.
Electrical Properties for the Solid Electrode Material
You must define the electric conductivity for the solid electrode portion of the analysis.
The electrical conductivity is defined as a function of the average normalized
concentration of the particle, (and, optionally, of temperature and user-defined field
variables).
Electrical and Ion Diffusion Properties for the Electrolyte Material
You define the electrolyte name and the ion charge number (z). The
charge number for the lithium ion is 1.0.
A label with the electrolyte name is used to identify all subsequent definitions required to
specify the electrical and ion diffusion properties of the electrolyte.
You define the electrical conductivity, diffusion coefficient, molar activity
coefficient, and the transference (migration). The electrical conductivity and diffusion
coefficients are functions of the lithium ion concentration and the electrolyte porosity
in the electrolyte. The molar activity coefficient and the transference are functions only
of the lithium ion concentration. In addition, all of the electrolyte properties can
depend on temperature and field variables. You can define Arrhenius temperature dependency
for the electrical conductivity and diffusivity. When you define temperature dependence
using a property table, the Arrhenius definition is ignored. For more details about the
Arrhenius dependency, see Arrhenius Temperature Dependency.
Defining Material Properties for the Porous Electrode and Microscale Particle
The following sections describe the material properties used to
define the multiscale nature of the electrodes and separator progressing from
the macroscale electrode level to the microscale particle level.
Electrode Definition
You must define various properties for the electrode at the macroscale, including a unique
region name; a region identifier characterizing the electrode as either an
anode, cathode, or separator; the volume fractions of the various phases in
the electrode; and the tortuosity factors in the three material directions.
You can specify a utilization fraction to account for inaccessible regions in
the electrode (see Ecker, 2015).
When lithium ions intercalate into a particle, the particle can swell
resulting in a convection of the electrolyte away from that region. If this
effect is significant and of interest to you, you can specify a factor of
one for the convection coefficient. By default, the value for the convection
coefficient is zero.
When particle swelling effects are included, the binder volume fraction, , and solid volume fraction, , are used to compute and update the electrolyte volume
fraction, .
Particle Definition
For each electrode you must define particle-specific properties such as the particle radius,
the initial concentration of the particle, and the contribution factor of
the particle when multiple particles are present. You can define multiple
particle types within the same electrode.
Particle Layer Definition
Each particle can have one layer or multiple concentric layers. All layers within a particle
must have the same material. For each layer, you must define a name for the
layer and a weight percentage per layer.
Particle Layer Mesh
The microscale particle is modeled as a sphere over which a one-dimensional finite element
mesh in spherical coordinates is created and the diffusion problem is
solved. You must specify the discretization () to use on each layer of the microscale particle. The
total number of elements in a particle can be between 1 and 100. Typically,
a mesh of 25 elements on a single particle gives good results. The meshing
can be uniform or biased toward the outside of the sphere using a quadratic
or equal volume approach.
Particle Layer Diffusion
You must define the theoretical maximum lithium capacity, , of each layer material in the particle and the diffusion
coefficient to use in the solution of the microscale lithium diffusion
equations. Several formats are available to define the diffusion coefficient
as a function of the normalized concentration, : tabular format, logarithmic tabular format, and spline
format.
For a tabular diffusion model, the diffusion coefficient is defined as a piecewise linear
function of the normalized particle layer concentration.
For the logarithmic definition of the layer diffusion, you specify the base 10 logarithm
of the diffusion coefficient as a tabular function of the normalized
particle layer concentration.
The spline input form for the diffusion curve requires the diffusion curve as a function
of the normalized concentration to be split into equal intervals
along the x-axis. Each interval is then curve
fit with a cubic spline equation. The parameter table has entries
ordered as the start and end of the interval range and the four
polynomial coefficients that correspond to the curve fit within that
range.
In addition, the diffusion coefficient can be a function of temperature. For
the spline format, only the Arrhenius form of temperature dependency is
supported.
Particle Swelling
The intercalation/deintercalation process of lithium atoms can lead to significant changes
in the volume of the spherical particles, which increases with the average
lithium concentration within the particle. Particle swelling also results in
the local displacement of the electrolyte (convection effects) to
accommodate the new particle volume. While elastic deformation is not
modeled explicitly by the thermal-electrochemical procedure, the effects of
particle swelling on the electrochemical process are accounted for by
considering its influence on the effective volume fractions of the solid and
liquid phases, as well as on the wetted particle surface area per unit
volume .
The volume change in the particle due to a change in average concentration is characterized
by a swelling coefficient, , as:
The solid volume fraction is defined as
where is the number of particles in the macroscopic volume, and is the particle initial volume.
Assuming that the local particle volume change is accommodated by displacing away the
electrolyte, with no overall macroscopic volume change, then
and, for the electrolyte
Similarly, the active surface area per volume, is also affected by swelling with
The electrochemical governing equations are modified accordingly to incorporate the effects
of particle swelling on , , and .
You can include the effects of particle swelling in the simulation by
specifying a swelling coefficient.
Butler-Volmer Definition
You must define the constants required to compute the current density using the
Butler-Volmer equation. You can choose whether the particle surface area per unit volume, , is computed within Abaqus or specified directly as a single value. By default, Abaqus computes the current exchange density, . You can specify in tabular form as a function of the normalized particle concentration
and temperature. is defined as follows for tabular definitions and does not include the
electrolyte concentration terms:
Specifying the Joule Heat Fraction in the Different Phases
You specify the contribution of the different phases to the overall joule heat fraction using
two parameter tables, one for the macroscale and one for the microscale. The
entropy table is used to define the derivative of the open circuit potential
with respect to temperature. This term is one of the five terms that contribute
to the overall joule heat fraction.
Arrhenius Temperature Dependency
An Arrhenius form of temperature dependence is also supported for the
definition of many of the material properties. The Arrhenius form of temperature
dependence provides a scaling of the corresponding material property given as:
where
is the activation energy,
is the universal gas constant,
is the absolute zero temperature,
is the reference temperature,
represents the value of material property at the reference
temperature, and
is the resulting value after accounting for Arrhenius
temperature dependency.
The Arrhenius form is ignored when you specify an explicit tabular temperature dependency in
the property table definitions for the different properties (using the
temperature parameter).
Universal Constants
You must define universal constants such as the Faraday number, gas constant, Avogadro's
number, Boltzmann number, absolute zero temperature, and the elementary charge
number.
Include File for Property and Parameter Table Definitions
The type definitions for the parameter and property tables described in the previous
sections are available as an include file (see Including Model or History Data from an External File). You can use the abaqus fetch
utility to obtain the file containing the type definitions of the parameter and property
tables used by the electrochemistry framework.
abaqus fetch job=ABQ_EchemPET_Types.inp
Solution-Dependent State Variables
The electrochemistry solution framework requires the use of solution-dependent state
variables at each integration point for the microscale problem at the anode and
cathode. Seventy-five solution-dependent state variables are required to
transfer information between the microscale and macroscale solutions. At the
microscale level, the number of solution-dependent state variables that are
required per particle is 822. Depending on the number of particles, , that are defined, the total number of state variables that
must be declared is . Only the anode and cathode regions of the battery require the
definition of solution-dependent state variables. Because the microscale problem
is not solved in the separator region, you are not required to define state
variables for the separator region of the battery.
Elements
The simultaneous solution in a coupled thermal-electrochemical analysis requires the use of
elements that have electric potential in the solid (degree of freedom 9),
temperature (degree of freedom 11), electric potential in the electrolyte (degree of
freedom 32), and ion concentration in the electrolyte (degree of freedom 33) as
nodal variables. The coupled thermal-electrochemical elements are available in Abaqus/Standard in three dimensions only (see Coupled Thermal-Electrochemical Elements).
Output
In addition to the output quantities available for the coupled thermal-electric procedure, you
can request the following output variables in a coupled thermal-electrochemical
analysis.
Nodal output variables:
EPOT
Electric potential in the solid phase.
EPOTE
Electric potential in the fluid
(electrolyte).
NNCE
Ion concentration in the
fluid (electrolyte).
RECURE
Reaction current in the
fluid (electrolyte).
RFLCE
Reaction ion
concentration in the fluid (electrolyte).
Element integration point variables:
ECDE
Magnitude and components
of electrical current density in the electrolyte.
MFLE
Ion flux rate in the
electrolyte.
CONCE
Ion concentration in the
electrolyte.
ELECPOTE
Electric potential in the
electrolyte.
ELECPOT
Electric potential in the
solid.
CONCGE
Magnitude and components
of ion concentration gradient in the electrolyte.
EPGE
Magnitude and components
of electrical potential gradient in the electrolyte.
For the microscale, the following subset of solution-dependent variables are most
relevant:
SDV1
Volume-averaged concentration on the surface of the microscale
based on the number of spherical particles defined in the
model.
SDV2
Volume-averaged Butler-Volmer current based on the number of
particles defined in the model.
Ebner, M., , and V. Wood, “Tool for Tortuosity Estimation in Lithium-Ion Battery Porous Electrodes,” Journal of the Electrochemical Society, vol. 162, no. 2, 2015.
Ecker, M., S. Kabitz, I. Laresgoiti, and D. Sauer, “Parameterization of a Physico-Chemical Model of a Lithium-Ion Battery II. Model Validation,” Journal of the Electrochemical Society, vol. 162, no. 9, 2015.
Ecker, M., T. K. D. Tran, P. Dechent, S. Kabitz, A. Warnecke, and D. Sauer, “Parameterization of a Physico-Chemical Model of a Lithium-Ion Battery I. Determination of Parameters,” Journal of the Electrochemical Society, vol. 162, no. 9, 2015.
Kumaresan, K., G. Sikha, and R E. White, “Thermal Model for a Li-Ion Cell,” Journal of the Electrochemical Society, vol. 155, no. 2, 2008.
Liu, B., X. Wang, H. S. Chen, S. Chen, H. Yang, J. Xu, H. Jiang, B. Liu, and D. N. Fang, “A Simultaneous Multiscale and Multiphysics Model and Numerical Implementation of a Core-Shell Model for Lithium-Ion Full-Cell Batteries,” Journal of Applied Mechanics, vol. 86, 2019.
Newman, J., , and K. E. Thomas-Alyea, “Electrochemical Systems,” Wiley-Interscience, Third Edition, 2004.
Wu, B., “Modeling and Design of Lithium-Ion Batteries: Mechanics and Electrochemistry,” Doctoral Dissertation, University of Michigan, 2019.