Thermodynamically Consistent Model for Space-Charge-Layer

Sep 4, 2015 - *Phone: +49 711 6862637. ... We derive a mathematical model for space-charge-layer formation in a solid electrolyte based on first princ...
0 downloads 0 Views 1MB Size
Article pubs.acs.org/JPCC

Thermodynamically Consistent Model for Space-Charge-Layer Formation in a Solid Electrolyte Stefanie Braun,*,†,‡ Chihiro Yada,§,⊥ and Arnulf Latz*,†,‡,∥ †

Institute of Engineering Thermodynamics, German Aerospace Center, Pfaffenwaldring 38-40, D-70569 Stuttgart, Germany Helmholtz Institute for Electrochemical Energy Storage Ulm, Helmholtzstraße 11, D-89081 Ulm, Germany § Toyota Motor Europe NV/SA, Hoge Wei 33A, B-1930 Zaventem, Belgium ∥ University of Ulm, Helmholtzstraße 11, D-89081 Ulm, Germany ‡

ABSTRACT: We derive a mathematical model for spacecharge-layer formation in a solid electrolyte based on first principles only. Consistent with the second law of thermodynamics, we employ mass, momentum, and energy conservation, supplemented with constitutive assumptions in the form of a Helmholtz free-energy functional. The resulting system of differential equations is solved semianalytically for a stationary 1D case, and the parametric dependencies of the space-charge layers forming at the boundaries under the influence of an external voltage are studied. We present results for different applied potentials, dielectric susceptibilities, and other parameters and compare our results with experiments. The predicted space-charge layers at the boundaries are in general not symmetric due to the restricted mobility of the anion lattice. Their size is found to be approximately 1 order of magnitude larger compared with liquid electrolytes, even if all macroscopic properties like mass density, dielectric constant, etc. are the same. Depending strongly on the dielectric properties of the material, typical widths of space-charge layers in some glass ceramics with very high dielectric susceptibility are predicted to be as large as several hundreds of nanometers, in qualitative agreement with experimental results in literature.

1. INTRODUCTION Batteries are likely to play a more and more prominent role within the world’s increasing efforts toward a sustainable energy supply based to a significant extent on renewables. Besides their use as storage devices for wind and solar power, batteries constitute an attractive candidate for automobile applications. With the first electric vehicles already on the market and in use, an increasing demand for greater energy storage densities and enhanced safety properties has inspired research to focus on developing new, improved battery concepts. Regarding safety issues, a particularly promising alternative to liquid electrolytes, which are often prone to thermal runaway, is to employ nonflammable solid electrolytes, i.e., fast ionic conductors. Other benefits of such solid-electrolyte-based systems include their favorable cycle and shelf lives1 and relatively low manufacturing costs. The most severe problem such systems have to face is their as yet rather limited power density, which is prescribed mainly to interfacial processes taking place at the electrode−electrolyte boundary, going along with a strongly increased interfacial resistance.2 The cause for this reduction of the effective ionic conductance is the formation of space-charge layers near the electrode, which can modify the dielectric properties of the material and are accompanied by strong electric fields impeding the transition of ions from the electrode to the electrolyte. A lot of important research has focused on the microscopic structure and properties of space-charge layers, in particular within the electrolyte itself at grain boundaries or © XXXX American Chemical Society

at artificial impurities, see, e.g., the numerous publications by Maier and co-workers.3−6 While those papers are mainly concerned with nanocrystalline structures, we want to derive the thermodynamically consistent transport equations for a homogeneous system to study the macroscopic implications of the space-charge layers at the electrodes. The inclusion of mesoscopic (not-point-like) impurities is not the focus of this paper, but can in principle be included later on via employing spatially inhomogeneous free-energy functionals. For homogeneous systems, continuum approaches are valid and have been employed successfully (and numerously) to describe liquid electrolytes (refs 7−9 as well as the review article by Bazant et al.10 and references therein). However, solid electrolytes differ from liquid electrolytes in several fundamental ways. Since solid electrolytes lack the additional neutral species constituted by the solvent in (diluted) liquid systems, the charge carriers are much more strongly concentrated in such systems. Furthermore, the fact that, in a “classical” ionic conductor (i.e., no electronic conduction), only one of the ion species (typically the cations) is mobile while the other species remains stationary changes the transport properties qualitatively. Existing theories for solid electrolytes typically model the ion migration in the solid by a hopping motion of the ions to nearby vacancies or Received: March 19, 2015 Revised: September 3, 2015

