Two States of the Cathode Catalyst Layer Operation in a PEM Fuel Cell

We measure impedance of a standard Pt/C–based PEM fuel cell for a series of current densities from 50 to 400 mA cm − 2 . Using our recent model for extraction of spatially–resolved data from impedance spectra, we plot the dependence of the oxygen diffusion coefﬁcient D ox in the cathode catalyst layer (CCL) and D b in the gas–diffusion layer (GDL) on the distance along the cathode channel. While the GDL oxygen diffusivity is fairly uniform over the cell active area, the shape of D ox indicates that the cell is separated into two domains with high and low water contents in the CCL. We attribute this effect to the positive feedback loop between the rates of oxygen transport and liquid water evaporation in the CCL, leading to local CCL ﬂooding. of

Impedance spectroscopy becomes increasingly popular method of PEM fuel cell characterization. The database of scientific publications Scopus shows that a number of publications on theory and applications of electrochemical impedance spectroscopy (EIS) in fuel cells exhibits exponential growth ( Figure 1). The characteristic time of the exponent is 7 years, which means that the number of publications on this topic doubles every 4 years and 10 months. This manifestation of the Moore's law is a signature of extremely powerful method.
A great advantage of EIS is the possibility to extract information on kinetic and transport coefficients in the cell components from the cell spectra. This, however, requires quite sophisticated physics-based modeling of the cell impedance. [1][2][3][4][5][6][7][8][9][10][11][12] Being fitted to an experimental spectrum, models of that type return oxygen diffusion coefficients in the porous layers, the cathode catalyst layer (CCL) proton conductivity, kinetic rate constants of the oxygen reduction reaction (ORR) etc.
Recently, we have shown that the spectrum of a fuel cell contains information on spatial distribution of the transport parameters over the cell active area. 13 To explain this effect, suppose that the cell is separated into N virtual segments. The impedance of each segment "knows" what is going on in the segments located upstream due to the oxygen flow in the cathode channel. This flow transports small-amplitude perturbations of the oxygen concentration from one segment to another, and in this way the total cell impedance contains information on the local transport and kinetic coefficients in each individual segment.
In this work, we extract the spatial shapes of the oxygen diffusion coefficient in the CCL D ox from the experimental spectra of a standard PEM fuel cell. These shapes show that the cell active area is split into two distinct domains with high and low D ox . We attribute this effect to the instability arising due to positive feedback loop between the rates of oxygen transport and liquid water evaporation in the CCL.

Experimental
Experimental work has been performed using a segmented cell system and a test station developed at Hawaii Natural Energy Institute ( Figure 2). The segmented cell setup consists of the cell hardware, current transducer system and data acquisition device. The current transducer system was custom designed and employed closed loop Hall sensors (Honeywell CSNN191) for current measurements. The segmented cell system is run as a single cell using the test station and standard operating protocols. This configuration reduces any interferences from the segmented cell setup and provides operating conditions identical to real ones, since only the overall cell load is governed by the test station. The whole system enables simultaneous measurements of spatial current, voltage and impedance. Details of the segmented cell system can be found in Refs. 14,15.
The segmented cell hardware contains standard non-segmented and segmented flow field plates and a 100 cm 2 MEA. The segmented flow field consists of 10 segments, each of them has an area of 7.6 cm 2 , its own current collector and gas diffusion layer (GDL). The inlet segment is segment 1, while segment 10 is the cell outlet. The system design allows either anode or cathode to be segmented. In this work, the segmentation was applied to the cathode. Both flow field plates have the same ten parallel channel serpentine design. The fuel and oxidant streams are arranged in co-flow configuration.
We used commercially available 100 cm 2 catalyst coated membranes (CCM) provided by Gore. Pt loading was 0.4 mg Pt cm −2 for anode and cathode with the catalyst layer thickness of 10-12 μm. The thickness of the reinforced membrane varied in the range of 16-18 μm. Sigracet 25 BC was used as GDLs for both electrodes. 25 BC consists of carbon paper and a microporous layer (MPL) with total thickness of 220-235 μm and porosity of 80%. Segmented GDLs were used on the cathode side, whereas a single and uniform GDL was applied at the anode. Teflon gaskets (125 μm) were employed for both electrodes to ensure required compression ratio.
The cell was operated with H 2 /air gas configuration at a temperature of 80 • C. The anode/cathode conditions were 2/9.5 stoichiometry, 100/50% relative humidity and 150/150 kPa absolute pressure. The air stoichiometry of 9.5 has been taken to minimize nonuniformity of the oxygen concentration over the cell active area. The impedance  measurements were performed under galvanostatic control of the cell over the frequency range from 0.05 Hz to 10 kHz with 11 steps per decade. The amplitude of the sinusoidal current perturbation signal was chosen to get a cell voltage response of 10 mV or lower. The impedance spectra have been recorded simultaneously from 10 segments and from the whole cell. In this work, we used only the whole cell spectra for processing. The cell operating and geometrical parameters are listed in Table I. Note that in our previous work, 13 we used a low-Pt loaded MEA with the Pt loading of 0.1 mg Pt cm −2 ; here, we perform experiments with the standard MEA from Gore having Pt loading of 0.4 mg Pt cm −2 .

