Modeling the impedance response of mixed-conducting thin film electrodes.

PCCP View Article Online

Published on 28 April 2014. Downloaded by Universitat Politècnica de València on 27/10/2014 12:43:18.

PAPER

Cite this: Phys. Chem. Chem. Phys., 2014, 16, 11573

View Journal | View Issue

Modeling the impedance response of mixed-conducting thin film electrodes Chi Chen,a Dengjie Chen,a William C. Chuehb and Francesco Ciucci*ac In this paper a novel numerical impedance model is developed for mixed-conducting thin films working as electrodes for solid oxide fuel cells. The relative importance of interfaces is considered by incorporating double layer contributions at the film/gas boundary. Simulations are performed on a model system,

Received 25th March 2014, Accepted 28th April 2014

namely doped ceria, in a symmetric cell configuration using geometrically well-defined patterned metal

DOI: 10.1039/c4cp01285b

can be extracted using the model. The impedance response depends strongly on the pattern spacing of

current collectors. Results reveal that experimentally consistent bulk impedances and surface capacitances the current collector, and is attributed to the electronic in-plane drift-diffusion as well as to the interplay

www.rsc.org/pccp

between the surface reaction resistance and the electronic/ionic bulk drift-diffusion resistance.

1 Introduction As electrodes of solid oxide fuel cells (SOFCs), mixed ionic– electronic conductors (MIECs) offer significant advantages over conventional composite ionic or electronic conductors, because the interfaces where reactions take place can extend beyond the triple-phase-boundary to the entire electrode|gas interface.1–3 Conventional porous MIEC electrodes have the advantage of a larger electrocatalytic area and consequently lower resistance. However, due to their complex morphology and microstructural features, it is difficult to map the measured electrochemical data using these electrodes to actual material properties. In particular, interpreting the mechanisms of electrochemical polarization poses a significant challenge. Patterned dense thin film electrodes can overcome this difficulty by excluding the geometrical uncertainties, thus allowing one to identify surface catalytic properties of the material as well as transport parameters.4–8 Electrochemical impedance spectroscopy (EIS)9 is one of the most commonly employed techniques for the characterization of MIEC thin films.10–17 The EIS responses of MIEC samples have been investigated on porous electrodes18–22 and on thin films,4,5,23 which even go down to the nanoscale.24 In a typical potentiostatic EIS experiment, a small sinusoidal voltage perturbation with an amplitude of a few mV drives the electrochemical system slightly beyond equilibrium, allowing the measurement of a small response current. From the amplitude ratio of the a

Department of Mechanical and Aerospace Engineering, The Hong Kong University of Science and Technology, Hong Kong SAR, China. E-mail: [email protected]; Tel: +852 2358 7187 b Department of Materials Science and Engineering, Stanford University, Stanford, USA c Department of Chemical and Biomolecular Engineering, The Hong Kong University of Science and Technology, Hong Kong SAR, China

This journal is © the Owner Societies 2014

excitation voltage and the measured current in addition to its phase angle difference, one obtains the impedance value at a certain frequency. This procedure is then repeated for a broad range of frequencies, leading to useful information that is typically obtained by interpreting the data using a suitable model. The tool of choice is a suitable equivalent circuit model (ECM).25–29 Aside from notable exceptions,25,27,30–34 the ECMs are the most commonly used models due to the ease of implementation. However, for given EIS data, the corresponding circuit may not be uniquely defined or even have physical meanings.35 For example, one specific arc in the Nyquist impedance plot may be related to more than one process.21 For complicated systems with multiple processes occurring simultaneously, the limitation of the ECM has to be considered. Compared to the ECMs, the numerical approaches offer the capability for virtual experiments, allowing one to assess readily the relative importance of different physical processes, such as surface reactions, bulk species transport and geometrical features, which are otherwise only assessable by advanced microfabrication and advanced probing techniques.36–39 An increasing number of studies have considered numerical modeling of the EIS response by solving partial differential equations (PDEs) via the finite element method (FEM).21,40–44 In previous studies, the EIS responses of MIEC thick samples were investigated numerically.45 The calculations were performed on a metal|MIEC|metal configuration. We stress that the thin film and the thick film configurations are significantly different, which would lead to quite different impedance responses even though the underlying PDEs describing the physical processes are similar.46 In thick samples, the resistance of the reaction and drift-diffusion dominates and the thickness-related area specific chemical capacitance is orders of magnitude higher than the interfacial capacitances.27 However, the latter may not be neglected in thin films since the contribution from the bulk chemical capacitance

Phys. Chem. Chem. Phys., 2014, 16, 11573--11583 | 11573

View Article Online


Paper

is small.47 Besides, in the presence of the electrolyte, which is the common case in SOFC application, the electron pathways are confined to the thin films and the electron transport may play an increasing role in this scenario. In this paper, an EIS numerical model is proposed which combines the solution of charge neutral Poisson–Nernst–Planck equations together with reaction boundary conditions. The model, in addition, incorporates a double layer treatment at the gas|film interface to study the peculiar characteristic behaviors of the thin films. The paper proceeds as follows: we first outline the impedance PDEs with suitable boundary conditions. Then, we solve the equations using FEM in the context of a symmetric cell configuration and we verify our approach by comparing the computed impedance and capacitance to the experimental data. Finally we study the impact of various geometrical and physical properties on the EIS responses, such as the current collector spacing and the surface reaction rate.

2 Cell geometry and assumptions The model is developed with the aim of understanding the electrochemical processes within the MIEC electrodes, specifically by linking them to the impedance responses. Our computations focus on the numerical analysis of a symmetric cell in the configuration shown in Fig. 1(a), where a dense thin MIEC layer is deposited on both sides of a solid electrolyte slab and a current collecting pattern of metal is present on top of the MIEC layers. The reactions occurring on the surface are schematically depicted in Fig. 1(b), where the entire surface provides reactive sites due to the mixed conducting properties of the film. The MIEC layer is thin, with the

PCCP