A

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C interstitials.11 However, the (few) existing treatments usually rely on the classical Nernst−Planck law, extended to include vacancies, which may have significant shortcomings with respect to the thermodynamic consistency.7 To avoid these inconsistencies, we carefully develop our model on the basis of first principles only, within the framework of rational thermodynamics.12−14 To this end, the conservation laws of continuum mechanics are augmented by constitutive relations, which are needed to close the model. These relations are derived from a free-energy functional that describes the properties of the system on the macroscopic scale. Further restrictions on the constitutive equations result from the requirement that the entropy inequality be fulfilled, and therefore, the thermodynamic consistency of the model is guaranteed by this approach. Finally, the electric properties of the system are described by Poisson’s equation. It is important to note that space-charge layers are intrinsically present in our model, in contrast to, e.g., Gouy−Chapman theory15,16 or Stern theory17 for fluids, in which the system has to be augmented by boundary-layer-like structures a posteriori. Consequently, a qualitative as well as quantitative analysis of the space-charge regions consistent with the adjacent bulk properties proves to be cumbersome in such a theory, whereas it is naturally possible within our new formalism. The remainder of the paper is organized as follows. In section 2, we introduce the mathematical model, including the constitutive assumptions. A stationary one-dimensional boundary value problem is solved semianalytically in section 3. In section 4, we evaluate this solution numerically and present results for different sets of parameters as well as a comparison with experiments. Finally, in section 5, we discuss the significance of these results and give a short summary in section 6.

itly. The continuity equations for the different particle species are19 ∂tρα + div(ρα vα) = 0 for α = a , c

(1)

and ∂tnα + div(nα vα) = 0 for α = v

(2)

where nα is the number density, ρα = mαnα the mass density, mα the atomic mass, and vα the velocity of particle species α. Since the anions are stationary, va = 0. The conservation equation for momentum can be written as14 ∂tρ v + div(ρ v ⊗ v) = div T + ρ b

(3)

with the total mass density ρ of the mixture and the barycentric velocity v defined as19 1 ρ ≡ ∑ ρα and v ≡ ∑ ρα vα ρ α (4) α respectively. T denotes the Cauchy stress tensor, and b represents the body forces acting on the electrolyte. Due to the conservation of angular momentum, the stress tensor must be symmetric. Note that, in a system under the influence of electromagnetic fields, the stress tensor consists of both mechanic and electromagnetic parts. 2.2. Poisson Equation. The electric behavior of the system is governed by the Poisson equation div D = n F

(5)

where n is the free charge density n ≡ ∑αzαnα and D = ϵ0E + P the electric displacement. Here, E denotes the electric field, P the electric polarization, ϵ0 the dielectric constant, and zα is the charge of particle species α. 2.3. Constitutive Relations. Constitutive relations for stress, chemical potentials, and polarizations are obtained from a free-energy function ψ, which describes the interactions and the entropy contributions. The nonequilibrium diffusive fluxes and electric current are derived from an entropy principle.7,14,18 We consider a solid electrolyte under isothermal conditions, which is subject to electric fields leading to a dielectric polarization P. For simplicity, we assume the susceptibility χ to be constant. Furthermore, the electrolyte is modeled as an isotropic elastic body with a constant bulk modulus K, where shear stresses are neglected. The corresponding free-energy functional reads7,11 F

F

2. MATHEMATICAL MODEL We consider a solid electrolyte clamped between two metal blocking electrodes to the left and right, which is subject to an external electric voltage. The electrolyte is modeled as a threecomponent system consisting of a stationary anion species (indexed a), mobile cations (indexed c), and vacancies (indexed v) in the form of Schottky-type defects. The vacancies do not carry any mass or charge. While the anions are assumed to be fixed on their respective lattice sites, the cations can diffuse through the system via hopping from vacancy to vacancy. Since, for pure defect-mediated transport, the sum of available defect sites and cations is constant, they can formally be described as an incompressible fluid, for which the sum of local volumetric concentrations of vacancies and cations is constant. At present, we do not include a distinction between defects of the anion sublattice and defects of the cation sublattice, which would lead to slightly different hopping rates regarding the respective defect sites. However, although relatively easy to include in the model, we do not expect this effect to play any significant role for the behavior of the cations. Under these assumptions, the equation of motion for the ions can be derived using rigorous methods from rational thermodynamics.7,13,14,18 The immobility of the anion lattice and the fact that vacancies are treated as massless particles lead to constitutive equations different from the one obtained in7 for liquid electrolytes. 2.1. Conservation Equations. Any proper thermodynamic model consists of the conservation equations for mass and momentum as well as the entropy inequality to ensure positive entropy production. Energy conservation is employed implic-

ρψ =

∑ nαψαR + ρψpol + ρψp + ρψmix α

(6)

where all quantities with a superscript R correspond to a reference configuration not subject to any external fields. The other contributions are 1 ρψpol = − ϵ0χE2 (7) 2 the polarization energy for a linear dielectric ⎛ n⎞ n ⎛n⎞ ρψp = (K − pR )⎜1 − R ⎟ + K R ln⎜ R ⎟ ⎝ ⎠ ⎝n ⎠ n n

(8)

