An Isothermal Model for Predicting Performance Loss in PEMFCs from Balance of Plant Leachates

A model is presented for the contamination of a proton exchange membrane fuel cell (PEMFC), including adsorption onto the Pt catalyst, absorption into the membrane, and ion exchange with ionomeric components. Model predictions for three sources of voltage loss account for the two-dimensional, time-dependent contamination along the channel and into the membrane. The model is developed by considering the well-known concepts of Langmuir adsorption, partition coefﬁcients, plug ﬂow reactors (PFRs), and dimensionless analysis. The phenomena are shown to be controlled by three important dimensionless groups: a Damk¨ohler number for the contamination reaction rate, a capacity ratio, and a coverage ratio for each contamination mechanism. These groups show how to scale ex situ equilibrium data for in situ predictions. The model predictions are shown to be reasonable when compared to in situ experiment data once ex situ data are used to provide reaction and equilibrium parameters. The predictions enable estimation of tolerance limits for leachates according to each mechanism. For typical parameters, the predicted voltage loss in the electrode ionomer by an ion-exchange mechanism shows slower reaction rates but greater performance losses than other mechanisms. © The Author(s) 2014. 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, provided the original work is properly cited. [DOI: 10.1149/2.0461414jes] All rights

Analysis of the potential for commercialization of proton exchange membrane fuel cells (PEMFCs) shows that the cost of balance of plant (BOP) materials can be as much as 30% of the system cost. 1 One opportunity to decrease these costs is to use off-the-shelf materials rather than custom-made materials, if leachates from less expensive materials would not affect performance and life-time. 2 This loss in performance is related to the contamination of the membrane 3 as well as the electrodes. Many previous studies for contamination [3][4][5][6][7][8][9][10][11][12][13] in PEMFCs have shown the sensitivity of performance to low levels of contamination. For example, cation leachates from gaskets or seals, from NH 3 in the fuel, or from catalyst metals [3][4][5][6][7][8][9] are well known contaminants of the membrane through an ion-exchange mechanism. Catalyst contamination by CO through an adsorption mechanism from feed gas contaminants is so well-known that those studies are too numerous to list here. Examples are also known for SO 2 , 10,11 for H 2 S, 12 and even for details of the reverse water-gas shift reaction of CO 2 . 13 The fundamental mechanism of contamination from organic contaminants such as aniline, ε-caprolactam, and diaminotoluene (DAT), [14][15][16][17][18] which have been identified as possible system contaminants, have been studied only briefly and recently, and it has been shown that their impact on PEMFC performance depends on their functional groups. For example, the ion-exchange reaction of the amine group in aniline with the Nafion membrane/ionomer causes a decrease in the number of available proton sites, thereby hindering the transport of protons and increasing ohmic losses. Aniline is also adsorbed on the carbonsupported platinum catalyst through interactions between the aromatic ring and the platinum metal and carbon supports. 16 A predictive and mechanistic model could aid in the interpretation of data and understanding of performance loss through contamination. Some of these mechanisms and the respective current responses were recently presented in Table II of Ref. 3, and this "library" of models is an excellent initial source for understanding in situ data exposed to changes in the contaminant concentration. As discussed there, 3 one goal for such a model would be the establishment of tolerance limits for exposure. If the model was based on the correct mechanisms, it would provide ideas for mitigation strategies. Another goal is to provide scaling functions that would allow translation of ex situ data into parameters useful for prediction of in situ losses from organic contaminants. This removes the need to use laboratory-scale in situ data for fitting the equilibrium or rate constants. Thus, we seek a contamination model which is widely applicable to the various types of contamination and which is readily scalable from bench-top to large-area cells. In this paper, we develop and verify a model based on fundamental contamination mechanisms of catalysts, electrode ionomers, and membranes in PEMFCs.
In previous PEMFC modeling studies of membrane/ionomer contamination, Okada et al. 19,20 developed a steady-state membrane model that correlated the cation impact on performance loss with water transport. Kienitz et al. 21 applied 1-D (i.e., anode to cathode) dilute solution theory to show how the cathode was affected by the mobility and distribution of cations in membranes as a function of current density, and Kienitz and Springer 22 extended the work by using concentrated solution theory and including the time dependency of this distribution. Weber and Delacourt 23 also developed a 1-D contamination model using concentrated solution theory, and thus both models required binary coefficients. Greszler et al. 24 included thermodynamic effects (such as the proton gradient between the anode and cathode) in their H 2 pump and H 2 /O 2 models. These previous models did not account for ion exchange physics in their models. Serincan et al. 25 included ion-exchange physics in their 1-D model with Langmuir isotherms. St-Pierre 3,26 proposed a time-dependent lumped capacitance ionomer contamination model that includes terms to characterize ion-exchange physics with a separation factor for the cation in the membrane. Note that the lumped capacitance models may a priori neglect the relationship of the distribution of the contaminant with the channel length, especially as area and channel dimensions change. In the model developed in this paper, we account for this distribution along the length of the channel that may be necessary as one changes from laboratoryscale lengths to full-scale cells with longer channels.
Models that address the absorption mechanism, which has been shown to exist for organics, 18 are non-existent in PEMFC literature. On the other hand, for the adsorption mechanism, many contamination models [27][28][29] have been developed and introduced. Typically, these consider only a single gaseous species, such as CO, SO 2 , NO 2 , or H 2 S, and the more detailed models are based on rate-determining reactions or steps (see Table II of Ref. 3) that produce unique in situ behaviors for step changes in the contaminant. For example, St-Pierre et al. 27 proposed surface-area-averaged or lumped-system adsorption contamination models that include kinetic and ohmic voltage losses. It should be noted that collecting data for experimental determination of the kinetic rate constants requires in situ experiments. Again, no distribution prediction for coverage by contaminant along the channel length was reported, and all of these models for contamination use an assumption of isothermal operation.
With this background, we propose below a model which can use equilibrium data from ex situ experiments, includes three fundamental contamination mechanisms, and considers a coverage or contamination distribution along the channel. In this paper, we focus on organic contaminants, which have been identified in leachates from BOP materials. These materials may be considered for use as structural materials or assembly aids in PEMFC systems. 14,15 The model was developed by considering the well-known chemical engineering concepts of Langmuir adsorption, partition coefficients, and plug flow reactors (PFRs). Finally, the model predictions are compared to in situ infusion experiment data by using parameters from ex situ data obtained with thin-film electrodes (TFE) and membrane ion-exchange and absorption experiments.

