PDF
Mass Transport and Electrochemical Reaction in a Fuel Cell Cathode
Introduction
One of the more important aspects of fuel cell modeling is the mass transport through the gas diffusion and reactive layers. Gas concentration gradients may often be quite large and are strongly coupled to the reactions that take place.
Figure 1 shows an example 3D geometry of a cathode from a fuel cell with perforated current collectors. This geometry configuration can be used for self-breathing cathodes or in small experimental cells. Due to the perforation layout, a 3D model is needed in the study of the mass transport, current, and reaction distributions.
Figure 1: A fuel cell cathode with a perforated current collector.
The model couples this mass transport to a concentration-dependent Butler–Volmer electrochemical kinetic expression in a porous gas diffusion electrode (the cathode). Darcy’s law is used to define the convective velocity in the porous gas diffusion electrode, whereas diffusion is modeled using the Maxwell–Stefan equations. A note here is that the molar fractions of the reactants and products, that is, oxygen and water vapor, are typically large (>10%), which makes Fickian diffusion an inappropriate assumption for modeling the diffusive mass transport.
For a more detailed description for how to build this model, including screen shots, see the Introduction to the Fuel Cell & Electrolyzer Module book.
The electrochemical reaction for a PEM fuel cell to produce electrical energy is given by:
where denotes the standard equilibrium potential of the cell reaction, assuming all reactants reacting in the gas phase at atmospheric pressure.
At the anode Hydrogen Oxidation Reaction (HOR) yield protons:
water is produced via Oxygen Reduction Reaction (ORR):
Model Definition
Figure 2 shows details for a unit cell, cut out from Figure 1. (In this case, the combination of a circular orifice and square unit cell eliminates the possibility of approximating the geometry with a rotationally symmetric model.) The circular hole in the collector acts as the inlet where the gas enters the modeling domain, and at this boundary the gas mixture composition and pressure is known. The upper and lower rectangular domains are the reaction-zone gas diffusion electrodes. They consist of a three-phase porous structure that contains the feed-gas mixture, an electronically-conducting material covered with an electrocatalyst, and an ionically-conducting electrolyte.
Figure 2: The modeled fuel cell unit cell. The quarter circle part of the top boundary is the surface of the cathode that is open to the feed gas inlet, while the rest of the top surface sits flush against a metal current collector. In the unit cell, the top domain is the porous cathode, the middle domain is the membrane, and the bottom domain is the porous anode.
The middle domain corresponds to a solid electrolyte membrane, ionically interconnecting the two electrodes of the fuel cell. No reaction takes place in this domain and the current is conducted ionically. In addition, there are no pores present to allow gas to flow, nor any material present for electronic current conduction.
The gas diffusion electrodes are 0.075 mm thick, as is the electrolyte layer. The unit cell is 1.5-by-1.5 mm in surface, and the gas inlet hole has a radius of 1.0 mm.
The Hydrogen Fuel Cell interface models the electronic and ionic current balances and solves for the potentials and in the electrode and electrolyte phases, respectively. The anode side of the cell is grounded, whereas the current collector boundary at the cathode is set to a cell potential value.
Mass transport and fluid flow are also modeled using the Hydrogen Fuel Cell interface. The species (mass) transport is modeled by the Maxwell–Stefan equations for the mass fractions of oxygen, water and nitrogen in the O2 gas phase. Mass transport is solved for in the cathode gas diffusion electrode domain only. Similarly, the pressure and the resulting velocity vector is solved for in the cathode gas diffusion electrode domain only using Darcy’s Law. As boundary conditions, inlet molar fractions are set for the three gas species corresponding to a humidified air mixture at 90% relative humidity at atmospheric pressure.
No mass or momentum transport effects are expected to occur at the hydrogen anode side. The partial pressure of hydrogen is set to be constant in the anode domain.
The cell operates at 70ºC. Reference equilibrium potentials for the higher temperature, the reference state, for each reaction are calculated automatically by the Hydrogen Fuel Cell interface from the standard free energies of formation (ΔH) and reaction entropies (ΔS) according to
where T denotes the operating temperature, n the number of electrons participating in the electrode reaction and F Faraday’s constant.
Generally, the equilibrium potentials of the electrode reactions will depend on the local partial pressures of the reacting species according to the Nernst Equation:
where νi are the stoichiometric coefficients of the reacting species.
The cathode electrode kinetics of the cathode are defined using a Butler–Volmer type of expression according to
where pi is the partial pressure of the reacting species, pref = 1 atm is the reference pressure and ηref, the overpotential with respect to the reference state, is defined as
.
The local current density expression in the cathode is multiplied by a specific area of 109 m2/m3 to create a volumetric current source term in the electrode domain. Assuming ideal kinetics according to the mass action law gives that .
For the anode domain, the kinetics is assumed to be so fast that a linearized Butler–Volmer expression may be used.
assuming . Also, the local current density expression in the anode is multiplied by a specific area of 109 m2/m3 to create a volumetric current source term in the electrode domain. The overpotential in the anode is defined as
.
In the first part of the model instructions below, a secondary (not concentration dependent) current distribution is modeled. In the second part, mass and momentum transport is incorporated in the O2 gas phase mixture (cathode domain), using Maxwell-Stefan diffusion and Darcy’s Law, respectively. In both parts of the tutorial, the model is solved for a range of cell potential values (0.5 V to 1 V in steps of 0.1 V) by the use of an auxiliary sweep in the stationary solver.
Results and Discussion
Figure 3: Polarization plot.
Figure 3 shows the polarization plot for the two scenarios investigated: limited and unlimited O2 gas phase transport. It can be seen that higher average cell current densities are achieved for the unlimited O2 gas phase transport scenario (that is, when no mass and momentum transport limitations are present).
Note that the plots and discussion in the rest of this section correspond to the limited O2 gas phase transport scenario, where diffusion and flow (in the cathode domain) has been considered, coupled to charge transport and the electrochemical reactions.
Figure 4: Mole fraction of oxygen at cell voltage of 0.7 V.
Figure 4 shows the oxygen mole fraction at a cell voltage of 0.7 V. The figure shows that mole fraction variations are small along the thickness of the cathode, while they are substantially larger along the electrode’s width.
Figure 5 shows the pressure and gas velocity streamlines in the porous cathode at the same cell voltage. There is a significant velocity peak at the edge of the inlet orifice. This is caused by the contributions of the reactive layer underneath the current collector because in this region the convective flux dominates the mass transport. The gas flows from the interior of the cell toward the circular hole. The reason for this is the oxygen reduction reaction, which creates two water gas molecules per oxygen molecule entering the cell, and this water must then leave the cell.
Figure 5: Pressure and velocity for the gas phase in the cathode’s porous reactive layer at cell voltage of 0.7 V.
The electrochemical reaction rate, represented by the local current density, is related to both the local overpotential and oxygen concentration in the cathode domain. Figure 6 depicts the local overpotential (at cell voltage of 0.7 V), which gets more negative toward the electrolyte domain.
The combination of the overpotential and oxygen concentration distributions will result in a highly uneven reaction rate in the reactive layer. One way to study the distribution of the reaction rate is to plot the ionic current density at the bottom boundary of the membrane layer. Figure 7 shows such a plot at cell voltage of 0.7 V.
Figure 6: Local overpotential in the cathode reactive layer at cell voltage of 0.7 V.
The current-density distribution shows that the variations are rather large. The reaction rate and the current production are higher beneath the orifice and decrease as the distance to the gas inlet increases. This means that the mass transport of reactant dictates the electrode’s efficiency for this design at these particular conditions.
Figure 7: Current density perpendicular to the lower membrane boundary at cell voltage of 0.7 V.
Figure 8 shows a plot of the cell-integrated power losses associated with the different processes in the cell. For this cell, the power losses due to kinetic activation of the oxygen reduction reaction dominate the losses, followed by the electrolyte transport losses.
By dividing the loss in power by the cell current, we can also compute the corresponding cell overpotentials, as shown in Figure 9. The oxygen reduction overpotential is highly nonlinear with respect to the cell current density. This stems from the Butler-Volmer kinetics deployed in the model. The electrolyte transport overpotential has a more linear behavior with respect to the cell current density.
Figure 8: Cell-integrated power losses
Figure 9: Cell overpotentials.
Application Library path: Fuel_Cell_and_Electrolyzer_Module/Fuel_Cells/fuel_cell_cathode
Modeling Instructions
From the File menu, choose New.
New
In the New window, click  Model Wizard.
Model Wizard
1
In the Model Wizard window, click  3D.
2
In the Select Physics tree, select Electrochemistry > Hydrogen Fuel Cells > Proton Exchange Membrane (fc).
3
Click Add.
4
Click  Study.
5
In the Select Study tree, select Preset Studies for Selected Physics Interfaces > Stationary with Initialization.
6
Click  Done.
Global Definitions
Parameters 1
Load some model parameters from a text file.
1
In the Model Builder window, under Global Definitions click Parameters 1.
2
In the Settings window for Parameters, locate the Parameters section.
3
Click  Load from File.
4
Browse to the model’s Application Libraries folder and double-click the file fuel_cell_cathode_parameters.txt.
Geometry 1
Now draw the model geometry. Use blocks to define the electrolyte and the porous electrode domains. Then use a work plane to draw the inlet hole at the top of the porous electrode. Facilitate geometry selection later (when setting up the physics interfaces) by enabling Resulting objects selection and renaming the geometry objects.
Membrane
1
In the Geometry toolbar, click  Block.
2
In the Settings window for Block, type Membrane in the Label text field.
3
Locate the Size and Shape section. In the Width text field, type S.
4
In the Depth text field, type S.
5
In the Height text field, type d.
6
Locate the Selections of Resulting Entities section. Select the Resulting objects selection checkbox.
By enabling Resulting objects selection here, the domain created by this rectangular block will also be available as a named domain option later on when setting up the physics.
Cathode Gas Diffusion Electrode
1
Right-click Membrane and choose Duplicate.
2
In the Settings window for Block, type Cathode Gas Diffusion Electrode in the Label text field.
3
Locate the Position section. In the z text field, type d.
Anode Gas Diffusion Electrode
1
Right-click Cathode Gas Diffusion Electrode and choose Duplicate.
2
In the Settings window for Block, type Anode Gas Diffusion Electrode in the Label text field.
3
Locate the Position section. In the z text field, type -d.
4
Click  Build Selected.
Your geometry should now look like this:
Inlet
Proceed to draw the inlet hole, placed at the top of the cathode gas diffusion electrode block.
1
In the Geometry toolbar, click  Work Plane.
2
In the Settings window for Work Plane, type Inlet in the Label text field.
3
Locate the Plane Definition section. In the z-coordinate text field, type 2*d.
4
Locate the Selections of Resulting Entities section. Select the Resulting objects selection checkbox.
Inlet (wp1) > Plane Geometry
In the Model Builder window, click Plane Geometry.
Inlet (wp1) > Circle 1 (c1)
1
In the Work Plane toolbar, click  Circle.
2
In the Settings window for Circle, locate the Size and Shape section.
3
In the Radius text field, type 1[mm].
4
In the Sector angle text field, type 90.
5
Locate the Position section. In the xw text field, type S.
6
In the yw text field, type S.
7
Locate the Rotation Angle section. In the Rotation text field, type 180.
8
Click  Build Selected.
Your work plane 2D geometry should now contain a quarter of a circle, looking as follows:
Form Union (fin)
1
In the Home toolbar, click  Build All.
2
Click the  Zoom Extents button in the Graphics toolbar.
The final 3D geometry should now look like this:
3
In the Model Builder window, under Component 1 (comp1) > Geometry 1 click Form Union (fin).
Hydrogen Fuel Cell (fc)
In the first part of the tutorial, a secondary (not concentration dependent) current distribution is modeled. Diffusion is hence disabled in the H2 and O2 gas phase mixtures. The default gas species are hydrogen and water on the anode side, and oxygen, nitrogen and water on the cathode side.
1
In the Model Builder window, under Component 1 (comp1) click Hydrogen Fuel Cell (fc).
2
In the Settings window for Hydrogen Fuel Cell, locate the H2 Gas Mixture section.
3
Find the Transport mechanisms subsection. Clear the Include gas phase diffusion checkbox.
4
Locate the O2 Gas Mixture section. Clear the Include gas phase diffusion checkbox.
A number of domain nodes, defining the different phases present in the model were added by default. The active selection of these nodes are locked, but may be controlled by adding additional domain nodes (such as Membrane etc). Start by adding these additional nodes, and make the corresponding selections on the geometry.
Membrane 1
1
In the Physics toolbar, click  Domains and choose Membrane.
2
In the Settings window for Membrane, locate the Domain Selection section.
3
From the Selection list, choose Membrane.
H2 Gas Diffusion Electrode 1
1
In the Physics toolbar, click  Domains and choose H2 Gas Diffusion Electrode.
2
In the Settings window for H2 Gas Diffusion Electrode, locate the Domain Selection section.
3
From the Selection list, choose Anode Gas Diffusion Electrode.
O2 Gas Diffusion Electrode 1
1
In the Physics toolbar, click  Domains and choose O2 Gas Diffusion Electrode.
2
In the Settings window for O2 Gas Diffusion Electrode, locate the Domain Selection section.
3
From the Selection list, choose Cathode Gas Diffusion Electrode.
Electrolyte Phase 1
The Electrolyte Phase node should now be active on all three domains. Define the conductivity in the Electrolyte Phase node by using the Fuel Cell and Electrolyzer material library, which contains conductivity data for some common electrolytes.
Add Material
1
In the Materials toolbar, click  Add Material to open the Add Material window.
2
Go to the Add Material window.
3
In the tree, select Fuel Cell and Electrolyzer > Polymer Electrolytes > Nafion®, EW 1100, Vapor Equilibrated, Protonated.
4
Right-click and choose Add to Component 1 (comp1).
5
In the Materials toolbar, click  Add Material to close the Add Material window.
Materials
Nafion®, EW 1100, Vapor Equilibrated, Protonated (mat1)
Note that the added material by default is active on all domains.
Hydrogen Fuel Cell (fc)
Electrolyte Phase 1
The polymer electrolyte conductivity depends on the temperature and the relative humidity. Specify the temperature globally in the Default Model Inputs node. The temperature defined in the Default Model Inputs node may be accessed by multiple physics nodes in the model (for instance in the Nernst and Butler-Volmer equations that will be set later in the Porous Electrode Reaction node). Specify the relative humidity for the membrane electrolyte in the Membrane node.
Global Definitions
Default Model Inputs
1
In the Model Builder window, under Global Definitions click Default Model Inputs.
2
In the Settings window for Default Model Inputs, locate the Browse Model Inputs section.
3
In the tree, select General > Temperature (K) - minput.T.
4
Find the Expression for remaining selection subsection. In the Temperature text field, type T.
Hydrogen Fuel Cell (fc)
Membrane 1
1
In the Model Builder window, under Component 1 (comp1) > Hydrogen Fuel Cell (fc) click Membrane 1.
2
In the Settings window for Membrane, locate the Electrolyte Water Activity for Material Model Input section.
3
In the aw text field, type RH.
Note that you only need to specify the water activity on the Membrane node. The water activity in the polymer of the gas diffusion electrodes is approximated to be in equilibrium with the adjacent gas phase in the pores, and is hence automatically set to equal the relative humidity in the gas phase in the Gas Diffusion Electrode nodes.
H2 Gas Phase 1
The H2 Gas Phase node should be active in domain 1 only.
Set up the composition of the H2 gas phase mixture using the Humidified mixture option.
1
In the Model Builder window, click H2 Gas Phase 1.
2
In the Settings window for H2 Gas Phase, locate the Composition section.
3
From the Mixture specification list, choose Humidified mixture.
4
In the RHhum text field, type RH.
5
In the Thum text field, type T.
O2 Gas Phase 1
The O2 Gas Phase node should be active in domain 3 only.
Similarly, set up the composition of the O2 gas phase mixture using the Humidified air option.
1
In the Model Builder window, click O2 Gas Phase 1.
2
In the Settings window for O2 Gas Phase, locate the Composition section.
3
From the Mixture specification list, choose Humidified air.
4
In the RHhum text field, type RH.
5
In the Thum text field, type T.
H2 Gas Diffusion Electrode 1
Next set up the properties of the H2 Gas Diffusion Electrode node. Note that the electrolyte volume fraction is used to calculate the effective electrolyte conductivity in the porous gas diffusion electrode.
1
In the Model Builder window, click H2 Gas Diffusion Electrode 1.
2
In the Settings window for H2 Gas Diffusion Electrode, locate the Electrode Charge Transport section.
3
In the σs text field, type sigma_s.
4
Locate the Effective Electrolyte Charge Transport section. In the εl text field, type eps_l.
H2 Gas Diffusion Electrode Reaction 1
The thermodynamics and kinetics of the hydrogen oxidation reaction are set in the child node that is added by default. Note that the reference equilibrium potential is calculated automatically when the default Built in option is used.
1
In the Model Builder window, click H2 Gas Diffusion Electrode Reaction 1.
2
In the Settings window for H2 Gas Diffusion Electrode Reaction, locate the Electrode Kinetics section.
3
From the Kinetics expression type list, choose Linearized Butler–Volmer.
4
In the i0,ref(T) text field, type i0_ref_H2.
5
Locate the Active Specific Surface Area section. In the av text field, type Av.
O2 Gas Diffusion Electrode 1
Set up the properties of the O2 Gas Diffusion Electrode node in the same way.
1
In the Model Builder window, under Component 1 (comp1) > Hydrogen Fuel Cell (fc) click O2 Gas Diffusion Electrode 1.
2
In the Settings window for O2 Gas Diffusion Electrode, locate the Electrode Charge Transport section.
3
In the σs text field, type sigma_s.
4
Locate the Effective Electrolyte Charge Transport section. In the εl text field, type eps_l.
O2 Gas Diffusion Electrode Reaction 1
The thermodynamics and kinetics of the oxygen reduction reaction are similarly set in the child node that is added by default. Note that the reference equilibrium potential is calculated automatically when the default Built in option is used.
1
In the Model Builder window, click O2 Gas Diffusion Electrode Reaction 1.
2
In the Settings window for O2 Gas Diffusion Electrode Reaction, locate the Electrode Kinetics section.
3
In the i0,ref(T) text field, type i0_ref_O2.
4
Locate the Active Specific Surface Area section. In the av text field, type Av.
Finalize the secondary current distribution model by setting up the boundary conditions for the potentials in the electronic conducting phase.
Electronic Conducting Phase 1
In the Model Builder window, under Component 1 (comp1) > Hydrogen Fuel Cell (fc) click Electronic Conducting Phase 1.
Electric Ground 1
1
In the Physics toolbar, click  Attributes and choose Electric Ground.
2
Select Boundary 3 only.
Electronic Conducting Phase 1
In the Model Builder window, click Electronic Conducting Phase 1.
Electric Potential 1
1
In the Physics toolbar, click  Attributes and choose Electric Potential.
2
Select Boundary 10 only.
3
In the Settings window for Electric Potential, locate the Electric Potential section.
4
In the ϕs,bnd text field, type E_cell.
Mesh 1
The default mesh that will be used automatically is fairly coarse, featuring only one or two mesh elements in the z direction. To improve accuracy of the results, the mesh needs to be refined. For this geometry a swept mesh can be used to get accurate control of the number of elements in the z direction. Inspect the default mesh before refining it.
1
In the Model Builder window, under Component 1 (comp1) right-click Mesh 1 and choose Build All.
2
In the Settings window for Mesh, locate the Physics-Controlled Mesh section.
3
From the Element size list, choose Fine.
Swept 1
1
In the Mesh toolbar, click  Swept.
2
In the Settings window for Swept, click to expand the Source Faces section.
3
Select Boundaries 10 and 14 only.
Distribution 1
1
Right-click Swept 1 and choose Distribution.
2
In the Settings window for Distribution, locate the Domain Selection section.
3
From the Selection list, choose Cathode Gas Diffusion Electrode.
4
Locate the Distribution section. From the Distribution type list, choose Predefined.
We will create a mesh with thinner elements in the cathode gas diffusion electrode domain toward the boundary facing the membrane domain.
5
In the Element ratio text field, type 5.
6
Select the Reverse direction checkbox.
7
Click  Build All.
The mesh should now look as follows:
Swept 1
To improve the resolution along the current collector-inlet hole edge, and the interior of the electrode, also add some Size nodes.
Size 1
1
In the Model Builder window, right-click Swept 1 and choose Size.
2
In the Settings window for Size, locate the Geometric Entity Selection section.
3
From the Geometric entity level list, choose Edge.
4
Select Edge 19 only.
5
Locate the Element Size section. From the Predefined list, choose Extra fine.
Size 2
1
Right-click Swept 1 and choose Size.
2
In the Settings window for Size, locate the Geometric Entity Selection section.
3
From the Geometric entity level list, choose Edge.
4
Select Edges 10 and 11 only.
5
Locate the Element Size section. From the Predefined list, choose Finer.
Size 3
1
Right-click Swept 1 and choose Size.
2
In the Settings window for Size, locate the Geometric Entity Selection section.
3
From the Geometric entity level list, choose Point.
4
Select Point 4 only.
5
Locate the Element Size section. From the Predefined list, choose Extra fine.
6
Click  Build All.
Study 1
The settings for the secondary current distribution model are now complete.
1
In the Study toolbar, click  Compute.
Results
Electrode Potential with Respect to Ground (fc)
Inspect the default plots.
1
In the Electrode Potential with Respect to Ground (fc) toolbar, click  Plot.
The electrode potential plot should look as follows:
Definitions
We will now compute and plot a polarization curve, that is, solve for a range of cell potentials, and plot these versus the average cell current density. First introduce a boundary probe for the average cell current density at the anode gas diffusion electrode boundary.
Boundary Probe 1 (bnd1)
1
In the Definitions toolbar, click  Probes and choose Boundary Probe.
2
In the Settings window for Boundary Probe, type I_avg in the Variable name text field.
3
Select Boundary 3 only.
4
Click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > fc.nIs - Normal electrode current density - A/m².
Change the sign of the expression.
5
Locate the Expression section. In the Expression text field, type -fc.nIs.
6
In the Table and plot unit field, type A/cm^2.
7
Select the Description checkbox. In the associated text field, type Average cell current density.
Study 1
Step 2: Stationary
Set up an auxiliary sweep to solve for a range of cell potential values and compute the model again.
1
In the Model Builder window, under Study 1 click Step 2: Stationary.
2
In the Settings window for Stationary, click to expand the Study Extensions section.
3
Select the Auxiliary sweep checkbox.
4
Click  Add.
5
In the table, enter the following settings:
Parameter name
Parameter value list
Parameter unit
E_cell (Cell voltage)
range(1,-0.1,0.5)
V
6
In the Study toolbar, click  Compute.
Results
Polarization Curve
A probe plot for the average current density was now created by default. Modify it as follows:
1
In the Model Builder window, under Results click Probe Plot Group 3.
2
In the Settings window for 1D Plot Group, type Polarization Curve in the Label text field.
3
Locate the Plot Settings section. Select the Flip the x- and y-axes checkbox.
4
Select the x-axis label checkbox. In the associated text field, type Average cell current density (A/cm<sup>2</sup>).
5
Select the y-axis label checkbox. In the associated text field, type Cell voltage (V).
6
In the Polarization Curve toolbar, click  Plot.
Hydrogen Fuel Cell (fc)
Now we start the second part of the tutorial to incorporate mass and momentum transport. Include mass transport using Maxwell-Stefan diffusion and momentum transport using Darcy’s Law in the O2 gas phase mixture.
1
In the Model Builder window, under Component 1 (comp1) click Hydrogen Fuel Cell (fc).
2
In the Settings window for Hydrogen Fuel Cell, locate the O2 Gas Mixture section.
3
Find the Transport mechanisms subsection. Select the Include gas phase diffusion checkbox.
4
Select the Use Darcy’s Law for momentum transport checkbox.
O2 Gas Phase 1
Inspect the settings of the O2 Gas Phase node. Note that since you are now including diffusion in the model, the composition values you specified earlier are no longer visible. (Settings for Diffusion are now present instead.) The initial and inlet conditions (composition and pressure) of the O2 gas phase mixture are now specified using child nodes. The gas composition is specified using the Humidified air option.
Set up the initial and inlet conditions (composition and pressure) of the O2 gas phase mixture.
Initial Values 1
1
In the Model Builder window, click Initial Values 1.
2
In the Settings window for Initial Values, locate the Initial Composition section.
3
From the Mixture specification list, choose Humidified air.
4
In the RHhum text field, type RH.
5
In the Thum text field, type T.
O2 Gas Phase 1
The same parameter values that were used for the initial values are used to specify the inlet mole fractions.
1
In the Model Builder window, click O2 Gas Phase 1.
O2 Inlet 1
1
In the Physics toolbar, click  Attributes and choose O2 Inlet.
2
In the Settings window for O2 Inlet, locate the Boundary Selection section.
3
From the Selection list, choose Inlet.
4
Locate the Inlet Flow Type section. From the Inlet flow type list, choose Mixture composition constraint.
O2 Gas Diffusion Electrode 1
Finally, set up the gas transport properties in the O2 Gas Diffusion Electrode node.
1
In the Model Builder window, under Component 1 (comp1) > Hydrogen Fuel Cell (fc) click O2 Gas Diffusion Electrode 1.
2
In the Settings window for O2 Gas Diffusion Electrode, locate the Gas Transport section.
3
In the εg text field, type eps_gas.
4
In the κg text field, type perm.
Note that the effect of varying concentration is automatically taken into account in the built-in thermodynamic and kinetic expressions of the oxygen reduction reaction in the O2 Gas Diffusion Electrode Reaction child node, and appropriate mass sources resulting from the electrochemical reaction in the O2 gas phases mixture are automatically defined. Hence, no additional settings are required in this node.
Results
Before proceeding to solve the model with transport effects, duplicate the probe table and rename the copy appropriately.
Unlimited O2 gas phase transport
1
In the Model Builder window, expand the Results > Tables node.
2
Right-click Probe Table 1 and choose Duplicate.
3
In the Settings window for Table, type Unlimited O2 gas phase transport in the Label text field.
Study 1
In order to generate new default plots related to the introduced mass transport, remove the old study sequence before recomputing.
Solver Configurations
1
In the Model Builder window, under Study 1 right-click Solver Configurations and choose Delete Configurations.
2
In the Home toolbar, click  Compute.
Results
Polarization Curve
Modify the Polarization Curve as follows to compare the concentration independent and concentration dependent solutions.
Probe Table Graph: Limited O2 gas phase transport
1
In the Model Builder window, expand the Polarization Curve node, then click Probe Table Graph 1.
2
In the Settings window for Table Graph, type Probe Table Graph: Limited O2 gas phase transport in the Label text field.
3
Click to expand the Legends section. From the Legends list, choose Manual.
4
In the table, enter the following settings:
Legends
Limited O2 gas phase transport
Probe Table Graph: Unlimited O2 gas phase transport
1
Right-click Probe Table Graph: Limited O2 gas phase transport and choose Duplicate.
2
In the Settings window for Table Graph, type Probe Table Graph: Unlimited O2 gas phase transport in the Label text field.
3
Locate the Data section. From the Table list, choose Unlimited O2 gas phase transport.
4
Locate the Legends section. In the table, enter the following settings:
Legends
Unlimited O2 gas phase transport
Polarization Curve
1
In the Model Builder window, click Polarization Curve.
2
In the Polarization Curve toolbar, click  Plot.
Mole Fraction, O2, Streamline (fc)
The mole fractions of the different species are plotted by default at the cell potential of 0.5 V. Modify the O2 plots as follows to plot at the cell potential of 0.7 V.
1
In the Model Builder window, click Mole Fraction, O2, Streamline (fc).
2
In the Settings window for 3D Plot Group, locate the Data section.
3
From the Parameter value (E_cell (V)) list, choose 0.7.
Note the direction of the arrows. Oxygen flows from the inlet hole into the porous cathode to react to form water.
Mole Fraction, O2, Surface (fc)
1
In the Model Builder window, click Mole Fraction, O2, Surface (fc).
2
In the Settings window for 3D Plot Group, locate the Data section.
3
From the Parameter value (E_cell (V)) list, choose 0.7.
The oxygen mole fraction gets low far away from the inlet hole.
Pressure (fc)
The Darcy pressure with velocity streamlines is also plotted by default at the cell potential of 0.5 V. Modify as follows to plot at the cell potential of 0.7 V.
1
In the Model Builder window, click Pressure (fc).
2
In the Settings window for 3D Plot Group, locate the Data section.
3
From the Parameter value (E_cell (V)) list, choose 0.7.
The direction of the net velocity is toward the inlet hole, that is, opposite to the oxygen flux. This a result of the production of two water molecules per consumed oxygen molecule in the cathode.
Also create some additional plots for the activation overpotential and local volumetric current density in the cathode gas diffusion electrode, and the current density at the anode gas diffusion electrode boundary.
Overpotential in Cathode
1
In the Results toolbar, click  3D Plot Group.
2
In the Settings window for 3D Plot Group, type Overpotential in Cathode in the Label text field.
3
Locate the Data section. From the Parameter value (E_cell (V)) list, choose 0.7.
Surface 1
1
Right-click Overpotential in Cathode and choose Surface.
2
In the Settings window for Surface, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Electrode kinetics > fc.eta_o2gder1 - Overpotential - V.
3
In the Overpotential in Cathode toolbar, click  Plot.
Generally the highest overpotentials (in magnitude) are found in the region facing the Membrane domain. Since the overpotential is the driving force for the electrochemical reactions, this is the region were we can expect higher reaction rates.
Local Volumetric Current Density in Cathode
1
In the Results toolbar, click  3D Plot Group.
2
In the Settings window for 3D Plot Group, type Local Volumetric Current Density in Cathode in the Label text field.
3
Locate the Data section. From the Parameter value (E_cell (V)) list, choose 0.7.
Surface 1
1
Right-click Local Volumetric Current Density in Cathode and choose Surface.
2
In the Settings window for Surface, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Electrode kinetics > fc.iv_o2gder1 - Local current source - A/m³.
3
In the Local Volumetric Current Density in Cathode toolbar, click  Plot.
As for the overpotentials, the highest current density magnitudes are found close to the Membrane domain.
Current Density at Anode Boundary
1
In the Results toolbar, click  3D Plot Group.
2
In the Settings window for 3D Plot Group, type Current Density at Anode Boundary in the Label text field.
3
Locate the Data section. From the Parameter value (E_cell (V)) list, choose 0.7.
4
Click to expand the Selection section. From the Geometric entity level list, choose Boundary.
5
Select Boundary 6 only.
Surface 1
1
Right-click Current Density at Anode Boundary and choose Surface.
2
In the Settings window for Surface, click Replace Expression in the upper-right corner of the Expression section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > fc.nIl - Normal electrolyte current density - A/m².
3
Locate the Expression section. In the Expression text field, type abs(fc.nIl).
The abs() is an operator which will return the absolute (positive) value of the argument.
4
Select the Description checkbox. In the associated text field, type Current density.
5
In the Current Density at Anode Boundary toolbar, click  Plot.
The region of highest current densities is located below the quarter circular edge of the inlet hole. In this area the combined effects of the ohmic and mass transfer losses in the gas diffusion electrode are at a minimum.
Power Losses
Finally, we will plot the different sources of power losses in the cell.
1
In the Home toolbar, click  Add Plot Group and choose 1D Plot Group.
2
In the Settings window for 1D Plot Group, type Power Losses in the Label text field.
Global 1
1
Right-click Power Losses and choose Global.
2
In the Settings window for Global, click Add Expression in the upper-right corner of the y-Axis Data section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Power losses > Feature-node integrated > fc.o2gde1.o2gder1.P_loss_act - Kinetic activation power loss - W.
3
Click Add Expression in the upper-right corner of the y-Axis Data section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Power losses > Cell integrated > fc.P_loss_l - Electrolyte transport power loss - W.
4
Click Add Expression in the upper-right corner of the y-Axis Data section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Power losses > Cell integrated > fc.P_loss_s - Electron conduction power loss - W.
5
Click Add Expression in the upper-right corner of the y-Axis Data section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Power losses > Cell integrated > fc.P_loss_gas - Gas transport power loss - W.
6
Click Add Expression in the upper-right corner of the y-Axis Data section. From the menu, choose Component 1 (comp1) > Hydrogen Fuel Cell > Power losses > Feature-node integrated > fc.h2gde1.h2gder1.P_loss_act - Kinetic activation power loss - W.
7
Locate the y-Axis Data section. In the table, enter the following settings:
Expression
Unit
Description
fc.o2gde1.o2gder1.P_loss_act
W
Kinetic activation, O2
fc.P_loss_l
W
Electrolyte transport
fc.P_loss_s
W
Electron conduction
fc.P_loss_gas
W
Gas transport
fc.h2gde1.h2gder1.P_loss_act
W
Kinetic activation, H2
8
Locate the x-Axis Data section. From the Parameter list, choose Expression.
9
In the Expression text field, type I_avg.
10
In the Unit field, type A/cm^2.
Power Losses
1
In the Model Builder window, click Power Losses.
2
In the Settings window for 1D Plot Group, locate the Plot Settings section.
3
Select the x-axis label checkbox. In the associated text field, type Average cell current density (A/cm<sup>2</sup>).
4
Select the y-axis label checkbox. In the associated text field, type Power Loss (W).
5
Locate the Legend section. From the Position list, choose Upper left.
6
In the Power Losses toolbar, click  Plot.
Cell-Averaged Overpotentials
By dividing the integrated power losses of the cell by the cell current we can also compute corresponding overpotentials for the whole cell.
1
Right-click Power Losses and choose Duplicate.
2
In the Settings window for 1D Plot Group, type Cell-Averaged Overpotentials in the Label text field.
Global 1
1
In the Model Builder window, expand the Cell-Averaged Overpotentials node, then click Global 1.
2
In the Settings window for Global, locate the y-Axis Data section.
3
In the table, enter the following settings:
Expression
Unit
Description
fc.o2gde1.o2gder1.P_loss_act/(I_avg*A_cell)
V
Kinetic activation, O2
fc.P_loss_l/(I_avg*A_cell)
V
Electrolyte transport
fc.P_loss_s/(I_avg*A_cell)
V
Electron conduction
fc.P_loss_gas/(I_avg*A_cell)
V
Gas transport
fc.h2gde1.h2gder1.P_loss_act/(I_avg*A_cell)
V
Kinetic activation, H2
Cell-Averaged Overpotentials
1
In the Model Builder window, click Cell-Averaged Overpotentials.
2
In the Settings window for 1D Plot Group, locate the Plot Settings section.
3
In the y-axis label text field, type Cell-averaged overpotential (V).
4
In the Cell-Averaged Overpotentials toolbar, click  Plot.