thickness ranging from 0.1 to 10 mm, while the electrolyte is on the order of 1 mm. The symmetric cell we consider here is put in a reducing atmosphere (H2–H2O–Ar) and at an intermediate temperature range from 550 1C to 650 1C. The thin film therefore works as an anode under similar working conditions as described in the experimental work by Chueh and Haile.16 Within the MIEC thin film layer, two mobile charge carriers are considered: the oxygen vacancies, denoted by the subscript ion, and the electrons, denoted by the subscript eon. Other charged dopants are assumed to be immobile and are named background charges. We assume that the background charge possesses one negative charge per particle and is characterized by a concentration B. We assume, in addition, that the oxygen vacancies flow perpendicularly to the MIEC|electrolyte boundary. From the above analysis, we show the electron and ion path lines we expect in the 2D cross-sectional system in Fig. 1(c). Moreover, the metal collector plays no role in the electrochemical reactions, as a general requirement in experiments for assessing the material properties. The interfaces MIEC|metal and MIEC| electrolyte are assumed to be both ideal. In light of the above analysis, we choose a model satisfying a number of assumptions. First, the YSZ slab behaviour is approximated to be onedimensional and the interface between YSZ and the MIEC under small bias conditions is reversible. This entails that when the system is perturbed slightly, as assumed here, the linear approximation holds and thus we can fix the potential to a constant. Second, the metal|MIEC interface is assumed to be reversible at a constant electric potential. Third, we assume that charge neutrality holds in the bulk material with the exception at the gas|MIEC boundary, where the capacitive responses that pertain to the gasinduced double layer may not be necessarily neglected.48

Fig. 1 (a) The symmetric cell configuration. A thick electrolyte is sandwiched by two thin MIEC layers, with metal current collectors placed atop. (b) Schematic depiction of the reactions occurring at the surface. (c) A cross-section of the cell symmetric unit, with intuitive electronic and ionic pathways within the material. (d) A cross sectional schematic of the corresponding 2D computation system.

11574 | Phys. Chem. Chem. Phys., 2014, 16, 11573--11583


View Article Online

PCCP

Paper

3 Numerical model and method


3.1 Governing equations from charge neutral Poisson–Nernst–Planck equations The electrons and vacancies are charge carriers, and their movements are subjected to the concentration gradient and also the external electric field. To incorporate the effects, the transport of charged vacancies and electrons can be described by the Poisson–Nernst–Planck equation set.26 Based on different physical processes that occur in the film, we divide the film into two regions: the first region is the bulk of the material that is away from the surface where we assume that the local charge is zero (i.e., the electroneutrality condition holds there); the other is the region near the surface where an electrical double layer is present and the non-zero local charge may occur. The species in the bulk material are subjected to a small sinusoidal perturbation as a consequence of the applied voltage. Under perturbation, the relevant physical quantities, such as the concentration of species and the electric potential, may be expressed as their equilibrium values plus a small perturbation correction. Consistent with previously published work,45 the equilibrium values are marked with a superscript (0) and the perturbations with a superscript (1), highlighting that only perturbations up to the first order are considered. The equilibrium concentrations of the vacancies and (0) electrons may be expressed as c(0) ion and ceon, respectively, and the perturbation counterparts are denoted by c(1) and c(1) , where the ion eon subscripts are chosen to be consistent with the work by Jamnik and Maier.25 Since the impedance is directly determined in the Fourier space, in the bulk region, the transport governing equations can be reduced to the following set of complex-valued PDEs by taking the Fourier transform:45 ˆ(1) Dxñ ˆ(1) = 0 iotn*n

(1a)

ˆ(1) Dx˜f^(1) = 0 iotf*n

(1b)

where i is the imaginary number, o is the angular frequency, tn* = (tn + nt % p/(4%p))/(1 + n/(4% % p)) and tf* = (tp tn)/(1 + 4%p/n) % are time constants from the derivation, with tn = lc2/Deon, tp = lc2/ Dion, n% = c(0) % = c(0) eon/B and p ion/B. In these formulae, lc is the characteristic length scale, Deon and Dion are diffusivities of the electrons and ions respectively. The variables n(1) = c(1) /c(0) , eon eon (1) (1) (0) (1) (1) p = cion/cion and f = ef0 /(kBT) are the normalized perturbed concentrations of electrons and ions, and electric potential respectively. The hat symbol in the equations denoting the ˆ(1) is the Fourier transformed n(1). kB is Fourier transform, e.g., n the Boltzmann constant, e is the elementary charge and T is the temperature. In addition, the charge neutrality condition gives ˆ(1) = n ˆ(1)n/(2% p % p). The reduced electrochemical potentials45 of the vacancies and electrons are computed by adopting the dilute solution ˆ(1)n/(2% approximation, which says (^ m(1) )* = n % p) + 2f^(1) and ion (1) (1) (1) ^ ˆ f . By solving the PDE set with suitable (^ meon)* = n boundary conditions given below, we can obtain the solution for the charge neutral model in the bulk region, which we


ˆ(1) ˆ(1) represent as (f^(1) EN, nEN, pEN), where the subscript EN is added to emphasize that the solution supposes charge neutrality. The equation we derived can be directly used to solve the bulk system since it is based on the charge neutrality assumption, however, in the region near the surface, charge neutrality may no longer hold due to the charge accumulation and eqn (1) needs to be corrected. One approach uses a linear superposition to represent the quantities in this region, which is expressed as: ^(1) f^(1) = f^(1) EN + fDL

(2a)

ˆ(1) ˆ(1) ˆ(1) = n n EN + nDL

(2b)

(1)

ˆ p

=

ˆ(1) p EN

+

ˆ(1) p DL

(2c)

It is important to note that the correction terms of the double layer (subscript DL) have been superimposed to the charge neutral solutions (subscript EN). We show in Appendix that the singular perturbation method can lead to analytical expressions for a capacitive current, which adds to the total current at the surface. The corresponding boundary regions labelled in Fig. 1(d) are discussed in more detail in the following section. 3.2

Boundary conditions