Model
The model is developed by considering Figure 1 and the concept of adsorption (or absorption or ion exchange). Figure 1 shows the channels, gas diffusion layer (GDL), and catalyst layer (CL) of a PEMFC. The dotted line represents the control volume of interest. Note that the control volume can be changed depending on the contamination mechanism. As shown in a cross-sectional view of the channel, contaminant A is convected along the channel in the x direction. The concentration along the channel decreases as contaminant A diffuses in the z direction (e.g., through the GDL) and ultimately adsorbs onto the CL or membrane.
The assumptions of the model are as follows: (1) the contaminant is dilute relative to the humid single phase gas stream; (2) only the x (along the channel) and the z (perpendicular to the membrane) directions are important; (3) negligible pressure drop along the channel; (4) isothermal conditions exist; (5) the effective diffusion coefficients in the porous GDL, ionomer, and membrane can be written using the MacMullin number; [30][31][32] (6) constant physical and transport parameters; (7) ideal gas law; (8) well-developed laminar flow in the channel; (9) the Butler-Volmer equation can be used to describe the cathode reaction with the Tafel assumption; (10) the rate of contamination is governed by the rate of mass transport of contaminant A by diffusion in liquid water, gas, or solid ionomer depending on the species and the saturation of the GDL, CL, or membrane; 3,26,33 (11) Langmuir adsorption isotherms apply for the adsorption behavior; (12) gas-phase transport of contaminants in the channel (although the contaminant is fed in the liquid state, it is nebulized prior to entering the cell in the experimental apparatus; thus, although salt cations are not present in the gaseous phase, their average bulk concentration in suspended droplets can be approximated as diluted by the gas stream and they can be transported with the gas stream (i.e., aerosol)) a ; (13) no interaction between contamination mechanisms; (14) neglect of migration and convection; (15) contaminant concentrations are equal in the liquid and gas phases because the concentration is very low; and (16) linear distributions of contaminants in the electrode and membrane.
Contamination by adsorption.-A reversible adsorption reaction of contaminant A onto the surface of the Pt catalyst (with S available sites) is expressed as A + Pt(S) ↔ A − Pt(AS). The material balance for contaminant A is made on a control volume in the channel with the above assumptions. The velocity of the fluid stream, u f , is in the x-direction through the channels and contaminant A diffuses in the z-direction (perpendicular to fluid flow). The rate of mass transfer by diffusion in the fluid phase is equal to the rate of reaction (adsorption) on the surface of Pt catalyst. Thus, the material balance on contaminant A is made on the volume (w + r ) · h · x over the time period from t to t + t. Application of the mean value theorems followed by the limiting conditions ( x and t→0) yields the following result: Note that assumption 2 neglects gradients in the y direction and that the ratio of rib width to GDL thickness is typically 2.5. This larger diffusion path in the y direction yields a smaller area correction than the (w+r)/w value of Eq. 1.
From the film diffusion model, 53 the rate of contaminant accumulation in the solid phase (on the Pt sites) is equal to that of contaminant transfer across the GDL according to the mass balance C * A is the concentration in the midpoint of electrode assumed to be in equilibrium with the Pt surface species, C AS , to be described in Eq. 5 below. Thus, the material balance equation for C A can be written as follows: The change in surface coverage by contaminant, C AS , with time is equal to the change in the equilibrium concentration of contaminant A, C * A , with time. Assumption 10 states that the rate of mass transfer by diffusion from the channel through the GDL is equal to the rate of reaction (i.e., adsorption) onto the surface of Pt catalyst and assumption 16 indicates that the surface concentration on the Pt catalyst can be represented by an average value at the midpoint of the catalyst layer (thus δ = z 1 +z 2 /2). Thus, As we mentioned before, the coverage and equilibrium concentration of contaminant A at the midpoint of the CL can be expressed by the Langmuir adsorption isotherm, Eq. 5.
Then, using the chain rule, a Since liquid droplets are governed by two-phase momentum equations, assumptions 1 and 12 may produce inaccuracies for the concentration and boundary condition at the GDL/channel interface. Note that assumption 10 requires the use of a diffusion coefficient for the appropriate phase or material.
we can obtain an expression for the time-dependent value of C * A that avoids the singularities that result from direct substitution of Eq. 5 into Eq. 4: [7] In summary, we have two coupled partial differential equations (PDEs) shown by Eq. 3 and 7, which require two initial conditions and one boundary condition. These could be written as follows: Defining dimensionless variables for distance (ξ), concentration ( A ), coverage by contaminant (θ), and time (τ) [10] and substituting into Eqs. 3 and 7 yield the dimensionless partial differential governing equations for contaminant A concentrations ( A , * A ). Solution of this and substitution of the dependent variables into Eq. 5 will yield the Pt coverage (θ).
Here, β is the geometric ratio of channel rib to width, and is the capacity ratio of Pt loading to inlet concentration of contaminant A. The Damköhler (Da) number is routinely used by chemical engineers as a measure of rate of reaction (i.e., to z-direction relative to convection to x-direction) at the inlet condition; and is the coverage ratio of contaminant A, which is associated with Langmuir parameters (B,K eq ).
= K eq C A0 , max = B K eq C A0 [13] The initial and boundary conditions are then defined again using the dimensionless variables. @ τ = 0, A = 0 and * Finally, using the dimensionless Langmuir isotherm Eq. 16 from Eq. 5, we can obtain the coverage change with time and distribution along the channel.
where z i is defined separately in the notation for the ionomer and the membrane.
In addition, the partition coefficient (K eq ) 33 and the coverage (χ) of acid sites by contaminant can be defined in terms of concentrations (e.g., C A , C AM , C H , and C HM ): Rearranging Eq. 19 by concentration gives us the following: Consistent with the derivation of the rate expression for the surface coverage, 34 we use the chain rule to remove C AM from the rate of change of C * A with time and write the following: [21] In summary, we have two coupled partial differential equations (PDEs), shown by Eqs. 3 and 21, which require the same two initial conditions and one boundary condition as Eqs. 8 and 9. Using the dimensionless groups above, we can write the equivalent dimensionless rate equation for Eq. 21: The right-hand side of Eq. 22 differs from Eq. 12 because the ion-exchange equilibrium is typically written from zero to one (i.e., equivalent molar fraction of cation) rather than from zero to a value of B for the coverage. 33 Here, β is again the geometric ratio of channel rib to width, is the molar ratio of acid sites to inlet concentration of contaminant A, Da is a measure of the rate of reaction (i.e., change in the z-direction relative to convection to x-direction) at the inlet condition, and is the partition coefficient for contaminant A. [23] Finally, using the dimensionless groups and rearranging Eqs. 18-20, we can obtain the coverage as a function of time and position along the channel: It may be useful to consider the diffusion length and the effective diffusion coefficient of contaminant A (D eff ) for the various penetration lengths of Figure 1. If we approximate those lengths by assuming series rather than parallel diffusion thorough the GDL, electrode, and membrane), we can obtain an effective diffusion coefficient: for membrane contamination [26] where δ = z G DL + z C L + z mem /2 for membrane contamination and δ = z G DL + z C L /2 for electrode contamination. This approximation assumes that the contaminant A diffuses through the GDL, then contaminates the ionomer in the electrode, and then diffuses through the membrane/catalyst layer interface, causing membrane contamination.
To reach an approximation of the diffusion through the ionomer of the electrodes, one could assume that the diffusion coefficient through the ionomer is equal to that through a Nafion membrane, which is generally measured by mass uptake, streaming potential, radiotracer or NMR methods. [35][36][37][38][39][40] Finally, an average coverage along the channel can be obtained by integration of the local coverage along the length of the channel, 0 < ξ < 1: Estimation of dilution by water generation/transport.-The model assumes isothermal conditions, as it is limited to an exit relative humidity of less than 100% and does not consider phase change. However as water is generated on the cathode by the ORR, one may be interested in estimating the extent to which the concentrations of contaminant A in the CL and channel change depending on operating conditions. If we assume (1) the current to be uniform over the length (i.e., x-direction of the channel) and (2) no water accumulation at the interface, the flux balance at the interface on the cathode for water vapor can be written as follows: Here, the α represents the ratio of net water flux to proton flux and combines the back diffusion from the cathode to the anode and the electro-osmotic drag (e.g., migration) from the anode to the cathode. 41 A negative α indicates that the back diffusion of water vapor from the cathode to the anode is dominant. As shown in Table A1, previous studies reported that the α is −0.053 for 25% RH and 0.28 for 50% RH with i = 0.4A/cm 2 for similar types of MEA. 42,43 From these values, one can estimate that a maximum dilution produces a 10% error in predictions using the simple model.
Predicting voltage change.-The total voltage change by contamination can be written as follows: Eq. 29 can be simplified if we assume no change in the OCV from the contamination, as observed in the experiments of Ref. 34, and negligible anodic overpotential change ( η a ). Thus, For membrane contamination, we can use Ohm's law (i.e., V = iR):.
Here, conductivity (κ) can be obtained from ex situ conductivity studies 46,47 as a function of coverage. Generally, the conductivity (κ) and coverage (χ) have the following exponential relationship at constant RHs: 34 ln κ = aχ + b, where a and b are constants [32] Thus, Eq. 31 can be written in terms of coverage (χ) by contaminant in the membrane for specified RHs: [33] To explain the change in cathodic overpotential ( η c ) with contamination, we used the Butler-Volmer equation for reaction j to explain the dependency between surface overpotential (η s ) and current density (i) on the electrode, 44 and we assume large cathodic overpotentials for the oxygen reduction reaction (ORR) so that For the ORR (i.e., 2H 2 O ↔ O 2 + 4H + + 4e − ), Eq. 34 can be written for pristine and contaminated electrodes as: where C H + and C O 2 are the surface concentrations of H + and O 2 , respectively, on the platinum. Finally, we can obtain an expression for the voltage loss by subtracting Eq. 36 from Eq. 35 to yield the following: Eq. 37 shows two contributions to voltage loss. One is from ECSA change by Pt contamination, and the other is from the proton concentration change by ionomer contamination. Note that q ij values of 1 for Pt contamination and of 1.7(±0.3) for ionomer contamination of the ORR have been measured directly for film electrodes (i.e., Pt/C + ionomer), and a previous report on Pt electrodes in an acid solution gave the reaction order for acid solutions as 1.5. 45 In addition, Tafel slopes of −60 mV/decade and −120 mV/decade for V(iR corrected)>0.85 V and V(iR corrected)<0.85 V, respectively, for the ORR on Pt/C electrode were used. 46 In summary, the total voltage loss from contamination has three contributions from the Pt, ionomer, and membrane. Figure 2 shows the predicted voltage change (i.e., loss) as a function of coverage of either the Pt (θ) or acid sites (χ). To calculate the V mem, one needs data for Eq. 32 at a given RH to calculate the line of (c) in Figure 2. Here, we used Na-exchange membrane data at 60% RH and 0.2A/cm 247 to obtain the constants a and b in Eqs. 32 and 33. Other values of a and b may yield a steeper slope for V mem versus χ, and constants a and b depend on the contaminant species and operating conditions (i.e., RH, I).