the linear elastic energy and ⎛ nc nv ⎞ nn ρψmix = kBT ⎜nc ln + nv ln ⎟+α c v nc + nv nc + nv ⎠ nc + nv ⎝ (9) B

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C the entropy of mixing of the constituents. Since the anions are completely stationary, they do not take part in the mixing process. The parameter α is a measure for the “nonideality” of the mixing, with α = 0 corresponding to an ideal mixing behavior. From a more fundamental point of view, the term containing α leads to a shift of the point at which the energy minimum is realized. Monte Carlo simulations have shown that such a simple analytical ansatz is in very good agreement with the numerical calculations of the Helmholtz-free-energy for yttria-stabilized zirconia.20 On the basis of this free energy the expression for the stress tensor including the pressure is very similar to the one for liquid electrolytes.7 Important differences appear for the fluxes due to the massless nature of the vacancies and the immobility of the anion lattice. The requirement that the entropy production rate ξ ξ=−

1 ϑ

∑ (∇μα − zα E)· α

∂tρc

nF

(10)

(11)

=0 (12)

α

holds; i.e., only two of the three fluxes are independent. Different from liquid electrolytes the neutral species in solid electrolytes are massless, i.e., Jv = 0. This leads to the additional constraint Ja = −Jc. On the other hand, Jv/mv = nv(vv − v) is finite, and we can use the continuity equation to relate it to Jc/ mc via J ⎛ m n + nv ⎞ = − c ⎜1 + c c ⎟+∇×k mv mc ⎝ ma na ⎠

(17)

λ 2∂zzϕ = −L2n F

(18)

⎞ mc ⎛ n + nv ∂zμv ⎟ ⎜∂zμa + c ma ⎝ na ⎠

⎛ m ⎞ + ⎜zc − c za⎟∂zϕ ma ⎠ ⎝

(19)

The constitutive relations for the chemical potentials are simplified by the assumption of incompressibility of the anion lattice and nc + nv = const, which follows directly from the definition of Schottky defects. The incompressibility limit can be performed as in ref 7 to obtain

(13)

which must hold for some vector k. Since, in equilibrium, all flows must vanish simultaneously, we can argue that k = 0. Therefore, Jc is the only independent variable of the three mass fluxes. Even the center of mass velocity v can be expressed by the flux Jc since the anions are not mobile and the vacancies are massless. From ρv = ρcvc = ρcv + Jc, we obtain (14)

From eq 10, assuming that the fluxes depend linearly on the thermodynamic potentials, we conclude that ⎛ m m ⎞ = − M(∇μc − c ∇μa + ⎜zc − c za⎟∇ϕ mc ma ma ⎠ ⎝ Jc

⎛ m n + nv ⎞ − ⎜1 + c c ⎟∇ μ ) ma na ⎠ v ⎝

(16)

a 2∂zp = −n F∂zϕ

0 = ∂z(μc − μv ) −

Jv

Jc = (ρ − ρc )v

−ϵ0(1 + χ )Δϕ

3. ONE-DIMENSIONAL EQUILIBRIUM PROBLEM To understand the nature of the lithium deficient layer perpendicular to some planar electrodes, we consider a problem, in which the cations move only in this direction; i.e., lateral concentration gradients parallel to the electrodes may be neglected. We want to focus on the equilibrium solution, which is characterized by vanishing time derivatives and vanishing mass fluxes, Jc = 0. Scaling the spatial coordinate x by the total length L of the system as x = Lz with z ∈ [0,1], we normalize all quantities with regard to their respective reference values [i.e., nα ↔ nα/nR, p ↔ p/pR, L ↔ L/LR, zα ↔ zα/e0, K ↔ K/pR, ϕ ↔ e0ϕ/kBT, mα ↔ mα/(mc + ma)]. Note that the densities and masses are normalized to the total reference density nR = nRc + nRa + nRv and the total mass of the system, respectively. The set of equations now reads

Note that, following directly from the definition, the relation

∑ Jα

=

together with the equation for the mass flux 15.

always be positive leads to a condition on the nonconvective mass fluxes Jα = ρα (vα − v)

⎛ ρ ⎞ −div⎜⎜ Jc ⎟⎟ ⎝ ρ − ρc ⎠

−n F ∇ϕ = ∂t(ρ v) + div(ρ v ⊗ v) + ∇p

Jα mα

=

(15)

for some positive mobility M, where we expressed the electric field E by the electrical potential ϕ, E = −∇ϕ, and the chemical potentials are obtained from the free-energy ψ as partial derivative, μα = ∂(ρ ψ)/∂(nα). This relation closes our system. In the absence of a body force and after applying all constitutive relations including the relation between the center of mass velocity and the flux Jc, the full set of equations reads

⎛ nc ⎞ μc = ψc* + a 2K (p − 1) + ln⎜ ⎟ ⎝ nc + nv ⎠

(20)

⎛ nv ⎞2 +α⎜ ⎟ ⎝ nc + nv ⎠

(21)

μa = ψa* + a 2K (p − 1)

(22)

⎛ nv ⎞ μv = ψv* + a 2K (p − 1) + ln⎜ ⎟ ⎝ nc + nv ⎠