Model
Model equations and impedance.-The experimental segmented cell with the meander-like channel ( Figure 2) is modeled below as a segmented cell with the single straight channel ( Figure 3). The experimental spectra have been fitted using the physics-based model for PEMFC impedance; 13,15 the model is based on the following idea.
The impedance model is obtained from the system of throughplane conservation equations in the individual segment; the boundary condition for the oxygen concentration perturbation in each segment is derived from the oxygen mass transport equation in the channel. In this way, individual segments are linked by the oxygen transport in the channel. The through-plane model includes the proton charge conservation Equation 1, the Ohm's law 2, and the oxygen mass transport equations in the CCL 3, and in the GDL 4: The oxygen transport in the channel is described by the mass transport equation assuming plug flow: [5] where is the local oxygen flux on the GDL side of the channel/GDL interface.
Here, c, c b , c h are the oxygen concentrations in the CCL, GDL and channel, respectively, x and z are the coordinates through the cathode and along the channel, respectively, c in h is the reference (inlet) oxygen concentration, C dl is the double volumetric layer capacitance, D ox , D b are the oxygen diffusion coefficients in the CCL and the GDL, respectively, h is the channel depth, i * is the ORR volumetric exchange current density, b is the ORR Tafel slope, σ 0 is the CCL proton conductivity at the membrane interface, β is the inverse characteristic scale of the proton conductivity decay toward the GDL, l t , l b are the CCL and the GDL thickness, respectively, L is the channel length, and v is the air flow velocity. Note that the MEA used exhibits a strong decay of the CCL proton conductivity with the distance through the CCL depth; this effect is taken into account by the exponential factor in Eq. 2 with β = 7.3. 17 Following the procedure of impedance measurements, we apply a small-amplitude perturbation of potential to all the segments. Equations 3-5 are nondimensionalized, linearized, and Fouriertransformed to obtain a linear system of equations for the small perturbation amplitudes of the overpotentialη 1 and the oxygen concentrationsc 1 ,c 1 b andc 1 h in the ω-space. Here, the superscript 1 marks the small-amplitude perturbations, and the sign tilde marks the dimensionless variables. The cell impedance is calculated in a following way.
First, we split the cell into N virtual segments. The variablec 1 h appears in the boundary condition for the through-plane system of equations along thex-coordinate. Suppose that the shape ofc 1 h along the channel is known; this allows us to calculate the impedance of each individual segment using the 1D through-plane impedance model. 18 The procedure returns also the shape of the oxygen flux perturbatioñ N 1 b at the channel/GDL interface alongz, which stands as the sink term in the linearized version of Eq. 5. Solving the latter equation we get the updated profile ofc 1 h , and repeat the procedure of calculation of segment impedances. These iterations are performed until the desired accuracy of the total cell impedance is achieved. Here,Z seg,n is the dimensionless impedance of the nth segment.
The key feature of the model is that each segment is allowed to have its own set of the transport and kinetic parameters. The local throughplane impedance model is based on Eqs. 1-4 and in each segment it includes five parameters: the ORR Tafel slope b, the CCL oxygen diffusivity D ox and proton conductivity σ p , the double layer capacitance C dl , and the oxygen diffusion coefficient D b in the GDL. Thus, for the cell divided into 10 virtual segments, we get a least-squares problem with 50 parameters. Thanks to extremely powerful fitting algorithms available nowadays, this problem can be solved yielding the shapes of the aforementioned parameters along the channel.
Numerical details.-Spectra fitting has been performed using a self-written Python code invoking the least_squares procedure from the SciPy library. The boundary-value problems for the perturbation amplitudes and for the static shapes have been solved using the SciPy solve_bvp procedure. The Cauchy problem for the oxygen perturbation amplitude in the channel has been solved with the SciPy odeint solver.
The fitting has been performed in two steps. First, the spectrum has been fitted assuming that all the fitting parameters (b, C dl , σ p , D ox and D b ) are uniform along the cathode channel. In the second step, the problem has been solved with variable along the z-coordinate parameters, taking the values obtained in the first step as initial conditions for least_squares internal iterations. To reduce the time of calculations, a parallel version of the code has been written, which employs a multi-core structure of modern processors for PC. On a 2.4 GHz PC with the quad-core processor fitting took about 3 hours per spectrum.

