Extending Newman’s Pseudo-Two-Dimensional Lithium-Ion Battery Impedance Simulation Approach to Include the Nonlinear Harmonic Response

The pseudo-two-dimensional (P2D) model of lithium-ion batteries couples a volume-averaged treatment of transport, reaction, and thermodynamics to solid-state lithium diffusion in electrode particles. Here we harness the linear and nonlinear physics of the P2D model to evaluate the fundamental (linear) and higher harmonic (nonlinear) response of a LiCoO 2 | LiC 6 cell subject to moderate-amplitude sinusoidal current modulations. An analytic-numeric approach allows the evaluation of the linearized frequency dispersion function that represents electrochemical impedance spectroscopy (EIS) and the higher harmonic dispersion functions we call nonlinear electrochemical impedance spectroscopy (NLEIS). Base case simulations show, for the ﬁrst time, the full spectrum second and third harmonic NLEIS response. The effect of kinetic, mass-transport, and thermodynamic parameters are explored. The nonlinear interactions that drive the harmonic response break some of the degeneracy found in linearized models. We show that the second harmonic is sensitive to the symmetry of the charge transfer reactions in the electrodes, whereas EIS is not. At low frequencies, NLEIS probes aspects of the cell thermodynamics that are not accessible with EIS. In short, NLEIS has the potential to increase the number of physicochemical parameters that can be assessed in experiments similar in complexity to standard EIS measurements. ©TheAuthor(s) 2017. Published by ECS. This is an open access article distributed under the terms of the Creative Commons Attribution 4.0 License (CC BY, http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse of the work in any medium,

Electrochemical impedance spectroscopy (EIS) is a commonly utilized tool for the noninvasive analysis of a wide range of electrochemical systems. 1,2 Using a sufficiently small perturbation amplitude simplifies the interpretation of the impedance response because the system under analysis can be approximated as linear. The condition of linearity and ease of fitting equivalent circuit analogs has led to the widespread adoption of EIS for both characterization [3][4][5] and prognostic [6][7][8] applications within lithium-ion battery research. In particular, EIS has been used to garner insight into important physicochemical processes including the growth of the solid electrolyte interphase (SEI) layer, 9 mass transfer and kinetics in electrode materials, 10 as well as other degradative and capacity loss mechanisms. 11 The challenge with EIS, however, is that the analysis is restricted to probing the linear regime, despite the inherently nonlinear nature of the electrochemical system. Fundamentally, linearization reduces the information content of a given EIS measurement. The consequence of linearization is model degeneracy, in that, one EIS dataset can be represented equally well by different physical models 12 and different circuit analogs, 13 because of the loss of information inherent with linearization. This degeneracy and loss of information may be an especially important limitation for the study of lithium-ion batteries if we hope to carefully dissect degradative interfacial processes. 14,15 The issue of linear EIS model degeneracy and information loss is particularly challenging for sealed batteries. Analytical solutions for the pseudo-two-dimensional (P2D) model impedance response under stationary conditions show that linearized EIS is unable to uniquely determine many of the model parameters needed for simulating chargedischarge curves. 16,17 For example, the anodic and cathodic transfer coefficients, α a and α c , respectively, cannot be uniquely determined; only the sum (α a + α c ) appears in the linearized solution. Changes to the battery potential (state of charge) do not alter the inaccessibility of unique transfer coefficients. Alternative ways to measure unique transfer coefficients using EIS, such as passing a mean current that drives the electrode kinetics into the Tafel regime, introduce the chal-lenge of maintaining a stationary state during EIS measurements in a sealed cell.
An approach that breaks some of the degeneracy of linearized EIS without violating the requirement for stationarity is the use of moderate-amplitude current modulations around a zero-mean current. Moderate-amplitude modulations are used to drive a weakly nonlinear response, as evidenced by higher harmonics in the frequency response of the system. Typically, higher harmonics in impedance experiments have been treated as noise, to be eliminated as a way of ensuring the linearity of EIS measurements. [18][19][20][21] However, measuring the nonlinear harmonic response of an electrochemical system has been used to study corrosion, [22][23][24] fuel cells, 12,25,26 and mass transport controlled hydrodynamic systems [27][28][29] as well as in chemical sensing applications. [30][31][32] The increased informational power of nonlinear EIS (NLEIS) measurements enabled Wilson et al. to distinguish between physical models which were degenerate in the linear regime, 12 while Nakata et al. used the sensitivity of the nonlinear response to distinguish between different chemical targets for improved sensor selectivity. 31 Extending EIS to probe the full-spectrum nonlinear response requires a more sophisticated approach than fitting lumped-parameter equivalent circuit analogs. Physics-based models of the full electrochemical system offer the flexibility and interpretability required for a meaningful analysis. Here we extend the P2D impedance model developed by Newman et al. 33,34 to include the linear and nonlinearlygenerated higher order harmonic response generated by a lithium-ion battery system. As with EIS measurements, NLEIS is a "local" probe of the battery around a single state of charge, so EIS and NLEIS experiments will typically involve a sequence of measurements at various open circuit potentials (states of charge) to parameterize the full battery behavior. Later papers show that the extra information available in an NLEIS experiment comes at a modest increase in experimental difficulty. To the best of our knowledge, full-spectrum nonlinear EIS has not been applied to the study of lithium-ion batteries. Yet, the ability to distinguish between different physical models, noninvasively probe nonlinear internal states, and investigate the symmetry and reversibility of degradative reactions make NLEIS a candidate for providing information-rich insights that can move the field forward.

The Mathematical Model
Pseudo-two-dimensional battery impedance model.-The widely adopted P2D battery model has been a powerful tool for understanding, 35 optimizing, 36,37 and controlling 38 lithium-ion batteries for many decades. The P2D model is a physics-based representation of the coupled nonlinear mass transport, thermodynamic, and reaction processes inside a lithium-ion battery. Porous electrode theory treatment of the electrodes requires volume averaging over the particles and electrolyte pores leading to two macrohomogeneous solid-and solution-phase continua superimposed in the x direction across the cell sandwich thickness. The x direction equations are coupled to solid-phase diffusion of lithium inside the electrode particles, governed by, where r is the radial distance from the center of the particle. Here we take the solid-phase diffusivity, D s , to be independent of concentration in the solid, c s . The concentration of lithium in the electrode particles is coupled to the interfacial flux from faradaic charge transfer by the boundary condition at the solid-solution interface, r = R p , where we assume Butler-Volmer kinetics, The total interfacial flux, j, also has a non-faradaic Helmholtz double-layer capacitive contribution, j dl , [4] where Here we follow the convention of Doyle et al. 33 for a simple model of interfacial processes where the double-layer capacitance, C dl , is dominated by cation rearrangement, and both C dl and the exchangecurrent density, i 0 , are independent of concentration. The current density in the electrolyte, i 2 , is related to interfacial processes by Faraday's law, where a is the specific interfacial area of the electrode. Binary electrolyte theory describes the concentration in the electrolyte, c, and solution-phase potential, φ, where is the void fraction, D ef f is the effective diffusivity, t 0 + is the transference number of Li + in the electrolyte, and κ ef f is the effective solution-phase conductivity. The potential in the solid, φ s , is governed by Ohms' law, [9] where σ ef f is the effective conductivity in the solid-phase and I is the total current density, −I = i 1 + i 2 .
Zero flux boundary conditions are taken at the positive and negative current-collectors, Current only flows in the solution-and solid-phases at the electrode-current collector and separator-electrode interfaces, respectively, i 2 = 0 and − σ ef f ∂φ s ∂ x = −I at x = 0 and x = L , [11][12] i 2 = −I and − σ ef f ∂φ s ∂ x = 0 at x = l neg and x = l neg + l sep [13][14] and the solution-phase potential is set to zero (φ = 0) at the negative current collector as a reference. The equations in the separator are given by Equations 6-8 with the interfacial flux set to zero ( j = 0). Continuity of the fluxes across the electrode/separator interfaces is maintained. In all cases, effective transport coefficients (diffusivities, conductivities) are modifications of the intrinsic transport coefficients using a Bruggeman-type tortuosity factor to account for the porous character of the medium. For example, D ef f = D Brugg for effective diffusivity, where D is the intrinsic diffusivity, and Brugg is the Bruggeman coefficient.
Frequency domain formulation.-The governing equations described above are both coupled and nonlinear. If we assume the system of equations is stable to finite amplitude current perturbations, then nonlinear theory suggests a pure single frequency sinusoidal input perturbation will generate an output response with spectral components at the excitation frequency (fundamental or 1 st harmonic) and integer higher harmonics (2 nd , 3 rd , and 4 th harmonic of the fundamental frequency). We have shown previously that a fastcomputing analytic-numeric method can be used for analyzing the linear and nonlinear response of many different electrochemical and transport-reaction systems when they are stable to finite amplitude perturbations. 12,28,29,39,40 Here we apply these methods to the analysis of the nonlinear impedance in a lithium-ion battery system using the above P2D equations. The key step in the fast analytic-numeric formulation is to expand all dependent variables in the governing equations and boundary condition as a complex Fourier series X (x, t; I, ω 1 ) = Re ∞ k=0X k (x; I, ω 1 ) exp ( jkω 1 t) [15] where X denotes any of the dependent variables (c s , c, i 2 , φ s , and φ), X k are the complex Fourier coefficients,X k =X k − jX k , for the respective dependent variables, x is the vector of spatial independent variables (x, r), t is time, I is the current modulation amplitude, ω 1 is the fundamental frequency of the perturbation, and k is an integer that indicates the harmonic number. The notation of X (a; b, c, . . .) is used to show X is the dependent variable, a is an independent variable, and b, c, . . . are parameters. For the weakly nonlinear, time-periodic (stationary) system studied here, the complex Fourier coefficients can be expressed by a Volterra series expansion in powers of the current modulation amplitude, The resulting double-indexed complex coefficientsX k, p (x; ω 1 ) within the summation of Equation 16 are the key frequency dispersion coefficients for the Fourier coefficients at harmonic number, k, and order of nonlinear coupling, p. When k = p = 1, the governing equations and boundary conditions represent the linear response of the system to infinitesimal current modulations (e.g., this is the condition Doyle et al. computed 33 to determine the traditional linear EIS of a half-cell). When k = p = 2, the governing equations and boundary conditions represent the leading order 2 nd harmonic response of the system of equations to a finite perturbation, and so on. Solutions for p > k, represent corrections to the leading order k = p solutions. Here, we ignore all correction terms and focus on the leading order (k = p) linear and nonlinear harmonics.
Computing the leading order linear and nonlinear spectra.-The resulting governing equations are shown in the Appendix. Table AI contains the governing equations and boundary conditions for simulating the linear impedance (k = p = 1) response. Tables AII and AIII contain the sets of equations for computing the second (k = p = 2) and third (k = p = 3) harmonic responses. The higher harmonic equations contain terms which depend on the lower harmonic solutions in a hierarchical manner such that the harmonics are solved sequentially (not simultaneously). The full hierarchical set of ordinary differential equations for the first three leading order impedance harmonics was solved using the generic Coefficient Form PDE physics module in COMSOL v4.4. We chose to mesh the solid-phase equations (A1.1, A2.1, and A3.1) instead of using an analytical solution in order to have a general approach for any arbitrary concentration dependent diffusivity parameters. The mesh nodes were spaced using a geometric ratio such that there existed a high concentration of nodes at the interfaces between the solid-and solution-phases (at r = R p,i ) and at the electrode/separator interfaces (at x = l neg and x = l neg + l sep ). Mesh refinement in the interfacial regions with high gradients was used so that the solution converged for all parameter and frequency combinations resulting in approximately 9,000 total elements. Further mesh refinement produced negligible changes in the computed results. The frequency parameter, ω 1 , was varied logarithmically from 10 −3 to 10 5 Hz with 10 points per decade. Each full spectrum (with 81 frequencies) required a computational time of around 25 minutes on a Dell Precision T1500 with an Intel Core i7 CPU @ 2.80 GHz and 8 GB RAM using a Windows 7 Professional 64-bit operating system. Parameters were updated and the results were saved to disk using the LiveLink for MATLAB version R2013b.
The leading order whole-cell impedance for the first harmonic (k = p = 1), the second harmonic (k = p = 2), and the third harmonic (k = p = 3) can be directly computed from the simulations asZ whereφ s k, p (x = L; ω 1 ) is the computed Fourier coefficient for the solid-state potential evaluated at the boundary of the positive current collector (x = L) andφ s k, p (x = 0; ω 1 ) is the computed potential at the negative current collector (x = 0). These impedance coefficients are complex with both a real and an imaginary component such thatZ k, p (ω 1 ) =Z k, p (ω 1 ) + jZ k, p (ω 1 ). The double indexed coeffi-cientZ 1,1 (ω 1 ) is identically equal to the normal (linear) electrochemical impedance of the system.Z 2,2 (ω 1 ) is the leading order nonlinear electrochemical impedance derived from the second harmonic of the voltage response to a pure cosine input current andZ 3,3 (ω 1 ) is the leading order nonlinear electrochemical impedance derived from the third harmonic response.

Results and Discussion
Simulations of a LiCoO 2 |LiC 6 cell.- Figure 1 shows Nyquist representations of simulated impedance spectra for a typical LiCoO 2 |LiC 6 cell using the parameters in Table I. The spectra cover a wide range of frequencies from 1 mHz to 100 kHz logarithmically spaced with 10 points per decade. The linear impedance response, Z 1,1 (ω 1 ), displays the expected characteristics of a lithium-ion battery with a purely real ohmic resistance at high frequencies, a depressed mid-frequency arc linked to kinetic processes, and a low-frequency mass transport and thermodynamic dominated tail. The high frequency resistance is given by the sum of the ohmic resistances in each region of the cell, [18] which evaluates to 0.467 m − m 2 for the parameters in Table I. The typical characteristic frequencies for the kinetics of each electrode are  Table I. The high frequency kinetics arc and low frequency mass transport tail are characteristic of a lithium-ion battery. The lack of a second harmonic response (b) in frequencies dominated by the reaction kinetics (1 Hz to 1 kHz) is related to the typical assumption of symmetric transfer coefficients, α a and α c .
given by, [19] which evaluate to 1284.5 H z and 1428.5 H z for the negative and positive electrodes, respectively. The two electrodes have a similar contribution to this kinetic arc for the base parameters (R ct,neg = 7.79 m − m 2 and R ct, pos = 7.0 m − m 2 ). The characteristic frequencies for the diffusional processes in the solution-and Anodic transfer coefficient   If the governing equations and boundary conditions representing the physics of a battery were linear, then the standard electrochemical impedance shown in Figure 1a would be the only response of the system to a sinusoidal current perturbation of any magnitude. However, because the P2D model is a coupled nonlinear set, it also produces a frequency-dependent higher harmonic response that we can compute. The five characteristic transport and reaction timescales, and the high frequency ohmic limit, described above, can be used to help interpret the higher order harmonic response. For instance, the second harmonic impedance,Z 2,2 (ω 1 ), shown in Figure 1b, only exhibits an appreciable signal in the characteristic frequency range for mass-transport and thermodynamics (ω 1 < 1 Hz). For the base case parameters used here, there is no appreciable second harmonic response at frequencies in the reaction kinetics regime (1 Hz − 10 kHz). As we will show in the next section, this lack of kinetically-driven second harmonic response is an artifact of the particular literature values used in the base case scenario (symmetric transfer coefficients on both electrodes); under more general conditions, the second harmonic is indeed sensitive to charge transfer.
The base case third harmonic spectrum, shown in Figure 1c, exhibits signal in the low-frequency mass-transport and thermodynamic regime as well as the mid-frequency regime dominated by the reaction kinetics, much like the standard linear electrochemical impedance spectrum shown in Figure 1a. For both the second and third harmonic, the nonlinear impedance asymptotes to zero in the limit of high frequencies. At high frequencies, the response of the battery is dominated by ohmic processes (that are governed by intrinsically linear physics). As a result, the nonlinear harmonics vanish at high frequencies whereas the linear impedance remains finite. Figure 2a shows that changing the positive electrode anodic transfer coefficients, α a, pos , in Equation 4, while keeping the sum α a + α c = 1, has no effect on the linear impedance,Z 1,1 (ω 1 ), response; all five curves lie directly on top of each other. The transfer coefficients, α a and α c , for a reaction are a measure of charge transfer symmetry on that electrode. The base case values of 0.5 for α a and α c used in Figure 1 for both electrodes means there is complete symmetry in oxidation and reduction on both electrodes; i.e. the magnitude of the kinetic overpotential is identical The linear response is insensitive to the transfer coefficients with all five curves superimposed upon one another. The second harmonic response is very sensitive to the asymmetry of the transfer coefficients with the presence of a kinetic arc only when α a = 0.5. While only the reaction symmetry of the positive electrode was varied here, a similar dependence on the negative electrode charge transfer symmetry also exists. The , , and symbols indicate 1 mHz, 1 Hz, and 1 kHz, respectively. during charge and discharge. There is no intrinsic reason to assume charge transfer symmetry on either electrode. A difference in α a and α c for either electrode would reflect asymmetry in the charge transfer processes on that electrode.

Effect of charge transfer symmetry.-
Mathematically, only the sum (α a + α c ) appears as a parameter in the linear k = p = 1 equations in the Appendix, and that sum is unity in all cases simulated, resulting in the linear impedance being insensitive to charge transfer symmetry and the overlapping spectra in Figure 2a. Figure 2b shows that symmetry in charge transfer is clearly important for the second harmonic response. The solid curve in Figure  2b is the same case as shown Figure 1b, with symmetric transfer coefficients in the battery (α a = α c = 0.5 for both electrodes). Mathematically, symmetric transfer coefficients nullify terms in the kinetic contributions to the second harmonic (k = p = 2) equations in the Appendix. When we break the symmetry of charge transfer in Figure 2b -making the positive electrode have more facile oxidation (α a > 0.5) or more facile reduction (α a < 0.5) -we see a second harmonic arc arise at frequencies that correspond to kinetic processes. Moreover, the direction of the kinetic arc depends on the direction of the charge transfer asymmetry. When the transfer coefficients promote more facile battery discharge (α a, pos < 0.5 and α a,neg > 0.5), the second harmonic kinetic arc moves toward the first quadrant in the Nyquist plot. In contrast, if the charge transfer asymmetry makes battery charging more facile (α a, pos > 0.5 and α a,neg < 0.5), the second harmonic response shifts to the third quadrant.
The mathematical origins of these second harmonic simulation results are given by the governing equations for k = p = 2. The real and imaginary oscillations in the second harmonic flux density at the solid/electrolyte interface are given by Equation A2.6r and A2.6i in the appendix, We find that terms involving oscillations in the overpotential at the fundamental frequency, η 1,1 , are significantly larger than the oscillating overpotential at the second harmonic, η 2,2 . As a result, the second term on the right-hand side (RHS) of both equations dominates the harmonic response in the second order flux density. The second term of the RHS disappears if α a = α c . This sensitivity to reaction asymmetry has been demonstrated experimentally by Xu and Riley 41 who studied the ferri-ferrocyanide redox couple as well as Heubner et al. 42 who used second harmonic measurements to determine the transfer coefficients in the charge transfer kinetics of LiFePO 4 .
The third harmonic shown in Figure 2c is also sensitive to charge transfer symmetry, but in an inverse sense from the second harmonic. The largest third harmonic kinetic loop is observed when there is symmetry, and then the kinetic loop shrinks as the symmetry is broken. It is also noteworthy that the shift in third harmonic is insensitive to whether the charge transfer asymmetry favors oxidation or reduction. Thus, Figure 2c produces overlapping curves for α a = 0.45 and α a = 0.55, and a different set of overlapping curves for α a = 0.35 and α a = 0.65.
The results of Figure 2 show that there is new information to be gleaned by measuring the second and third harmonics. In the case of the transfer coefficients, the standard linear EIS spectra (Figure 2a) are insensitive to the specific values, whereas the second harmonic provides information about the magnitude of the asymmetry and direction (charge or discharge) that is kinetically more facile. The third harmonic just gives insights into the magnitude of the asymmetry. Figure 3 shows the effect of the exchange current densities, i 0,i , and double-layer capacitances, C dl,i , for the positive and negative electrodes. An asymmetric set of transfer coefficients (α a,neg = 0.45 and α a, pos = 0.55) has been assumed to ensure there is a kinetic arc in the second harmonic. The remaining discussion of individual parameter effects on the linear and nonlinear harmonic response will focus on only the first and second harmonic, for clarity.

Effect of kinetic parameters.-
The linear impedance spectra in Figures 3a and 3c show the sensitivity of the kinetic arc width to the changes in exchange current  Table I except α a,neg = 0.45 and α a, pos = 0.55. The , , and symbols indicate 1 mHz, 1 Hz, and 1 kHz, respectively. density on each electrode. The increasing size of the kinetic arc is normally interpreted through the charge transfer resistances, R ct , with smaller exchange current densities resulting in larger R ct . Similar to the linear response, the second harmonic spectra in Figures 3b and 3d show a growing kinetic arc (at frequencies in excess of 1 Hz) as the exchange-current density is decreased. The width of the second harmonic kinetic arc is dependent on both the asymmetry of the charge transfer as well as on the value of the exchange current density.
Varying the double-layer capacitance, C dl , changes the characteristic frequency for each electrode without changing the overall charge transfer resistance of the electrode. Figures 3e-3h show that independently varying the capacitances by four orders of magnitude  Table I  separates the timescales for each electrode kinetic process and introduces distinct arcs for the two electrodes. In many ways, the linear and second harmonic show similar behaviors, namely, they go from being one distinctive arc to two clear arcs. Because the timescales between the kinetic, transport, and thermodynamic regimes are well separated for the base case parameters, changing either the exchange-current densities or the double-layer capacitances has little effect on the low frequency response of any of the harmonics for the values shown here.
Effect of mass transport.-In equivalent circuit approaches to impedance analysis, the low frequency tail of the linear EIS response is often used to extract information about the diffusional processes. In their original description of the P2D impedance model on which this work builds, Doyle et al. 33 noted that lithium-ion batteries have many different diffusional impedances, and for simple analysis, the low frequency spectrum must be dominated by a transport process in order for Warburg plots to be accurate. The overlapping diffusional time constants remain a challenge for extracting individual coefficients from impedance spectra. Here we explore if the nonlinear second harmonic spectra may offer more sensitivity. Figure 4 shows the effects of varying the solution-and solid-phase diffusion coefficients over several orders of magnitude. The characteristic "wagging" of the low frequency linear impedance tail shown in Figures 4a, 4c, 4e makes clear the difficulty in extracting a single diffusion coefficient from a Warburg analysis, since all three diffusional impedance responses manifest similarly. Figures 4b, 4d, 4f shows the effect of the mass-transport parameters on the second harmonic response. Similarly to the linear impedance response, the high frequency arc remains constant due to the well separated timescales between the kinetic and diffusion regimes. Figures 4b and 4d show that varying either the negative or positive solid-phase diffusion coefficients has the effect of changing the magnitude and phase of the low frequency response like the linear impedance. For particularly low values of the solution-phase diffusivity, however, Figure 4f shows that the structure of the low frequency second harmonic rearranges as this diffusional impedance dominates the response. In nearly all cases, the low frequency changes in second harmonic impedance are more dramatic than in the linear response, suggesting a full physics analysis of experimental linear and second harmonic data may result in a much greater ability to fit low frequency processes. We are exploring this currently.
Effect of thermodynamic parameters.-The thermodynamics of the cell determine the impedance response as the frequency goes to zero. The important thermodynamic parameters in the low frequency range are the derivatives of the open circuit potential, U (c s ), with respect to the concentration of lithium in the electrodes. Mathematically, the linear response shows a sensitivity to only the first derivative, dU dc s , while the second and higher harmonics depend on the second, d 2 U dc s 2 , and higher order derivatives, respectively. Each of these derivatives are evaluated at the lithium concentration in the solid present at the particular state of charge being evaluated, c s 0 . The sensitivity to second and higher order derivatives could be useful for applications such as estimation of state of charge or probing fundamental thermodynamic models predictions for lithium intercalation. The sign of the second derivative (the curvature of U (c s )) changes in between different phase-change plateaus in the open circuit potential. Figure 5 shows the dependence of these open circuit potential derivatives around the base case parameters on the linear and second harmonic spectra. As with the mass transfer case, one sees similar effects between the linear and second harmonic impedance spectra, though with exaggerated effects in the second harmonic, when the first derivative of the open circuit potential is varied (Figures 4a-4d). In contrast, the first harmonic is completely insensitive to the second derivative of open circuit potential, whereas the second harmonic is sensitive at low frequencies.

Conclusions and Future Implications
Physics-based impedance models are an important component of understanding and interpreting impedance spectra due to their flexibility and physically meaningful parameters. Here, we presented an extension of Newman's P2D impedance model for a lithium-ion battery cell sandwich and demonstrated the simulated linear and harmonic response for a LiCoO 2 |LiC 6 cell. The strong dependence of the second harmonic on the symmetry of the reaction kinetics for each electrode was demonstrated and the individual effects of many kinetic, mass transport, and thermodynamic parameters were shown.
Introducing higher harmonic spectra into the analysis of an electrochemical system provides additional informational content over the linear response alone. Expanding the model to introduce additional nonlinear contributions to the physics (concentration dependence of parameters like the diffusion coefficients, exchange-current densities, and double-layer capacitances, for example) could lead to the generation of additional insight into the P2D model and its parameters. Quantifying the sensitivity of the model to the many different physicochemical and geometric parameters could provide insight into the experimentally measurable information in the linear and harmonic spectra. Preliminary experimental data indicate that these higher harmonic signatures are easily measurable for commercial lithium-ion  Table I except α a,neg = 0.45 and α a, pos = 0.55. The , , and symbols indicate 1 mHz, 1 Hz, and 1 kHz, respectively. cells using standard potentiostats and frequency response analyzers. These results will be explored in the near future.