(23)

⎛ nc ⎞2 +α⎜ ⎟ ⎝ nc + nv ⎠

(24)

The independent variables of the system are ϕ, nc, and the pressure p. The two dimensionless parameters a and λ are defined as a2 = C

pR R

n kBT

λ2 =

kBT ϵ0(1 + χ ) e02n R (L R )2

(25)

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C ∂zϕ = F(n F)∂zn F

They define the scale of the pressure relative to the ideal gas pressure at some reference density nR and temperature T as well as the spatial scale relative to the length LR of the electrolyte. We see that the spatial scale λ behaves similarly to the Debye length in liquid electrolytes. λ2 is especially proportional to the dielectric permittivity ε = 1 + χ and inversely proportional to the density of charge carriers. The final set of equations to be solved as a result of the constitutive equations for μα and the conservation equations is 2

F

a ∂zp = −n ∂zϕ F

(note that ∂zn = zc∂znc), we obtain the expression γ α F (n F ) = 2 2 + ∑ F i νzc n + βi i = 1,2

−n = λ ∂zzϕ

(34)

with γ1 = −

(26)

2

(33)

F

1 zc

γ2 =

1 zc

β1 = zcna

β2 = −zc(1 − 2na) (35)

Here, we employed mc + ma = 1, and we assumed that za = −zc. Before we proceed, we define

(27)

⎛ ⎛n ⎞ n ⎞⎞ ⎛ m ⎞ ⎛ 0 = ∂z⎜ln⎜ c ⎟ + α⎜1 − 2 c ⎟⎟ + ⎜zc − c za⎟∂zϕ ⎝ ⎝ ⎝ν⎠ ν ⎠⎠ ⎝ ma ⎠

G (n F ) ≡

∫ F(nF) dnF = ν2zα2 nF + ∑ γi ln|nF + βi| c

2 ⎞⎞ m ⎛⎛ ν⎞ ν ⎛ ⎛ ν − nc ⎟⎞ ⎜⎛ nc ⎟⎞ ⎟⎟ − c ⎜⎜⎜1 + ⎟a 2∂zp + ∂z⎜⎜ln⎜ +α ma ⎝⎝ na ⎠ na ⎝ ⎝ ν ⎠ ⎝ ν ⎠ ⎟⎠⎟⎠

i

(36)

and H (n F ) ≡

(28)

F 2

∫ nFF(nF) dnF = ν2zα2 (n2)



c

with ν = nc + nv = const and L = 1. 3.1. Boundary Value Problem. We consider the situation of a solid electrolyte subject to an external voltage applied at the electrodes as depicted in Figure 1.

∑ γβi i ln|nF + βi | i

(37)

which will be used excessively in the following. Equation 33 can be integrated up to obtain an expression for ϕ(nF)

ϕ(n F ) = G(n F ) + cϕ

(38)

and with the boundary conditions in eq 29, this expression becomes ϕ(n F) = G(n F) +

∫0

(39)

(40)

which can be integrated up once to obtain (29)

∂zn F =

for some given values of ϕL at the left electrode and ϕR at the right electrode. The pressure is fixed at one boundary, p(z = 0) = p0. In addition, we impose electroneutrality on the system 1

G(x0) + G(x1) 2

λ 2(F(n F)∂zzn F + F ′(n F)(∂zn F)2 ) = −n F

The corresponding boundary conditions are ϕ(z = 1) = ϕR

2



with x0 ≡ nF(z = 0) and x1 ≡ nF(z = 1). The combination of eqs 27 and 33 leads to an equation for nF

Figure 1. Sketch of the solid electrolyte.

ϕ(z = 0) = ϕL

ϕL + ϕR

1 2 cnF − 2 H(n F) F λ F (n )

(41)

This equation is separable and is solved by x

n F (z ) d z = 0



(30)

Note that this condition does not enforce electroneutrality pointwise, but only on the electrolyte as a whole. Together with these conditions, the system in eqs 26−28 forms a well-posed problem. 3.2. Semianalytical Solution. Equation 26 can be integrated up easily to obtain an expression p(ϕ)

x0

cnF −

2 H (n F ) λ2

=z (42)

cnF

The integration constant can be determined from the constraint that z ∈ [0,1], leading to 2 = λ

1 λ2 (∂zϕ)2 + cp (31) 2 a2 where cp can be determined from p0 once ϕ(z) is known. Thus, the pressure can be eliminated from the system, and we only solve the remaining two equations for ϕ and nc. Equation 28 can be rewritten as p=

ν⎛ 1 α⎞ ∂zϕ = − ⎜ − 2 2⎟ zc ⎝ nc(ν − nc) ν ⎠

F (n F ) d n F

∫x

x1

F (x ) d x cnF ̃ − H (x )

0

(43)

cnF

which defines uniquely once x0 and x1, i.e., the charge densities at the electrodes, are known. These are set by the boundary conditions on ϕ and the electroneutrality constraint, which lead to the following system of equations