Results and Discussion
The cell geometrical and operating parameters used in the calculations are listed in Table I Most interesting is that at the current densities of 100 to 400 mA cm −2 , the cell active area is split into two well-separated domains with high and low D ox (Figures 5b-5e). In the following discussion we will refer to the variant with 10 virtual segments. The last three segments 8-10 operate at D ox , which is several times less than D ox in the segments 1 to 7. At the cell current of 200 mA cm −2 , the ratio of the high to low D ox is nearly 10 (Figure 5c). The only plausible explanation of this effect is flooding: the CCL in the last three segments works under severe flooding conditions. Amazingly, the GDL in the segments 8 to 10 is not flooded and the GDL oxygen diffusivity does not change with the cell current density ( Figure 6). Note that the local cell current density is nearly uniform along the cathode channel; the decay of D ox in the last three segments does not lower the local current there (Figure 7). Furthermore, the distribution of the high-frequency resistance also does not show any significant decay in the last segments ( Figure 8). Since the largest contribution to the HFR gives the membrane resistance, we conclude that the membrane water content exhibits only a minor variation along the channel. Flooding in the segments 8 to 10, thus, affects only the catalyst layer.  To explain this effect we note that the key parameter describing operation of a partially flooded CCL is the reaction penetration depth l ox due to oxygen transport 19 where c 1 is the oxygen concentration at the CCL/GDL interface, and j 0 is the local cell current density. The regime of the CCL operation is determined by relation between l ox and the CCL thickness l t . If l ox l t , oxygen easily penetrates through the CCL and the shape of the ORR rate through the CCL depth is determined by proton transport. If the proton transport is fast, the ORR rate is uniform through the CCL depth, otherwise this rate peaks at the membrane surface.
In the opposite limit of l ox < l t , the ORR runs mainly at the CCL/GDL interface, as the oxygen does not penetrate deep into the CCL. This regime is sketched in Figure 9 for the case of l ox /l t 0.4. As can be seen, about 60% of the CCL thickness adjacent to the membrane is inactive: in this domain, the ORR rate is close to zero and hence the proton current density does not change (Figure 9). The effective CCL volume is, thus, reduced by the factor of 0.4.
On the other hand, liquid water produced in the ORR must be evaporated in-place to provide "clean" pathway for gaseous transport of oxygen. The total rate of water evaporation in the CCL is proportional to the active CCL volume (Eq. 5 in Refs. 20 and 21). Estimates of Eikerling 20 show that at a normal operating regime with l ox > l t , a standard 10-μm CCL is able to evaporate all liquid water produced in the ORR up to the cell current density on the order of 1 A cm −2 . However, reduction of the effective CCL depth lowers the CCL evaporation efficiency, which may lead to catalyst layer flooding. In terms of Figure 9. Schematic of the CCL operation at low reaction penetration depth due to poor oxygen diffusion, Eq. 8. Shown are the dimensionless shapes of the proton current densityj, the local ORR overpotentialη, the oxygen concentrationc and the ORR rateS. Fast proton transport is assumed, which translates into uniform profile ofη. The ratio l ox /l t 0.4 is assumed. Figure 9, the remaining in operation 40%-part of the CCL is not able to evaporate all liquid water produced in the ORR and it is flooded. (It is worth noting that the CCL with reduced catalyst loading and thereby lower depth is more prone to flooding.)