Model Predictions
To solve the two hyperbolic coupled PDEs, for each contamination mechanism, we used a ( x 2 ) finite difference method for the length direction and the Gear algorithm for the time derivative set of discretized ODEs. The results reported below are independent of the size of both the time step and the spatial discretization.
Values for dimensionless numbers.-Before proceeding to the predictions of the model, we will discuss the selection of the dimensionless numbers (i.e., Da, ψ, ). We selected the dimensionless numbers by considering the general operating variables (T, i, RH, backpressure, stoic, C A0 ), the system geometry (w, r, h), and the physical constants (N m , δ, D, K, B) used for in situ infusion experiments. 2,[14][15][16][17][18][30][31][32]34,52 Tables I to III list the variables for each contamination source (i.e., Pt = 1, ionomer = 2, membrane = 3). For the purpose of exploring the parameter space for the model predictions, we choose a 2 2 factorial design for each source in terms of the dimensionless numbers (i.e., high and low) as summarized in Table IV for each mechanism. For 1, 2, and 3 (capacity ratios), the Pt catalyst and H + ionomer and membrane concentrations change with feed contaminant concentrations for a specified membrane electrode assembly (MEA). For the feed concentrations, we selected values from the analytical study (i.e., GC and ICP) of BOP leachates. 14,15 The Pt surface concentration can be fixed for the same Pt loading and electrochemical surface area (ECSA = 68 m 2 /g pt ). The ion exchange capacity (IEC = 0.97 mmeq/g) of Nafion was used for the calculation of the H + concentration in the ionomer and membrane. The Da number (ratio of reaction rate to convection transport down the channel) depends on the effective diffusion rate toward the electrode or membrane (in the z-direction, per Figure 1) and the dimensions and geometry of the system components (e.g., 50 cm 2 single-cell hardware which was assembled with bipolar plates with triple serpentine gas channels). Finally, the (coverage ratio) can be obtained from ex situ adsorption and ion-exchange isotherms of Pt and the Nafion membrane for selected compounds. [16][17][18] Details of the calculations for Table IV can be found in Ref. 34.
Pt contamination.- Figure 3a shows the predictions for the dimensionless concentrations of the contaminant at the entrance (ξ = 0) of the channel ( A ) and at the midpoint of the CL ( A * ). These are predicted for the case where 1 = 6.0, B = 0.6, and Da 1 and 1 are varied as shown in Table IV. The predicted value of the concentration of contaminant in the channel has a step change from 0 to 1 at ξ = 0 according to the initial and boundary conditions. A * approaches 1.0 from an initial value of 0.0 at different rates for different combinations of Da 1 and 1 . The difference between these concentrations  is the driving force in the PDEs. The fastest concentration change is for (Da 1 , 1 ) = HH, the change is moderate for HL (or LH), and the change is negligible for the LL values over the dimensionless time shown. That is, the driving force in Eqs. 11 and 12 decreases rapidly for a high diffusion of A through the GDL relative to convection down the channel, and it decreases rapidly for a large value of initial moles of contaminant relative to the active sites of Pt (i.e., HH in Table IV and Figure 3a). Thus, steady-state is approached quickly for the HH conditions. This approach to steady-state is slower with LH and HL conditions and slowest for LL conditions. Again, these are all for ξ = 0 (at the entrance). Figures 3b and 3c show the average coverage and voltage change corresponding to the change in driving force shown in 3a as a function of Da and . Again, the fast decay of Pt availability is observed for the high dimensionless numbers of Da 1 and 1 . It is interesting to note that all cases reached the same coverage loss at steady-state because 1 is fixed for all combinations. Figure 4a shows the effect of varying 1 , which represents the dimensionless concentration of contaminant A at ξ = 0 in the channel ( A ) and in the CL ( A * ) predicted by the model for the case where Da 1, 1 , and B are fixed  Table IV. The average coverage and voltage change are also shown in Figures 4b and 4c. The higher 1 shows severe loss of Pt availability and voltage due to the larger coverage ratio. Here, the slopes of the predicted voltage curves for each case are the same, as the contamination reaction rates are equal because Da 1 and 1 are fixed. Note that the maximum coverage (B) in Figures 3b and 4b depends on the specific contaminant, as described below for the ex situ isotherms.
In summary, the analysis with dimensionless numbers for Pt contamination reveals that a high Da 1 and 1 yield a fast contamination reaction, and a high 1 causes a greater coverage ratio. If one uses these predictions to establish tolerance levels, then the coverage (i.e., ECSA loss) should be less than 0.35 if one chooses to maintain a voltage loss from the beginning of life (BOL) voltage of less than 50 mV, and then only if the Pt contamination alone is relevant (see Figure 2). We specifically mention BOL because the contamination must be addressed separately from the long-term degradation of the Pt/C. Thus, only the (Da 1 , 1 ) = LL case satisfies these tolerance limits until the dimensionless time of 4 × 10 5 (≈50 h). Clearly, lower Da 1 and 1 are required to mitigate the impact of Pt contamination. Highcurrent and high-RH conditions are preferred because of the high flow rate/dilution effect. Other changes in geometry and capacity, such as higher loadings or more moles of ionomer, that change Da 1 and 1 may not yield desirable higher economic tolerances.
Ionomer contamination.-Next, we consider ionomer contamination and voltage change resulting from the ion-exchange reaction. Figure 5a shows the prediction for the dimensionless concentration   Table IV. Figures 5a and 5b show predictions for average coverage and voltage change at various values of Da 2 and 2 corresponding to those used in Figure 5a. The trends for this ionomer contamination are similar to Pt contamination, where greater values of Da 2 and 2 have a significant impact on contamination. It is interesting to note that the voltage loss caused by ionomer contamination is more sensitive to coverage than that caused by Pt contamination, as shown by the exponent for H + in the ionomer of the electrode being 1.7, as obtained from experimental data. Recall that the exponent of 1.5 was obtained for a crystalline Pt crystalline electrode. The maximum coverage for the ion-exchange mechanism is 1.0 (= fully exchanged) because it depends on the molar ratio of contaminants to acid sites. 33 However, the maximum Pt coverage is typically less than 1.0, corresponding to the Langmuir behavior. 56 It is also interesting to note that, as shown in Figure 6, affects the rate coverage for the membrane and ionomer, whereas the value of B affects the final value for Pt coverage (see Figure 5). This is because the maximum coverage value is 1.0 for the ion exchange process and typically <1.0 for adsorption.
In summary, the analysis of dimensionless numbers for ionomer contamination also reveals that higher Da 2, 2 , and 2 lead to rapid contamination and voltage change, and that the reaction continues until the coverage reaches 1.0. These predictions show that coverage of less than 0.2 is required for the BOL voltage loss to be less than 50 mV for ionomer contamination, if only the ionomer contamination is relevant (see Figure 2) and all other contaminations and interactions are negligible. Another possible ionomer contamination mechanism is absorption of contaminant A into the ionomer in the electrode. In this case, we can use the same design equations and analysis as with Pt adsorption cases, so that the relevant isotherms are the Langmuir absorption isotherms, and we substitute these in the design equations to predict the impact of ionomer contamination. To quantify the voltage loss from proton concentration change by ionomer contamination, ex situ conductivity data for absorbed compounds are needed to compare to those of ion-exchange compounds (e.g., Na + ). For example, conductivity losses for diethylene glycol monoethyl ether (DEGEE) absorbed into NRE211 are equal to 30% of those of Na-exchanged NRE211. Thus, maximum coverage (i.e., proton concentration change in Eq. 37) of DEGEE for the absorption mechanism can be estimated as 0.3 to quantify the voltage loss of the ionomer in the electrode by the model. 18,34 Membrane contamination.-In a similar fashion to that described above, the analysis of dimensionless numbers for the membrane contamination by ion-exchange reaction also shows that higher Da 3, 3 and 3 yield the fastest contamination and voltage loss, as shown in Figure 6d. The interesting aspect of membrane contamination is that 3 is a small number due to the large number of moles of ionomer in a typical membrane. Thus, the contamination and voltage loss are slower for the membrane at a given feed rate than for the catalyst Table V. Effect of Da 1 and Da 2 on coverage and V distributions for the a) adsorption mechanism with 1 = 0.14 and 1 = 6, b) ion-exchange mechanism with 2 = 0.003 and 2 = 1.0. or ionomer in the electrode. There is a similar slowing of membrane contamination relative to the electrode contamination because the diffusion lengths are increased as the reaction plane moves toward the center of the membrane (i.e., Da 1,2 > Da 3 ). The membrane contamination by ion-exchange reaction continues until the acid sites are fully exchanged, which leads to a total loss of conductivity.
The predicted results show that membrane contamination is more tolerant than ionomer contamination, and is the last to be affected by the voltage loss. However, a more complex model is necessary, as we believe that the two mechanisms interact in that a greater capacity in the membrane can decrease the effect of ionomer contamination because the membrane will withdraw contaminant from the ionomer. The equilibrium coverage would then be divided between the region of less capacity (i.e., the ionomer) and the region of greater capacity (i.e., the membrane) according to the diffusion coefficient of contaminant. Note that the predicted results in Figures 3 to 6 are prepared by assuming non-interaction between contamination sinks (i.e., ionomer and membrane). We will discuss in the next section this interaction of contamination between the ionomer and the membrane.
Distributions along the channel.- Table V shows the local coverage and voltage distributions within various fractions of the dimensionless length (ξ) for three values of increasing Da for the adsorption and ion exchange mechanisms. For both mechanisms, Table V shows that the concentration distribution becomes more important as the Da number increases. Thus, the concentration distribution can become important for thin GDL, CL, or membrane thicknesses, when there is less ionomer content in the electrode, or when the contaminant has a large diffusion coefficient. Table VI shows the local coverage and voltage distributions for contaminant A across the dimensionless length (ξ) as increases for each mechanism (i.e., adsorption and ion exchange). Like the Da, both mechanisms show that the distribution becomes more important with an increase to (i.e., 0.4 to 0.05 to 0.01 mg/cm 2 Pt loadings). Note that the prediction for the local voltage for the ion exchange takes the ionomer in the electrode into consideration. Thus, the concentration distribution can become important for large feed concentrations, a low number of Pt and acid sites, or a thin membrane thickness. For example, cases (a) through (c) for the adsorption mechanism have identical leaching conditions of 1280 ppm of DEGEE with a flow rate of 0.23 ccm, which leads to a contaminant concentration of 1.42 × 10 −8 mol/cm 3 in the channel, when the cell operates at T = 80 • C, RH = 32%, stoics = 2.0, and back pressure = 150 kPa. Only Pt catalyst loadings are assumed to be different for the cases (i.e., 0.4 mg/cm 2 , 0.05 mg/cm 2 , 0.01 mg/cm 2 ). Again, the concentration distribution can become significant for ultra-low loading of Pt, and it will reduce the durability of PEMFCs.