(32)

G(x1) − G(x0) = ϕR − ϕL

(44)

H(x1) − H(x0) = 0

(45)

The resulting values for x0 and x1 are uniquely defined for nF ∈ [−na, 1 − 2na] (corresponding to the situations of all cations

and using nF instead of nc and defining F(nF) via D

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

much smaller value λ ∼ 10−5. For a Li1+x+3zAlx(Ti,Ge)2−xSi3zP3−zO12 glass ceramic, which is the material closest to the one used for the experiments, we compare our theory with that for which we were able to obtain the relevant material parameters, and values for the dielectric susceptibility can become as high as 4.84 × 105. With a specific density of 3.05 g/ cm3 (material parameters kindly provided by Dr. Y. Inda, Ohara Inc.), the resulting value for λ is as large as ∼1.5 × 10−3. We will study the effect of varying λ in the following; for all other figures where λ is constant, we use λ = 1.5 × 10−3. Figure 2 shows the evolution of the normalized cation density nc/n in the solid electrolyte. As anticipated, the bulk

sites vacant and all vacant sites filled, respectively) and can be obtained numerically.

4. RESULTS AND COMPARISON WITH EXPERIMENTS We solve the system in eqs 44 and 45 as well as eq 43 numerically using Mathematica 10. Since the parameter λ, corresponding to what would be the Debye length divided by the system length in a fluid, is rather small, we expect the integration constant c̃nF to lie very close to the maximum of H(x), which is realized at x = 0. Therefore, we decompose c̃nF = ϵ + H(0) and solve for ϵ instead. In all calculations, we use the following set of parameters: zc = 1, za = −1, mc = 0.3, na = 0.4. Note that mc is normalized to the total mass of the system. The other parameters (Δϕ, α, λ) are varied in order to study the properties of the space-charge layers under different conditions. The parameters for our study were, as far as possible, set to match realistic experimental conditions. Unless indicated otherwise, in the following we will use a value for the cell voltage of Δϕ = 4.0 V. The parameter α is a measure for the “nonideality” of the mixing behavior of cations and vacancies. Since the focus of this work is to shed light on the general parametric behavior of the space-charge layers and not to analyze mixing properties, in order to avoid distortions of the results of interest we assume ideal mixing behavior (α = 0) except for the case where we study directly the α dependence. There is some question as to how large, or rather how small, the parameter λ (∼(1 + χ)/nR) should be. While it is no problem to determine the reference density of a given material, the choice of the appropriate permittivity is more subtle. From general linear response theory we know that in stationary situations the induced polarization charge is given by the zero frequency limit of the dielectric permittivity P = ε∞ + ε′(ω = 0) E. This treatment is valid, as long as there is a time scale separation between transport phenomena and relaxation phenomena which give rise to the frequency dependence of the dielectric permittivity. Since the dielectric excitations in the Nasicon structures are typically on the order of milliseconds or less (see the experimental results of Bucheli et al.21), this assumption is justified for the static case and for typical transport phenomena except for frequencies obtained in impedance spectroscopy. Contrary to the high frequency dielectric permittivity ε∞, the zero frequency limit of the real part of the frequency-dependent permittivity depends strongly on temperature and the number of mobile carriers. (e.g., ref 22). For the two fast-ion conductors La0.5Li0.5TiO3 perovskite and Li1.2Ti1.8Al0.2(PO4)3 Nasicon, Bucheli et al. measured the frequency dependence of the dielectric permittivity.21 Their results reveal a very strong variation of the dielectric permittivity from the order of 101 to 105 when going from high frequency to low at approximately room temperature, with the transition occurring at frequencies around 106−107 rad/s. This indicates that the reorientation of the microscopic dipoles, which is the origin for the increased dielectric response, is much faster than any of the transport processes occurring inside the battery, which happen on much longer time scales. Therefore, for the purpose of our model, the dipole reorientation can be viewed as occurring instantaneously, and consequently, the zero-frequency limit of the dielectric permittivity is the relevant value for determining the parameter λ. For La0.5Li0.5TiO3, with an average density of 4 g/cm3,23 a molar mass of 168.822 g/ mol, and a value for the permittivity of 105, we find that λ = 3.16 × 10−4. In the case of LiPO3 LIPON, the dielectric permittivity is several orders of magnitude lower, resulting in a

Figure 2. Fraction of cations vs normalized spatial coordinate z for λ = 1.5 × 10−3, na = 0.4, Δϕ = 4.0 V, α = 0.

concentration of the cations is very close to that of the anions, so that electroneutrality is maintained throughout almost the entire length of the electrolyte. Only near the electrodes do space-charge layers form to shield the applied potential difference, which are on the order of up to 400−600 nm for the glass ceramic we use. Note that we do not face the problem of overshooting cation density at the electrodes with systems relying on the Poisson−Boltzmann equation encounter,7 but the cation density stays within the physically allowed regime due to the coupling of the chemical potential to the pressure, which is dramatically enhanced by the strong electric fields in the space-charge layer. The corresponding electric potential distribution for an applied voltage difference of 4 V is shown in Figure 3. Clearly,