F825
The data in Table I and in Figures 5, 6 allow us to estimate l ox in the flooded domain. The oxygen concentration at the CCL/GDL interface is c 1 = c in h (1 − j 0 /j lim ), where j lim = 4F D b c in h /l b is the limiting current density due to oxygen transport in the GDL. With the data in Table I, and D b 0.02 cm 2 s −1 (Figure 6), for j 0 = 0.2 A cm −2 we get c 1 8.5 · 10 −6 mol cm −3 . With this c 1 and D ox 0.7 · 10 −4 cm 2 s −1 (Figure 5c, segments 8 -10), from Eq. 8 for the same current density we obtain l ox = 0.0012 cm. The CCL thus operates in the transition regime with l ox l t , which means quite a noticeable decay of the oxygen concentration c through the CCL depth: assuming exponential decay of c toward the membrane, c ∼ exp(−x/l ox ), (Eq. (22) in Ref. 19) c 0 at the membrane interface is nearly 3 times less than c 1 at the CCL/GDL interface. This decay shifts the ORR peak to the CCL/GDL interface and it lowers the effective volume for the liquid water evaporation.
This estimate does not take into account possible lowering of local CCL temperature close to the channel outlet. This may further contribute to lowering of CCL evaporation performance in the last segments, which appears to be insufficient already at j 0 = 100 mA cm −2 (Figure 5b). Note that in the non-flooded domain of the cell (segments 1 to 7), l ox is five to ten times larger and hence in this domain the relation l ox l t holds. The positive feedback loop between the CCL oxygen diffusivity D ox and the evaporation rate could, therefore, lead to instability of the CCL operation. Indeed, lowering of D ox reduces the oxygen penetration depth into the electrode, thereby reducing the effective volume of liquid water evaporation, which leads to higher liquid water saturation in the CCL. This, in turn, translates into lower D ox , which further reduces the water evaporation rate. As a result, the cell active area appears to be split into two domains with the low and high water contents in the catalyst layer. Moreover, the experiment shows that this two-domain configuration is stable with respect to variation of the cell current density ( Figure 5). It should be emphasized that by instability we mean the mechanism leading to the two separated in space stable domains of the CCL operation. We have no direct experimental evidences for transient functioning of this mechanism; however, the fact that the cell active area is separated into two domains with the non-flooded and flooded CCL is a clear consequence of this transient.
The instability discussed is a critical phenomenon: it does not happen when l ox l t , and it starts to develop as soon as the relation l ox l t holds. It is, thus, important to keep l ox two to three times larger than l t during the cell operation. Interestingly, the process of CCL flooding does not affect the GDL transport properties ( Figure 6).

Conclusions
We use our recent model to extract spatially resolved oxygen diffusion coefficients in the catalyst D ox and gas-diffusion D b layers from the impedance spectra of a PEM fuel cell measured at the cell current densities in the range of 50 to 400 mA cm −2 . The shape of D ox along the oxygen channel shows that the cell active area is split into two domains with high and low D ox . We attribute this effect to insufficient capability of CCL to evaporate liquid water produced in the ORR in the domain close to the channel outlet. Flooding of this domain of the cell occurs due to positive feedback loop between the oxygen diffusivity and the rate of liquid water evaporation in the CCL. Decrease in D ox reduces the active volume of the CCL, which, in turn, lowers the total rate of water evaporation in the layer. In the flooded domain, the reaction penetration depth due to oxygen transport in the CCL appears to be close to the CCL thickness, which provides conditions for the development of oxygen transport -water evaporation instability.

Acknowledgments
T. Reshetenko gratefully acknowledges funding from US Office of Naval Research (N00014-15-1-0028) and US Army Research Office (W911NF-15-1-0188). The authors are thankful to Günter Randolf and Jack Huizingh for valuable help in the system operation and the Hawaiian Electric Company for ongoing support of the Hawaii Sustainable Energy Research Facility.