Comparison to Experimental Data
In this section, the contamination model presented above is compared with in situ experimental data. To make these comparisons, we demonstrate the impact of the parameters of the contaminant on the dimensionless numbers. Thus, experimental designs similar to those shown in Figures 3-6 may not be possible. The parameter values used were obtained from the literature and from separate ex situ experiments. [16][17][18]34,55 For the comparisons, we chose data for organic compounds, DEGEE and 2,6-diaminotoluene (2,6-DAT), which shows all three contamination mechanisms (adsorption onto the Pt, absorption into the ionomer, and ion-exchange reaction), and for sodium (Na + ), which only drives the ion-exchange phenomenon.
The conditions of the DEGEE, 2,6-DAT and Na + infusion experiment are listed in Table VII and detailed in Refs. 34 and 52. The calculated values of Da and for each mechanism are shown in Tables VIII and IX with subscripts corresponding to 1 = Pt, 2 = ionomer, and 3 = membrane. For the calculation, ECSA = 68.0 m 2 /g pt is used for each case, and the effective diffusion coefficient is obtained from Eq. 25 and 26 assuming that half of the ionomer thickness is 2 μm. (see Appendix). The diffusion coefficient of water vapor through the GDL was used for gas (i.e., water vapor+air+contaminant) diffusion through GDL, and the diffusion coefficients of methanol, ammonium and Na + in Nafion membrane were used to represent the diffusion of DEGEE, 2,6-DAT, and Na + in the ionomer and membrane, respectively, because of the lack of data. [35][36][37][38][39][40][49][50][51] Scaling from ex situ experiments to in situ predictions.-Our goal for the model is to predict voltage losses using ex situ measurements. For the Pt coverage, our data indicates that equilibrium is obtained   55 This is consistent with the model formulation using the equilibrium constants. However, for the model developed here (in situ PEMFCs), we are concerned with transport of contaminants through the gas phase, 34,52 per assumption (12). Because the concentrations are very low in the gas phase, the liquid-vapor equilibrium curve extrapolates to satisfy assumption (15), that vapor and liquid phase concentrations are equal. Furthermore, to establish the relationship between the in situ gas-phase coverage and the coverage measured on an ex situ TFE using contaminants in the liquid phase, we match the dimensionless numbers corresponding to the capacity ratios as shown in Eq. 38 to scale the concentrations: In Eq. 38, the in situ moles available are calculated using the volume of the channel which corresponds to the product of the electrode area (A in ) and the channel height (h in ). This corresponds to the moles available for contamination since the gas diffusion coefficients are large and yield well-mixed conditions in the channel, so that gradients only exist in the GDL. For the ex situ moles available, we use the product of the liquid phase contaminant concentration, thin film electrode geometric area (A ex ), and the diffusion layer (h ex = y) corresponding to mass transfer from an infinite stagnant medium to a plane electrode. Note that since there is a concentration gradient in the diffusion layer, we use C A0,ex = C A0,bulk,RDE /2. Then, to determine the diffusion layer, y, of the TFE, we use Eq. 39, 53 where ς = dimensionless diffusion layer, D R = diffusion coefficient of reactant, and t = time. Thus, the thickness of the diffusion layer of the TFE is 0.04 cm for our ex situ condition. From Eq. 38, we can relate the in situ and ex situ concentrations with a scaling factor of 40 for a Pt loading of 0.4 mg/cm 2 ). Details of the calculations are shown in the Appendix. This scaling factor indicates that the concentrations of the adsorption isotherm shown in Figure 7 must be divided by 40 to apply to the in situ data for calculation of the coverage (or, conversely, the in situ data, C A0, in, can be multiplied by 40 to use the ex situ adsorption isotherm).
C A0,in = C A0,ex /40 or C A0,ex = 40C A0,in [40] For example, consider a contaminant concentration in the channel, C A0, in = 1.66 × 10 −9 mol/cm 3 , which corresponds to a molar feed rate to the cathode of 8.89 × 10 −9 mol/s for conditions of 10 A, stoic 2.0, 32% inlet RH, and 80 • C in a 50 cm 2 cell. The scaled value is 40 C A0, in , or 6.63 × 10 −8 mol/cm 3 , which can be combined with the ex situ adsorption isotherm to obtain a surface coverage of 0.42. Note that the predicted coverage for the scaled and unscaled concentrations is significantly different, as shown in Figure 7. The predicted coverage for the unscaled C A0 (i.e., 1.66 × 10 −9 mol/cm 3 ) is 0.18, and the coverage for the scaled C A0 (i.e., 6.63 × 10 −8 mol/cm 3 ) is 0.42. Alternatively, one could adjust the equilibrium constant. This calculation corresponds to case (c) in Table VII. We will discuss more results and analysis for Table VII in the next sections. For the ion-exchange mechanism, the scaling factor for the electrode and the membrane can be determined without considering the diffusion layer by Eq. 41, since the thickness of each membrane section is known.
Adsorption and absorption mechanisms.- Figure 8 shows the prediction results as lines for comparison with in situ experiment results 34 for the three cases in Table VII corresponding to changes in feed concentration and current density. Comparison of the model prediction with in situ experiment data for all cases shows that the predicted change in coverage, ECSA, and voltage loss (iR corrected) are reasonable (see Table VII and Figure 8). Table VIII shows that the 1 and 2 for cases (a) and (b) are almost equal because the feed concentrations in the channel are close. Thus, the model predicts the same loss of ECSA (i.e., coverage). Prediction results for coverage seem reasonable on the basis of experimental measurements of in situ cyclic voltammetry (CV) at the end of infusion. Additionally, the Da 1 and Da 2 for cases (b) and (c) in Table VIII are equal because they use the same species, so that D eff is constant, and operating conditions, such as RH, i, T, and stoics. Note that for the case (a) in Table VIII, i = 1.0A/cm 2 was higher than that of cases (b) and (c) (i.e., i = 0.2A/cm 2 ); and for a fixed stoic of 2.0, the u f for case (a) is larger than fort cases (b) and (c). The larger u f decreases the Da number. Likewise, the model predictions with the higher Da and yield the fastest contamination. Also, the higher causes a higher coverage ratio, resulting in a greater voltage loss at steady-state.    For the case (a) in Figure 8, the contribution from each of the contamination locations (i.e., Pt, ionomer, and membrane) is shown in Figure 9. The voltage losses by Pt (30%) and membrane contamination (10%) have relatively small impacts on the total loss, whereas the ionomer contamination (60%) has the most significant impact. As discussed earlier, the predicted voltage loss for each source was different even when the coverage was equal (see Figure 2). The 80 mV of total voltage loss that occurs within 5 h would be unacceptable. However, we could observe that the ionomer contamination by the absorption of DEGEE is recoverable, because the DEGEE does not have any cationic functional groups. Hence, there should be no ion exchange with the acid groups in the ionomer. This assumption was confirmed by in situ experiments based on the HFR change during contamination and during a recovery experiment (shown in Ref. 34).
Ion-exchange mechanism.- Figure 10a shows the prediction results as lines compared to experiment results shown as dots for an NaOH infusion of 30 ppm at a flow rate of 0.03 cm 3 /min. The com-  parison of the model prediction with in situ experiment data shows that the predicted voltage loss is reasonable for the initial infusion region (0-70 h) and that the voltage loss by ionomer contamination is significantly larger than that by the other sources. It is noteworthy that the Pt contamination can be neglected for Na + contamination. A previous study 5 also reported that cation contamination on the electrode did not affect the Pt catalyst but speculated that there was an effect in the ionomer. In situ CV measurements during the 80 h of NaOH infusion (see the vertical data at 21 h, 45 h, and 78 h in Figure 10a) support the hypothesis of no change in available Pt concentration (i.e., negligible ECSA changes <5%), as shown in Figure 10b (See reference 34 for experimental details). The membrane contamination also shows a relatively small impact on PEMFC performance relative to ionomer contamination. In addition, we believe that the contamination by the ion exchange cannot be recovered fully, as the Na + strongly interacts with the sulfonic acid groups in ionomeric components (e.g., Nafion).
For the 2,6-DAT contamination, we found that the adsorption on Pt and ion-exchange reactions in the ionomeric components are possible mechanisms considering its functional groups (i.e., aromatic benzene ring and amine groups) from the ex situ studies, while Na + only showed the ion-exchange mechanism, as discussed above. Figure 11 shows the prediction results as lines compared to the experimental results shown as dots for an infusion of 128 ppm 2,6-DAT with a flow rate of 0.03 cm 3 /min (=1.82 × 10 −10 mol/cm 3 ). The comparison of the model prediction with in situ experiment data shows that the predicted voltage loss is reasonable for the 100 h infusion region. Again, the voltage loss by ionomer contamination is significantly larger than that from other sources (i.e., Pt and membrane), but membrane contamination becomes important as 2,6-DAT contamination continues (i.e., t = 20-100 h). It should be noted that the point at t = 20 h indicates the reaction plane movement toward the membrane from the ionomer, resulting in the slope change. For example, a discrete peak point at 20 h rather than the smooth curve represents the time at which the governing mechanism changes from the ionomer to the membrane. The reason for this is that as the contamination reaction plane moves toward the membrane from the electrode, the effective diffusion length is increased, lowering the Da number. We believe that the initial contamination is governed by ionomer contamination, and that membrane contamination becomes more important as contamination continues. In situ CV measurement at 100 h of 2,6-DAT infusion shows good agreement with the prediction (i.e, −18.0% and −16.1%), as shown in Table VII. In addition, the contamination by ion exchange cannot be fully recovered after operating the cell with the contaminant-free normal operating conditions (i.e., t = 100-120 h).