Figure 3. Electric potential ϕ vs normalized spatial coordinate z for λ = 1.5 × 10−3, na = 0.4, Δϕ = 4.0 V, α = 0. E

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C the range over which a potential variation occurs coincides with the range of the space-charge regions. In the bulk region, the electric potential is constant, since no currents are applied. In Figure 4, we compare the normalized cation densities for

Figure 5. Fraction of cations for different values of α and Δϕ = 4.0 V, na = 0.4, λ = 1.5 × 10−3: dashed, α = 0.0; dotted, α = 1.0; dash-dotted, α = 1.5; top, cathode; bottom, anode. Figure 4. Fraction of cations for different values of Δϕ and λ = 1.5 × 10−3, na = 0.4, α = 0: dashed, Δϕ = 2.0 V; dotted, Δϕ = 3.0 V; dashdotted, Δϕ = 4.0 V; top, cathode; bottom, anode.

different values of the external potential difference Δϕ. It is obvious that the larger Δϕ is, the larger the overall width of the space-charge layers is. Looking more closely at these regions, we see that they consist of a region where either all lattice sites are filled with ions (on the anode side) or the electrolyte is completely void of cations (on the cathode side). Between the fully filled or fully void regions and the bulk, there is an intermediate region where significant amounts of both cations and vacancies are present. In the bulk region, electroneutrality holds, and therefore, the number of vacancies approximately equals the number of defects found in an electrolyte not subject to an electric potential. Note that, with increasing Δϕ, only the fully stripped and fully filled regions of the space-charge layers increase in size, but the transition regions toward the bulk have the same thickness. In Figure 5, we study the normalized cation densities when the parameter α, corresponding to the “nonideality” of the mixing behavior of the components, is varied. It is easy to see that although the ratio between the fully filled/stripped region and the intermediate region is slightly different, there are no significant changes on the overall width of the space-charge layers. This also justifies our assumption of vanishing α in the preceding diagrams. Finally, Figure 6 shows the fraction of cations when the parameter λ is varied. Since the dielectric permittivity depends very strongly on the material as well as the operating conditions (remember also the discussion above) and can vary over several orders of magnitude, the parameter λ can vary significantly as well (λ is proportional to the dielectric

Figure 6. Fraction of cations for different values of λ and Δϕ = 4.0 V, na = 0.4, α = 0: dashed, λ = 1.5 × 10−5; dotted, λ = 1.5 × 10−4; dashdotted, λ = 1.5 × 10−3; top, cathode; bottom, anode.

permittivity). As we can see from the diagram, this may lead to substantial changes in the thickness of the space-charge layers. We observe that the widths of the space-charge regions are proportional to ∼10λ. It was found in literature that, for liquid F

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C

distances, since effects taking place on the atomic scale (e.g., layering effects) cannot be captured by such an approach. In the free energy which we employ, perfect isotropy is implicitly assumed; i.e., the ions can hop in every direction, which is a reasonable approximation on length scales above ≳10 Å. On smaller length scales, however, the hopping motion is severely restricted by the crystal structure and the existence of neighboring vacancies, which are not included in the model. Therefore, strictly speaking, our model is only valid for parameters resulting in space-charge layers above ≳10 Å. However, it has been demonstrated in other fields of physics that continuum-based models, while failing to provide exact quantitative results on small length scales, are still able to predict qualitative features as well as rough quantitative estimates correctly even well outside their “classical” range of validity if supplemented with appropriate constitutive relations.28,29 Therefore, even for smaller values of λ, we assume our model to still predict parametric dependencies accurately, although in such boundary cases one should not rely too heavily on exact numbers. For not too small values of λ resulting in space-charge layers larger than a few interatomic distances, to the best of our knowledge, all approximations we made are within their respective range of validity.