3.2.1 Symmetric and ionic/electronic blocking boundary conditions. The boundaries G2 and G3 are symmetric boundaries, implying that no flux may occur at the boundaries, ˆ(1)|G2,G3 = @ x˜f^(1)|G2,G3 = 0. leading to @ xñ At the boundary G1, the electrons are blocked but the ions can flow through. Consequently, for the electrons, the normal gradient of the corresponding chemical potential vanishes, ˆ(1) @ ˜yf^(1))|G1 = 0. which gives (@ ˜yn Since the ionic flow across G1 is reversible, we may assign the electric potential to a reference value f^(1)|G1 = 0. Alternatively we can pin the electrochemical potential of the electrons and ions simultaneously.45 We assume that ions cannot penetrate the metal|MIEC interface. Thus, the boundary conditions may be set to ˆ(1)n/(2% (@ ˜yn % p) + 2@ ˜yf^(1))|G4 = 0 and f^(1)|G4 = 1. 3.2.2 Reaction at the MIEC|gas interface. At the MIEC|gas interface, denoted by G5 in Fig. 1(d), gaseous hydrogen is adsorbed on the material and subsequently reacts with oxygen anions. The details of this mechanism are rather complicated.1,49 In this work we treat the surface reaction as a single step reaction with the reaction formula described by the following mechanism:45 0 H2 ðgasÞ þ O O Ð H2 OðgasÞ þ VO þ 2e

(3)

The products of the forward reaction, i.e., the ions and electrons, are injected into the system through the reaction boundary. Under the small perturbation assumption, the perturbed production rates for the electrons and ions in the Fourier space may be calculated based on a linear reaction (0) ˜ ˆ(1) and o ^ Sion = ^ Seon = 2kf (1 + c(0) kinetics, leading to45 o eon/4cion)pH2n S ^ eon/2 for the electrons and ions, respectively, where o . kf ¼ 2Dion k~0 p~b lc , and the pressures are normalized by the f O2


View Article Online

Paper

PCCP

standard condition, 1 atm. In the following, we define k~f ¼ k~0 p~b # particles per m3 for convenience and choose b = 0.25

Table 1

Parameter

Definition

Default value

to match experimental observations. 3.2.3 Double layer correction at the reaction surface. One may apply singular perturbations in order to capture the double layer effect at the surface. This leads to the following boundary condition at the gas|MIEC interface G5,

W1 W2 lc B

Current collector width MIEC exposed width Characteristic length Background charge concentration Temperature Normalized reaction constant Normalized oxygen partial pressure

1.5 mm 2.5 mm 10 mm 3.47 1027 # particle per m3 650 1C 1031 2 1024


f O2

@ n^ð1Þ ¼ An n^ð1Þ þ Bn @ y~

(4a)

ð1Þ

^ @f ¼ Af n^ð1Þ þ Bf @ y~

(4b)

where the relevant parameters are shown in Appendix. 3.3

Computational methods

The PDEs (1) with boundary conditions described above were solved using an h-adapted FEM on an anisotropic triangular mesh. The equations were solved using the free and open source package FreeFem++.50 The computational mesh was refined 11 times for each frequency. The method ensures that the meshes are finer in the regions where the gradients are steeper and coarser in other areas, representing faithfully the exact solution. In the computation, we choose 15% samarium doped ceria as the MIEC material,51 because this material has been studied frequently and published experimental data are available.

4 Results and discussion 4.1

Model verification

4.1.1 Interpreting experimental impedance data. The impedance values are computed by calculating the ratio of electrochemical potential drop of electrons across the film and the current. In Fourier space, the perturbed current density ð1Þ the film is computed using45 j^cross-plane ¼ . ð0Þ Ð ^ð1Þ Deon eceon G4 r m ey dx ðW1 þ W2 Þ, where ey is the outward eon

T ˜k0 f ˜ O2 p

Data for the domain geometry and default parameters

circuit featuring a resistor and a capacitor (RC) in parallel. As indicated in the figure, increasing oxygen partial pressures leads to larger impedances. This can be explained by the fact that in the reducing atmosphere, as in our case, the concentration of the electrons within the film scales linearly with 1=4 p~O2 .16 With fewer electrons, the current reduces, resulting in

the larger impedance. In addition, for higher oxygen content, the reaction of producing oxygen vacancies, eqn (3), is not as fast as in the lower oxygen content environment, which is another aspect that leads to the increase of impedances. As shown in the figure, given the same operating conditions, the numerical data capture the experimental response well, suggesting the validity of the impedance model. 4.1.2 The interplay between surface and chemical capacitances. The typical Nyquist plot of the impedance, as shown previously, exhibits an ideal semicircle. The overall capacitance of the cell can be extracted on the basis of RC responses, given the characteristic time scale. Two capacitive contributions are expected: one is from the bulk chemical capacitance and the other from the surface contribution. The bulk chemical capacitance26,27 is associated with the local concentration of species and is an extensive quantity that scales with the volume. The surface capacitance Csurf, on the other hand, does not depend on the volume, but rather on the active area surface extension.

across

unit normal of the boundary, W1 and W2 are the widths of the current collector and the exposed area in the 2D domain, respectively, as shown in Fig. 1(d). The electrochemical ð1Þ potential drop of the electrons across the MIEC film is V^ ¼ . Ð . Ð UT G4 m ^ð1Þ dx W1 G1 m ^ð1Þ dx ðW1 þ W2 Þ , where eon eon

UT = kBT/e is the thermal voltage. The film impedance is therefore calculated using the following expression: Z¼

ð1Þ V^ ð1Þ j^cross-plane

(5)

Unless otherwise specified, the parameters are provided in Table 1. The material properties and the parameters that we use in the simulation have been provided in previous work.45 Numerical impedances with experimental data are shown in Fig. 2, where the semicircles indicate the same responses of a


Fig. 2 Numerical impedance responses (solid line with k˜0f = 6 1030) of a 195 nm-thick thin film are compared to the experimental data (markers) at 650 1C. The experimental impedance data are taken from the paper written by Chueh and Haile.16


View Article Online


PCCP

In the derivation of our model, the surface capacitance is incorporated by a double layer correction, which results in an expression that contains the surface constant C0, a quantity introduced for simplifying the notation, as we show in Appendix. The link between the surface constant C0 and the surface capacitance Csurf is explored by changing C0 and subsequently measuring the change in Csurf. The Csurf is extracted by extrapolating a thickness-capacitance curve to zero thickness, since the surface capacitance is expected to be independent of the thickness. The overall cell capacitances with varying film thicknesses are shown in Fig. 3, where several C0s are used to reproduce the experimental data (the hollow dots) and the correct fitting of the slope was done with ˜k0f = 5 1031 at the experimental oxygen pressure. When C0 = 0, the extrapolated line passes through the origin, indicating vanishing surface capacitance, which is clearly not applicable in this case. Increasing C0 shifts the line upwards, resulting in the increase of Csurf correspondingly. Their relationship is shown in the inset of Fig. 3, where a linear line is obtained. The shift of the lines by varying C0 also suggests that the surface capacitance might be additive to the bulk capacitance, since no change in the line slope can be detected. Thus, one may expect that the surface capacitor is parallel to the thickness-dependent chemical capacitor in the corresponding ECM. In addition, our model is consistent with physical intuition as well as experimental evidence since the chemical capacitance grows linearly with the thickness.26 The experimental values are fitted with C0 = 0.005 F cm2 and the corresponding surface capacitance is determined to be Csurf E 0.00154 F cm2. Two more experimental data sets with different oxygen partial pressures are compared to computational EIS using the same C0 = 0.005 F cm2 and are presented in Fig. 4. The

Paper

Fig. 4 Capacitances versus the thin film thickness with various oxygen partial pressures. The solid line is from numerical calculation, and the markers are the experimental results.16 The dashed lines are extrapolated based on the numerical values.

three extrapolated lines appear to intersect at the same point on the ordinate, suggesting that the surface capacitances do not change with oxygen partial pressures. However, the thicknessrelated chemical capacitances behave quite differently: the chemical capacitances grow faster with thickness at a smaller oxygen partial pressure. This is caused by the aforementioned fact that the migrating species concentration decreases with increasing oxygen pressure, which results in a higher chemical capacitance at lower oxygen pressures.16 4.2 Interplay between surface reaction and in-plane electronic drift-diffusion

Fig. 3 The overall capacitances with different film thicknesses at T = 25 600 1C and p ˜ O2 = 2 10 . The solid lines are numerical capacitances with different C0s, with the extended dashed parts denoting the extrapolated values, whose intersections with the ordinate represent the surface capacitances Csurf. The inset shows that the actual surface capacitance Csurf is linear with C0.


4.2.1 Vacancy and electron trajectories. In Fig. 1(c), we have illustrated the expected pathway of the migrating species within the bulk material. In the numerical experiments, we can further verify our intuitive understanding by computing the electron and vacancy path lines at zero frequency. The numerical results show that the electrons produced by reaction at the gas|ceria interface are collected by the metal collector, and none of them pass through the interface of ceria|YSZ due to the electronic blocking. Similarly, the flux of ions produced by the reaction converges at the ceria|YSZ interface, penetrating across the interface, as shown in Fig. 5. This is consistent with the physical understanding that the metal current collector is impenetrable to vacancies and the electrolyte is impenetrable to the electrons. The penetration depth of electrons into the film is also of interest, since it reflects the surface catalytic effect on the bulk material. Similar to the thick sample case, we discover that the electron penetration depth increases with thickness gradually until reaching a stable value, but the ionic pathways are barely affected. Under the mentioned conditions and the physical parameters, the maximum penetration depth is of the order


View Article Online


Paper

PCCP

important role in determining the shape of the Nyquist plot, shown in Fig. 6. In the simulations we kept W1 = W2, ensuring an equal proportion of the current collector and the exposed surface. If W1 = W2 is sufficiently small then impedance is in essence a semicircle as reported by the experimental results used in this work,16 suggesting the sole contribution from the surface reactions. This is reasonable since for small W1 the electron in-plane transport is small. In the middle range, from 10 mm to 1000 mm, the shape changes from a semicircle to an asymmetric depressed arc. This confirms that two processes play the role. Given the analysis above we can safely assume that one process corresponds to the surface reaction and the other corresponds to the in-plane transport. If we increase W1 = W2 further, the shape does not change with an additional increase in the spacing, indicating that this is a limiting phenomenon. The phenomenological shape changes are, of course, concurrent with the change in resistance values. Herein ? we define the total surface impedance as45 Zion ¼ Ð . . . ð1Þ Ð ð1Þ ð1Þ ^ îon dx W2 G4 m êon dx W1 j in-plane , where UT G5 m ˆj (1) is the in-plane current from the electron transport, in-plane

i.e., the current collected by the metal. The total surface impedance consists of both the contribution from the surface Ð Ð ð1Þ . ð1Þ ð1Þ îon dx G5 m êon dx reaction Zsurf ¼ UT G5 m W2 jîn-plane

Fig. 5 Numerical flow lines of electrons (left) and vacancies (right) of different film thickness at excitation frequency f = 0.

of 2–3 mm. This value is consistent with previous work.52,53 We stress that it is not obvious that those two values are similar due to the distinctly different boundary conditions in the thin film configuration and in the thick film sample. 4.2.2 Influence of the current collector distance on the EIS. The current collector spacings, namely W1 and W2, play an

and that from the bulk in-plane electronic transport Zbulk. We further define the R> as the total surface resistance, which is ion the value of Z> at zero frequency, and the surface resistance ion Rsurf, the value of Zsurf at zero frequency. As an example, the change in surface resistance R> with increasing spacing is ion discussed, as shown in Fig. 7. By definition here, the surface resistance is the sum of both the reaction resistance and the inplane transport resistance. Its value, therefore, reflects the overall resistance within the cell when the bulk ionic transport is diminished and is of great importance in the thin film. In the low spacing region, e.g., below 10 mm, the R> ion is insensitive to the change in spacing, which proves the sole role of the surface reaction resistance. In the intermediate region, from 10 mm to

Fig. 6 Normalized impedances of a 1 mm-thick film with changing current collector spacing. The spacing changes from (a) W1 = W2 = 1 mm to 10 mm, (b) W1 = W2 = 10 mm to 1000 mm and (c) W1 = W2 = 1000 mm to 10 000 mm.



View Article Online


PCCP

Fig. 7 Total surface resistance upon changing the spacing.

1000 mm, the surface resistance increases from the low constant region to an exponentially increasing region. In the high spacing region, however, the total surface resistance increases 2 exponentially with the spacing, i.e., R> ion p W1 according to the fitting, and the slope in the logarithm scale was demonstrated using a triangle. In conjunction with Fig. 6, the change in the physical phenomenon and the rate limiting steps is well elucidated. The transition from reaction controlled resistance to the transport controlled resistance is explained previously. However, the quantitative description of the interplay has not been drawn yet. Thus, as a further step, the relative importance of the surface reaction resistance and the in-plane electronic transport resistance is studied, as shown in Fig. 7. The relative importance ratio is defined as fsurf = Rsurf/R> , ion where Rsurf is the contribution from the surface reaction. The change in the total surface resistance is achieved by altering the surface reaction rate constant and meanwhile keeping the collector and exposed area width identical, i.e., W1 = W2. The numerical results show that for the small spacing, e.g., W1 = W2 = 2.5 mm, the surface reaction dominates the overall resistance, rendering a fsurf value close to unity, as shown in Fig. 8. Instead, at larger spacing, the electronic in-plane diffusion starts to play the main role, as indicated in the line with W1 = W2 = 40 mm, where fsurf sharply decreases with decreasing R> . ion Furthermore one notices that if the fsurf factor is identical, lower spacing generally implies a lower overall resistance R> ion. These results are again consistent with previous work.45,52 In all, the current collector spacing is an important quantity that affects the impedance responses. It also suggests that the studies of electrodes using geometrically ill-defined porous current collectors may lack consistencies in different reports, since the feature of the current collector can hardly be controlled, which however would affect the final impedance results greatly.


Paper

Fig. 8 The ratio of the reaction resistance to the total surface resistance, fsurf, varies with the total surface resistance R> and with different collector ion pattern distances W1 and W2. The thickness t = 1 mm.

5 Conclusions A general 2D numerical work for impedance spectroscopy responses of a dense MIEC thin film in a symmetric cell configuration was developed. The model assumes charge neutrality in the bulk material and employs a double layer correction at the reaction surface to capture the boundary layer capacitive behaviour, an important feature that could dominate the capacitive responses in the case of thin films. Boundary chemical reactions as well as the bulk transport phenomenon are incorporated in our model to explain their interplay within the material. As an example, we performed FEM simulations that model single crystal doped ceria thin film electrodes in a symmetric cell configuration in a reducing atmosphere. We compared our simulation results with the experimental data that are available in literature. The validity of our model is demonstrated by showing that both the numerical resistance and capacitance have captured the experimental results well. The overall capacitance derived from the impedance arcs is fitted by taking into account the surface capacitance, which is shown to be additive to the bulk chemical capacitance. We also investigated various conditions that could affect the impedance response. First the oxygen vacancy and electron fluxes are computed and are shown to agree well with intuition. The electronic and ionic rails suggest the existence of a maximum electronic penetration depth for the thin film. We further find that increasing the current collector pattern spacing changes the impedance Nyquist plot from a semicircle to an asymmetric depressed shape, accompanied by a dramatic increase in the resistance values, indicating that multiple processes may be identified simultaneously, most notably the surface reaction and the electronic in-plane transport. Lastly, the film resistances with respect to reaction rate constants and


View Article Online

Paper

PCCP


the collector pattern distances are presented. A quantitative understanding of the transition from reaction rate controlled resistance at small pattern spacing to electronic in-plane diffusion controlled resistance at large pattern spacing highlights the importance of the current collector spacing, especially for the thin films, and the need for multidimensional modeling under those conditions.

List of symbols cm c(0) m c(1) m f f(1) B kB m~eon

m~ion

m~(1) eon

m~(1) ion

Concentration of species m, m = eon for electrons and m = ion for ions Equilibrium concentration of species m First order perturbation concentration of species m Electric potential First order perturbed electric potential Background charge concentration Boltzmann constant Electrochemical potential of electrons, defined as ceon ~eon ¼ m0eon þ kB T ln ð0Þ ef m ceon Electrochemical potential of ions, defined as cion ~ion ¼ m0ion þ kB T ln ð0Þ þ 2ef m cion First order perturbed electrochemical potential of ð1Þ ceon ð1Þ ~eon ¼ kB T ð0Þ efð1Þ electrons m ceon First order perturbed electrochemical potential of ions ð1Þ

ð1Þ

~ion ¼ kB T m (~ m(1) )* m

ð1Þ

~m m

ð1Þ

n(1) p(1) f~(1) ˜m p Kg

Kr

F ˆ(1) n ˆ(1) p f^(1) e lD lc

cion ð0Þ

ceon

þ 2efð1Þ

ð1Þ

¼

~m m zm kB T

ceon ð0Þ

ceon ð1Þ cion ð0Þ

cion efð1Þ kB T Pressure for gas m normalized by 1 atm Equilibrium constant for oxygen and hydrogen gas reaction, p~2H2 O Kg ¼ p~O2 p~H22 Equilibrium constant for oxygen defect reaction, ! ð0Þ 2 ð0Þ ceon cion 1=2 Kr ¼ p~ B B O2 Fourier transform operator F(n(1)) F( p(1)) F(f^(1)) rffiffiffiffiffiffiffiffiffiffiffi Dielectric constant ekB T Debye length, defined as e2 B Characteristic length


l

n%

Inverse of the normalized Debye length, lc defined as lD Normalized equilibrium electronic concentration, ð0Þ

ceon B Normalized equilibrium ionic concentration, defined ð0Þ c as ion B Mobility of species m um kB T Diffusivity of species m, computed by Dm ¼ jzej

defined as p%

um Dm tn

tp tn* ti* kB UT Z> ion Zsurf R> ion Rsurf

Characteristic time for electronic diffusion, defined as lc2 Deon lc2 Characteristic time for ionic diffusion, defined as Dion n n 1þ tn þ tp 4 p 4 p 4 p 1þ tp tn n Boltzmann constant Thermal potential, defined as kBT/e The total surface impedance, incorporating the surface reaction impedance and the bulk transport impedance The surface reaction impedance The value of Z> at zero frequency ion The value of Zsurf at zero frequency

Appendix The rationale behind boundary correction is illustrated in Fig. 9. At the reaction boundary, charge neutrality of perturbed system fails, and full equations should be considered. The full Poisson–Nernst–Planck (PNP) model at the boundary is the sum of contributions due to different models, one characterized by charge neutrality and one corresponding to surface/double layer effects. The full model is derived by first linearizing the Poisson–Nernst–Planck equations and by subsequently taking its Fourier transform, thus any solution to this model may be written as the superposition of two sets of solutions, namely an

Fig. 9 A schematic of ionic and electronic fluxes at the reaction boundary. The total flux of each species consists of the capacitive current and the reaction faradaic current.


View Article Online

PCCP

Paper

electrically neutral (EN) and a double layer (DL) solution. This may be written as follows ˆ(1) ˆ(1) ˆ(1) ˆ(1) + f^(1) ) = l2(%n(n +n ) 2%p(p +p )) Dx˜(f^(1) EN DL EN DL EN DL

(6a)

ˆ(1) ˆ(1) ˆ(1) ˆ(1) iotn(n +n ) + Dx˜(f^(1) + f^(1) ) Dx˜(n +n )=0 EN DL EN DL EN DL

(6b)


ˆ(1) iotp(p EN

ˆ(1) p ) DL

2Dx˜(f^(1) EN

f^(1) ) DL

ˆ(1) Dx˜(p EN

ˆ(1) p ) DL

+ + = 0 (6c) pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi where l ¼ lc e2 B=ekB T is the characteristic length lc divided pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi by the Debye length lD ¼ ekB T=e2 B and e is the dielectric constant. ˆ(1) ˆ(1) We can plug in the charge neutral solution (f^(1) EN, nEN, pEN) into the corresponding DL solution and this will yield +

ˆ(1) ˆ(1) + iotn*n = l2(nn 2%pp ) Dx˜f^(1) % ˆ(1) DL EN DL DL ˆ(1) iotnn DL

+

Dx˜f^(1) DL

ˆ(1) Dxñ DL

=0

ˆ(1) ˆ(1) iotpp 2Dx˜f^(1) Dx˜p =0 DL DL DL

(7a)

We can compute the charge excess for unit area ð1 lD eB ð1Þ npffiffiffiffiffiffiffiffiffiffiffiffiffiffiADL ¼ e nB^ nDL dy ¼ q^ð1Þ n n þ 4 p 0 ð1 lD eB ð1Þ ppffiffiffiffiffiffiffiffiffiffiffiffiffiffiADL pB^ pDL dy ¼ 4 q^ð1Þ p ¼ 2e n þ 4 p 0

io

(7c)

otn ð1Þ 1 @2 @ 2 ^ ð1Þ 1 @2 @2 ð1Þ ^ f n^DL ¼ 0 þ þ þ n DL DL l2 l2 @xþ2 @yþ2 l2 @xþ2 @yþ2 (8b)

otp ð1Þ 1 @2 @ 2 ^ ð1Þ 1 @2 @2 ð1Þ fDL 2 p^DL ¼ 0 þ þ i 2 p^DL 2 2 l l @xþ2 @yþ2 l @xþ2 @yþ2 (8c)

@ 2 ^ ð1Þ ð1Þ ð1Þ f ¼ nn^DL 2 pp^DL @yþ2 DL

(9a)

@ 2 ^ ð1Þ @ 2 ð1Þ n^ ¼ 0 fDL 2 @yþ @yþ2 DL

(9b)

2

(13a)

(13b)

q^ð1Þ n ^ Seon jêon jDL!BULK ¼ 0 þo e

(14a)

q^ð1Þ p ^ Sion jîon jDL!BULK ¼ 0 þo 2e

(14b)

io

^ Sion are the Fourier transformed electron and ion ^ Seon o where o production rates, respectively, due to the surface reaction. ! ð0Þ Dion ~0 1=4 ceon S ^ eon ¼ 4 k p~ 1 þ ð0Þ p~H2 nð1Þ (15) o lc f O2 4cion . . ð0Þ ð0Þ If we define k ¼ 2Dion k~0f p~bO2 1 þ ceon 4cion lc and 1=4 k~f ¼ k~0 p~ , the mass balance equations at the boundary become: f O2

eDeon cð0Þ eon

ð0Þ

It is reasonable to assume that for sufficiently low frequency, we have otn/l2, otp/l2 { 1 and consequently, otn*/l2 { 1. Hence O(1/l2) terms can be neglected and the PNP equations become:

(12)

where ADL is an integration constant. If we integrate the Nernst– Planck equation in the double layer for both electrons and ions and Fourier transform them, we can obtain:

(7b)

One notes that if l is large, then we can solve this problem by stretching the corresponding coordinate and transforming the coordinate from x˜ to x+, where (x+,y+) = (x˜,y˜l). The PNP model in the new coordinates reads: 1 @2 @ 2 ^ ð1Þ otn ð1Þ ð1Þ ð1Þ fDL þ i 2 nÊN ¼ nn^DL 2 þ pp^DL (8a) 2 @x 2 2 @yþ l l þ i

which leads to the solution for f^(1) DL pffiffiffiffiffiffiffiffiffiffiffiffiffiffi y ^ ð1Þ ¼ ADL exp 4 p þ n f DL lD

2eDion cion

^ ð1Þ @ n^ð1Þ @ f @ y~ @ y~ ^ ð1Þ @ p^ð1Þ @f þ2 @ y~ @ y~

! ¼ 2elc k~ pH2 n^ð1Þ þ io^ qð1Þ n lc

!

(16a)

p ¼ 2elc k~ pH2 n^ð1Þ þ io4 q^ð1Þ lc (16b) n n

These two flux conditions define the boundary condition at the reaction interface. By further simplifying the equations we get the corrected boundary condition at the reaction surface: @ n^ð1Þ ¼ An n^ð1Þ þ Bn @ y~

(17a)

ð1Þ

^ @f ¼ Af n^ð1Þ þ Bf @ y~

(17b)

2

@ ^ ð1Þ @ ð1Þ f þ 2 p^ ¼ 0 @yþ2 DL @yþ2 DL

(9c)

From (9b) and (9c) and their corresponding boundary conˆ(1) ˆ(1) ditions, i.e., |f^(1) n | | o N, f^(1) n | = 0, DL DL y+-N DL DL BULK (1) (1) (1) (1) ˆDL|y+-N| o N and 2f^DL + p ˆDL|BULK = 0 we can derive that |2f^DL + p ˆ(1) =n f^(1) DL DL

(10a)

ˆ(1) = 2f^(1) p DL DL

(10b)

Thus, the following boundary value problem is obtained. @ 2 ^ ð1Þ ^ ð1Þ p Þf f ¼ ðn þ 4 DL @yþ2 DL


(11)

where An ¼ k~f

p~H2 ð0Þ

cion

ð0Þ

1þ4

Dion cion

!

ð0Þ

Deon ceon

1 1 Deon Dion lD lc pffiffiffiffiffiffiffiffiffiffiffiffiffiffiADL Bn ¼ 4io ð0Þ n þ 4 p ceon 4 þ ð0Þ cion Af ¼ k~f

p~H2 ð0Þ

cion

1

Dion Deon

(18a)

(18b)

(18c)


View Article Online

Paper

PCCP ð0Þ

Bf ¼ io

4 ceon 1 þ Dion cð0Þ Deon ion

ð0Þ

ceon ð0Þ

cion

þ4

lD lc pffiffiffiffiffiffiffiffiffiffiffiffiffi ffiADL n þ 4 p

(18d)


We define the surface constant: ð0Þ

4lD ADL c e2 ffi C0 ¼ pffiffiffiffiffiffiffiffiffiffiffiffiffiion n þ 4 p kB T

(19)

Bn and Bf can be rewritten as: 1 1 D Dion kB T C0 Bn ¼ iolc eon 2 ð0Þ ð0Þ 4cion þ ceon e

(20a)

ð0Þ

Bf ¼ iolc

1 ceon 1 þ Dion cð0Þ 4Deon kB T ion

ð0Þ

ð0Þ

ceon þ 4cion

e2

C0

(20b)

Acknowledgements The authors gratefully acknowledge HKUST for providing startup funds, support from the University Grants Committee of Hong Kong through the projects DAG12EG07-12, DAG12EG06, and ECS 639713. The authors thank Prof. S. M. Haile (California Institute of Technology) for fruitful discussions and significant suggestions.

Notes and references 1 S. B. Adler, Chem. Rev., 2004, 104, 4791–4843. 2 W. C. Chueh, Y. Hao, W. Jung and S. M. Haile, Nat. Mater., 2012, 11, 155–161. 3 Z. P. Shao and S. M. Haile, Nature, 2004, 431, 170–173. 4 F. S. Baumann, J. Maier and J. Fleig, Solid State Ionics, 2008, 179, 1198–1204. 5 J. Fleig, F. S. Baumann, V. Brichzin, H. R. Kim, J. Jamnik, G. Cristiani, H. U. Habermeier and J. Maier, Fuel Cells, 2006, 6, 284–292. 6 W. Jung and H. L. Tuller, Adv. Energy Mater., 2011, 1, 1184–1191. 7 S. J. Litzelman, J. L. Hertz, W. Jung and H. L. Tuller, Fuel Cells, 2008, 8, 294–302. 8 Y. L. Yang, C. L. Chen, S. Y. Chen, C. W. Chu and A. J. Jacobson, J. Electrochem. Soc., 2000, 147, 4001–4007. 9 J. R. Macdonald, Ann. Biomed. Eng., 1992, 20, 289–305. 10 E. Koep, D. S. Mebane, R. Das, C. Compson and M. Liu, Electrochem. Solid-State Lett., 2005, 8, A592–A595. 11 G. T. Kim, S. Wang, A. J. Jacobson, Z. Yuan and C. Chen, J. Mater. Chem., 2007, 17, 1316–1320. 12 F. Baumann, J. Maier and J. Fleig, Solid State Ionics, 2008, 179, 1198–1204. 13 Y. Yang, C. Chen, S. Chen, C. Chu and A. Jacobson, J. Electrochem. Soc., 2000, 147, 4001–4007.


14 W. Jung and H. L. Tuller, Solid State Ionics, 2009, 180, 843–847. 15 H.-S. Noh, J.-W. Son, H. Lee, H.-S. Song, H.-W. Lee and J.-H. Lee, J. Electrochem. Soc., 2009, 156, B1484–B1490. 16 W. C. Chueh and S. M. Haile, Phys. Chem. Chem. Phys., 2009, 11, 8144–8148. 17 J. Nielsen and J. Hjelm, Electrochim. Acta, 2014, 115, 31–45. 18 D. Chen, R. Ran and Z. Shao, J. Power Sources, 2010, 195, 7187–7195. 19 D. Chen, R. Ran, K. Zhang, J. Wang and Z. Shao, J. Power Sources, 2009, 188, 96–105. 20 V. Dusastre and J. A. Kilner, Solid State Ionics, 1999, 126, 163–174. 21 J. Fleig, J. Electroceram., 2004, 13, 637–644. 22 S. McIntosh, J. Vohs and R. Gorte, J. Electrochem. Soc., 2003, 150, A1305–A1312. 23 K. Katsiev, B. Yildiz, K. Balasubramaniam and P. A. Salvador, Appl. Phys. Lett., 2009, 95, 092106. 24 M. W. Louie, A. Hightower and S. M. Haile, ACS Nano, 2010, 4, 2811–2821. 25 J. Jamnik and J. Maier, J. Electrochem. Soc., 1999, 146, 4183–4188. 26 J. Jamnik and J. Maier, Phys. Chem. Chem. Phys., 2001, 3, 1668–1678. 27 W. Lai and S. M. Haile, J. Am. Ceram. Soc., 2005, 88, 2979–2997. 28 M. E. Lynch and M. L. Liu, J. Power Sources, 2010, 195, 5155–5166. 29 D. D. Macdonald, Electrochim. Acta, 2006, 51, 1376–1388. 30 W. Lai and S. M. Haile, Phys. Chem. Chem. Phys., 2008, 10, 865–883. 31 M. L. Liu and J. Winnick, Solid State Ionics, 1999, 118, 11–21. 32 H. C. Chang and G. Jaffe, J. Chem. Phys., 1952, 20, 1071. 33 S. Hershkovitz, S. Baltianski and Y. Tsur, Solid State Ionics, 2011, 188, 104–109. 34 Q. A. Huang, R. Hui, B. W. Wang and H. J. Zhang, Electrochim. Acta, 2007, 52, 8144–8164. 35 M. E. Orazem and B. Tribollet, Electrochemical impedance spectroscopy, Wiley. com., 2011. 36 A. N. Morozovska, E. A. Eliseev, N. Balke and S. V. Kalinin, J. Appl. Phys., 2010, 108, 053712. 37 A. N. Morozovska, E. A. Eliseev, S. L. Bravina, F. Ciucci, G. S. Svechnikov, L. Q. Chen and S. V. Kalinin, J. Appl. Phys., 2012, 111, 014107. 38 S. Jesse, A. Kumar, T. M. Arruda, Y. Kim, S. V. Kalinin and F. Ciucci, MRS Bull., 2012, 37, 651–658. 39 Y.-M. Kim, J. He, M. D. Biegalski, H. Ambaye, V. Lauter, H. M. Christen, S. T. Pantelides, S. J. Pennycook, S. V. Kalinin and A. Y. Borisevich, Nat. Mater., 2012, 11, 888–894. 40 S. B. Adler, Solid State Ionics, 2000, 135, 603–612. 41 F. Ciucci and D. G. Goodwin, ECS Trans., 2007, 7, 2075–2082. é, Electro42 A. Kromp, H. Geisler, A. Weber and E. Ivers-Tiffe chim. Acta, 2013, 418–424. é, 43 H. Zhu, A. Kromp, A. Leonide, E. Ivers-Tiffe O. Deutschmann and R. J. Kee, J. Electrochem. Soc., 2012, 159, F255–F266.


View Article Online


PCCP

¨ffelin, J. Joos, M. Ender, A. Weber and E. Ivers-Tiffe é, 44 A. Ha J. Electrochem. Soc., 2013, 160, F867–F876. 45 F. Ciucci, Y. Hao and D. G. Goodwin, Phys. Chem. Chem. Phys., 2009, 11, 11243–11257. 46 M. E. Lynch, D. S. Mebane, Y. Liu and M. Liu, J. Electrochem. Soc., 2008, 155, B635–B643. 47 D. Chen, S. R. Bishop and H. L. Tuller, J. Electroceram., 2012, 28, 62–69. 48 J. Fleig, Phys. Chem. Chem. Phys., 2005, 7, 2027–2037.


Paper

49 S. Adler and W. Bessler, Handbook of Fuel Cells–Fundamentals, Technology and Applications, 2009, pp. 441–462. 50 F. Hecht, J. Numer. Math., 2012, 20, 251–266. 51 W. C. Chueh and S. M. Haile, Philos. Trans. R. Soc., A, 2010, 368, 3269–3294. 52 F. Ciucci, W. C. Chueh, D. G. Goodwin and S. M. Haile, Phys. Chem. Chem. Phys., 2011, 13, 2121–2135. 53 A. Tarancon, M. Burriel, J. Santiso, S. J. Skinner and J. A. Kilner, J. Mater. Chem., 2010, 20, 3799–3813.


[Thin-film electrodes for medical diagnosis].

Modeling thin-film piezoelectric polymer ultrasonic sensors.

Evaluation of back contact in spray deposited SnS thin film solar cells by impedance analysis.

Engineering interface-type resistive switching in BiFeO3 thin film switches by Ti implantation of bottom electrodes.

Titanyl phthalocyanine ambipolar thin film transistors making use of carbon nanotube electrodes.

Absorptive carbon nanotube electrodes: consequences of optical interference loss in thin film solar cells.

Rapid prototyping of all-solution-processed multi-lengthscale electrodes using polymer-induced thin film wrinkling.

Thin-film electrodes for high-capacity lithium-ion batteries: influence of phase transformations on stress.

Aging-induced chemical and morphological modifications of thin film iron oxide electrodes for lithium-ion batteries.

Fabrication and Characterization of Thin Film Nickel Hydroxide Electrodes for Micropower Applications.

Clean graphene electrodes on organic thin-film devices via orthogonal fluorinated chemistry.

Transparent and flexible supercapacitors with single walled carbon nanotube thin film electrodes.

Carbon nanotube and graphene hybrid thin film for transparent electrodes and field effect transistors.

Graphene-based electrodes for enhanced organic thin film transistors based on pentacene.

High performance semiconducting enriched carbon nanotube thin film transistors using metallic carbon nanotubes as electrodes.

Probing the morphological influence on solid electrolyte interphase and impedance response in intercalation electrodes.

Thin film characterization for modeling and optimization of silver-dielectric color filters.

Applications of thin-film sandwich crystallization platforms.

Interlaboratory Comparison of Magnetic Thin Film Measurements.

Multiscale modeling of the electrostatic impact of self-assembled monolayers used as gate dielectric treatment in organic thin-film transistors.

Flexible logic circuits based on top-gate thin film transistors with printed semiconductor carbon nanotubes and top electrodes.

SnO2 electrodes for flexible and transparent thin film heaters.

Materials analyses and electrochemical impedance of implantable metal electrodes.

Comparison of impedance cardiographic measurements using band and spot electrodes.