Conclusions
A simple isothermal time-dependent 2-D model that accounts for three contamination mechanisms in PEMFCs in terms of measurable properties has been presented. The model equations enable prediction of performance loss and distribution (i.e., voltage distribution) for PEMFCs. The contamination phenomena are shown to be controlled by three important dimensionless groups: a Damköhler number for the contamination reaction rate (Da), a capacity ratio ( ), and a coverage ratio ( ) for each mechanism. The experiments and model predictions are compared in order to evaluate the model, and the predictions are shown to be reasonable. A method for scaling the equilibrium constants between liquid phase ex situ and gas phase in situ has been presented. The model can be used to provide tolerance limits for contamination and to predict voltage losses attributable to adsorption on the Pt/C catalysts and absorption into the ionomer. The reaction plane concept does not adequately account for the interaction of the acid sites in the ionomer of the electrode with the membrane in this model, but future work will address their interaction. Table A1 shows the initial (i.e., feed) and the change of total moles of species on the cathode. The maximum RH at the inlet and outlet of the channels is also listed in Table A1 at I = 10 A, P = 150 kPa, stoic = 2.0. The calculation shows that the dilution effect on the contaminant by the water can be negligible (i.e., n T n T 0 = 0.1), and that the water saturation at the interface of the membranes on the cathode may start over RH = 50% (inlet).The water vapor concentration and RH can be defined thus: The molar flow rate for each species at the inlet of the channel can be defined thus: From Eq. 28, the difference of water vapor concentrations between the interface and the channel showed 10% (i.e., C * w −Cw Cw × 100) with selected α (shown in Table A1). In other words, we can speculate that the 10% of RH difference exists between the interface (electrode and membrane) and the channel. Thus, water saturation may start at the interface of the membrane over RH = 90% of the channel condition. Thus, we can assume that the presented model is applicable for feed RHs below 50% under the above conditions.  where h ex = 0.04 cm, h in = 0.084 cm, A ex = 0.2475 cm 2 , A in = 50 cm 2 , Pt loading of ex-situ = 0.02 mg/cm 2 , in-situ = 0.4 mg/cm 2 For the ex situ available moles, we use the product of the liquid phase contaminant concentration, thin film electrode geometric area, and the diffusion boundary layer. Note that since there is a concentration gradient in the diffusion layer, we use C A0,ex /2. Thus, the scaling factor is obtained as 40 (i.e., 20 × 2) by Eq. B3. This scaling factor indicates that the concentrations of adsorption isotherm shown in Figure 7 must either be divided by 40 to be applied to the in situ data or else multiplied by the 40 to use ex situ isotherms.

B. Scaling Factor Calculation
The scaling factor is not dependent on the reaction order of 1.7. The scaling factor was obtained by the comparison between the ex-situ thin film electrode experiment and the in-situ single cell infusion experiment.

C. Ionomer Thickness
We assume and use the following variables to calculate the ionomer thickness: carbon supported 45.5 wt% Pt catalyst, Pt loading is 0.4 mg/cm 2 , a 50 cm 2 single cell, ionomer to carbon (I/C) ratio is 1.0, ionomer density (ρ) is 2.01 g/cm 3