electrolytes, the space-charge region is typically of the order of λ.7 Therefore, relative to the liquid electrolyte, the space-charge layer is approximately 1 order of magnitude larger in a solid electrolyte. In contrast to the case of varying Δϕ, the length scale on which the charge adjusts itself from electroneutrality in the bulk to the nonzero values near the electrodes depends strongly on λ. Thus, not only the fully filled/stripped regions become larger when λ is larger but also the transition regions. Finally, we want to point out that the positive and negative space-charge layers are in general not symmetric due to the fixed anion lattice. The free-charge density nF varies on the interval [−na, nmax − na] = [−na, ν − na], and this interval is c only symmetric if ν = 2na, corresponding to a perfect crystal without any vacancies. Experiments resolving the potential and ion concentrations in solid electrolytes on the nanometer scale around the electrodes have proven difficult to carry out, and literature on this topic is scarce. Yamamoto et al. succeeded in visualizing the potential distribution in the vicinity of the electrodes of a Li1+x+yAlyTi2−ySixP3−xO12 ceramic.24 Their results suggest very large space-charge regions on the order of 1−2 μm. Spacecharge layer of this size seem to be very unusual, especially since there are experimental data for other solid electrolytes25,26 which indicate layers more on the order of 10 nm. Even within our theory sizes of 2 μm are difficult to justify. On the basis of the data of ref 21, a layer on the order of about 400−600 nm is feasible. To explain the remaining gap requires more detailed knowledge about the material properties of the solid electrolyte used in ref 24. We have seen that our model reacts very sensitively to changes in our dimensionless spatial scale parameter λ (which is proportional to the permittivity) and that the zero frequency limit of the dielectric permittivity is very sensitive to temperatures and material properties as, e.g., the availability of mobile charge carriers.22 In addition, also the frequency scale of dielectric relaxation can be very different for different solid electrolytes (compare, e.g., the data for Nasicon structures in ref 21 and for LiPO3 in ref 22). Differences in the size of the space-charge layers can therefore be explained within our theory by the different frequency-dependent dielectric response of solid electrolytes. It would be interesting to study this correlation further experimentally, but this is beyond the scope of this paper. Support for the correlation between dielectric properties and the size of space-charge layers comes also from experiments on dielectric interface modifications in solid electrolytes.26,27 In a forthcoming paper, we will show that our model is able to explain variation of space-charge layers and strong reduction of interfacial resistance by introducing dielectric modifications of the interface similar to the experimental findings in refs 26,27.

6. CONCLUSIONS We have presented a new, mathematically rigorous model for ionic transport in a solid electrolyte with special emphasis put on the thermodynamic consistency. We were able to identify a semianalytical solution for the one-dimensional stationary case, which we used to study the parametric dependencies of the cation distribution on different quantities. In equilibrium, the width of the space-charge region is mainly determined by the dielectric properties of the material, whereas other parameters, such as the external potential difference or the mixing behavior, only play a minor role. Typical widths for the space-charge layers of ionic conductors with high dielectric permittivity are found to be on the order of up to several hundreds of nanometers, which are much larger than values predicted by standard Gouy−Chapman theory (a few angstroms) for liquid electrolytes, but only slightly (a factor of 2−3) lower than values suggested by experiments for a similar material (up to a micrometer). If the discrepancy between our theory and said experiment can be explained by the different material parameters or if some other physical mechanism is responsible can only be clarified in a direct comparison with exactly matching conditions. The fact that it is possible to study the system semianalytically greatly eases the interpretation of the results in comparison with a full numerical simulation, in which parametric dependencies cannot be tracked as effortlessly as in an analytical context. Furthermore, our solution can provide an excellent benchmark for future simulations of the full, timedependent problem, which will be an important step toward understanding and improving the performance of all-solid-state batteries.

5. DISCUSSION The widths of the space-charge layers depend strongly on the conditions and the specific material parameters. For the set of parameters used in the diagrams above (corresponding to La0.5Li0.5TiO3), the space-charge layers we obtain from our model are of the order of ∼10−100 nm. For materials with different dielectric properties (smaller dielectric permittivity, resulting in smaller λ ≃ 1× 10−6), the resulting space-charge layers may be significantly smaller, down to the order of a few angstroms (in agreement with the Gouy−Chapman theory). Continuum models as the one underlying our calculation are typically only valid on length scales above a few interatomic



AUTHOR INFORMATION

Corresponding Authors

*Phone: +49 711 6862637. Fax: +49 731 5 03 4011. E-mail: [email protected]. *E-mail: [email protected]. G

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry C Present Address

synthesized by wet chemical methods and the effect of Ce, Zr substitution at Ti site. J. Mater. Sci. 2005, 40, 4737−4748. (24) Yamamoto, K.; Iriyama, Y.; Asaka, T.; Hirayama, T.; Fujita, H.; Fisher, C. a. J.; Nonaka, K.; Sugita, Y.; Ogumi, Z. Dynamic visualization of the electric potential in an all-solid-state rechargeable lithium battery. Angew. Chem., Int. Ed. 2010, 49, 4414−7. (25) Sata, N.; Eberman, K.; Eberl, K.; Maier, J. Mesoscopic fast ion conduction in nanometre-scale planar heterostructures. Nature 2000, 408, 946−949. (26) Takada, K.; Ohta, N.; Zhang, L.; Xu, X.; Hang, B. T.; Ohnishi, T.; Osada, M.; Sasaki, T. Interfacial phenomena in solid-state lithium battery with sulfide solid electrolyte. Solid State Ionics 2012, 225, 594− 597. (27) Yada, C.; Ohmori, A.; Ide, K.; Yamasaki, H.; Kato, T.; Saito, T.; Sagane, F.; Iriyama, Y. Dielectric modification of 5V-class cathodes for high-voltage all-solid-state lithium batteries. Adv. Energy Mater. 2014, 4, 1301416. (28) Cleri, F.; Philippot, S.; Wolf, D.; Yip, S. Atomistic simulations of materials fracture and the link between atomic and continuum length scales. J. Am. Ceram. Soc. 1998, 81, 501−516. (29) Ippolito, M.; Mattoni, A.; Colombo, L.; Pugno, N. Role of lattice discreteness on brittle fracture: Atomistic simulations versus analytical models. Phys. Rev. B: Condens. Matter Mater. Phys. 2006, 73, 104111.



Battery Research Division, Toyota Motor Cooperation, Mishuku 1200, Susono, Shizuoka, Japan.

Notes

The authors declare no competing financial interest.



REFERENCES

(1) Bates, J. B.; Dudney, N. J.; Neudecker, B.; Ueda, A.; Evans, C. D. Thin-film lithium and lithium-ion batteries. Solid State Ionics 2000, 135, 33−45. (2) Ohta, N.; Takada, K.; Zhang, L.; Ma, R.; Osada, M.; Sasaki, T. Enhancement of the high-rate capability of solid-state lithium batteries by nanoscale interfacial modification. Adv. Mater. 2006, 18, 2226− 2229. (3) Maier, J. Ionic conduction in space charge regions. Prog. Solid State Chem. 1995, 23, 171−263. (4) Maier, J. Thermodynamic aspects and morphology of nanostructured ion conductors aspects of nano-ionics part I. Solid State Ionics 2002, 154−155, 291−301. (5) Maier, J. Nano-sized mixed conductors (aspects of nano-ionics. Part III). Solid State Ionics 2002, 148, 367−374. (6) Jamnik, J.; Maier, J. Nanocrystallinity effects in lithium battery materials. Phys. Chem. Chem. Phys. 2003, 5, 5215. (7) Dreyer, W.; Guhlke, C.; Müller, R. Overcoming the shortcomings of the Nernst-Planck model. Phys. Chem. Chem. Phys. 2013, 15, 7075− 86. (8) Latz, A.; Zausch, J. Thermodynamic consistent transport theory of lithium-ion batteries. J. Power Sources 2011, 196, 3296−3302. (9) Less, G. B.; Seo, J. H.; Han, S.; Sastry, a. M.; Zausch, J.; Latz, A.; Schmidt, S.; Wieser, C.; Kehrwald, D.; Fell, S. Micro-scale modeling of Li-ion batteries: parameterization and validation. J. Electrochem. Soc. 2012, 159, A697−A704. (10) Bazant, M. Z.; Thornton, K.; Ajdari, A. Diffuse-charge dynamics in electrochemical systems. Phys. Rev. E 2004, 70, 021506. (11) Landstorfer, M.; Funken, S.; Jacob, T. An advanced model framework for solid electrolyte intercalation batteries. Phys. Chem. Chem. Phys. 2011, 13, 12817−12825. (12) Coleman, B. D.; Noll, W. The thermodynamics of elastic materials with heat conduction and viscosity. Arch. Ration. Mech. Anal. 1963, 13, 167−178. (13) Müller, I. Thermodynamics; Pitman Publishing Limited: London, 1985. (14) Kovetz, A. Electromagnetic Theory; Oxford University Press: Oxford, U.K., 2006. (15) Gouy, G. Constitution of the electric charge at the surface of an electrolyte. J. Phys. Theor. Appl. 1910, 9, 457. (16) Chapman, D. L. A contribution to the theory of electrocapillarity. Philos. Mag. 1913, 25, 475−481. (17) Stern, O. Zur Theorie der elektrolytischen Doppelschicht. Z. Elektrochem 1924, 30, 508−516. (18) Latz, A.; Zausch, J. Multiscale modeling of lithium ion batteries: thermal aspects. Beilstein J. Nanotechnol. 2015, 6, 987−1007. (19) deGroot, S.; Mazur, P. Non-Equilibrium Thermodynamics; NorthHolland: Amsterdam, 1969. (20) Horton, R.; Haslam, A.; Galindo, A.; Jackson, G.; Finnis, M. New methods for calculating the free energy of charged defects in solid electrolytes. J. Phys.: Condens. Matter 2013, 25, 395001. (21) Bucheli, W.; Arbi, K.; Sanz, J.; Nuzhnyy, D.; Kamba, S.; Várez, A.; Jimenez, R. Near constant loss regime in fast ionic conductors analyzed by impedance and NMR spectroscopies. Phys. Chem. Chem. Phys. 2014, 16, 15346−54. (22) Martin, S. W.; Yao, W.; Berg, K. Space charge polarization measurements as a method to determine the temperature dependence of the number density of mobile cations in ion conducting glasses. Z. Phys. Chem. 2009, 223, 1379−1393. (23) Jena, H.; Kutty, K. V. G.; Kutty, T. R. N. Studies on the ionic transport and structural investigations of La 0.5Li0.5TiO3 perovskite H

DOI: 10.1021/acs.jpcc.5b02679 J. Phys. Chem. C XXXX, XXX, XXX−XXX