Theory of Multicomponent Phenomena in Cation-Exchange Membranes: Part III. Transport in Vanadium Redox-Flow-Battery Separators

Transport through vanadium redox- ﬂ ow-battery membranes strongly in ﬂ uences cell performance. In this work, we use a multicomponent concentrated-solution model of transport and thermodynamics in phase-separated cation-exchange membranes, the most common separator type, to develop structure-performance relationships. The model incorporates species partitioning into the membrane, thermodynamic nonidealities, and Stefan-Maxwell-Onsager frictions between species. Molecular-thermodynamics and -transport theories parameterize the model. We validate the calculations against measured Coulombic and voltage ef ﬁ ciencies of a vanadium ﬂ ow battery as a function of current density. Our model shows that species transport is the result of collective interactions between all species present in the system. The magnitude of coupling suggests that predictions made using dilute- solution theory for transport in these systems will be misleading in many situations. As a demonstration of the capabilities of the model, we predict cell performance, incorporating these interactions, as a function of electrolyte concentration and composition and membrane equivalent weight and backbone modulus. We ﬁ nd that electrolytes with high sulfuric acid concentrations provide the greatest cell performance (quanti ﬁ ed by maximizing power density at a target energy ef ﬁ ciency). In the case of membrane properties, low equivalent-weight polymers perform better; at high equivalent weights, a low membrane modulus is preferred.

Membrane separator properties critically impact the performance of vanadium redox-flow batteries (VRFB). [1][2][3][4][5][6][7][8][9] These separators, which are typically polymer membranes, facilitate ionic current between the positive and negative electrodes while limiting shorting and self-discharge due to crossover of vanadium active species. 3,5 To optimize VRFB performance, there is an optimal design window between conductivity and crossover, which are at odds, as they necessitate different membrane morphologies and intrinsic properties. 1,3,5 Thus, determining structure-performance relationships for membranes is key to successful deployment of VRFBs. As a result, the electrochemistry community has researched these transport properties using measurements and models for a variety of separator materials and operating conditions, as numerous reviews and articles outline. [1][2][3][4][5][6][7][8][9][10][11][12][13][14][15] Despite these efforts, transport in VRFB separators is still a poorly understood process due to the numerous species and modalities involved, as displayed in Fig. 1. 12,13,[15][16][17][18][19] Such transport involves both vanadium partitioning into the separator from the external electrolyte solution, 10,17,20,21 followed by ion and solvent transport across the separator driven by concentration gradients (i.e. diffusion) and the electric field (i.e. migration). 1,10,13,16,19,[22][23][24][25][26] The impact of different driving forces on species flux is particularly complicated due to the high concentration of multiple ionic species present. 16,27 High concentrations create strongly nonideal thermodynamics (e.g. large excess chemical potentials) and frictional interactions between species couple transport (e.g. concentration gradients of species i cause transport of species j). 16,23,[27][28][29] As a result of these complex conditions, transport and partitioning coefficients are functions of both the composition of electrolyte in the electrodes and membrane chemistry. 6,9,17,20,21,29,30 Mathematical models for transport in membranes deconvolute these various effects. 1,4,10,12,19,[25][26][27]29,31 Previously developed models have been invaluable for understanding cell performance across multiple cycles and how electrolyte transports through the cell. 1,4,10,12,19,[25][26][27]29,31 However, many of these efforts are restricted to dilute-solution approximations that do not consider coupled transport and/or account for concentration-dependent transport properties. The validity of these assumptions under flow-battery operating conditions is not known. Unfortunately, approaches that use concentrated-solution theory 27,29,31 are often intractable because they require numerous transport parameters, of which the few that are measured have large associated experimental error. 31 To overcome these challenges, Part I and II of this series developed and validated molecular models for coupled, multi-ion thermodynamic and transport properties in phase-separated cation-exchange membranes. 32,33 This approach completely specifies the concentrated-solution transport properties as a function of concentration and membrane chemistry. Using these models, this paper explores the role of transport coupling and links molecular-scale behavior to macroscale performance of VRFB membranes, providing design criteria and guidelines for both cell developers and membrane chemists.
The outline of this paper is as follows. In the theory section, we summarize the relevant thermodynamic and transport equations. We then examine the extent and nature of transport coupling and the impact on membrane properties. A 1-D, quasi-steady-state model of a VRFB membrane uses the transport parameters to show how vanadium crossover is detrimental. Finally, we examine how the electrolyte composition and the structure and chemistry of the membrane impact the cell performance. electrode electrolytes does not change significantly during a cycle. This work compliments previous VRFB modeling efforts that have focused on overall cell operation (see Fig. 1 top). 10,12,13,25,26,[34][35][36][37][38][39][40] This section first summarizes the VRFB system, including relevant electrochemical and chemical reactions. We then show how ohmic overpotentials and crossover of active species impact cell operation. The final sections formulate the transport and thermodynamic equations governing species flux and chemical potential drop across the membrane.

Vanadium Flow-Battery System
This study examines an all-vanadium redox-flow battery, pictured in Fig. 1. The system consists of electrolyte-filled porous electrodes that contain aqueous vanadium sulfate salts with sulfuricacid supporting electrolyte (water, H 2 O, protons, H + , and sulfate, SO 4 2− , are abbreviated 0, H, and SO 4 , respectively). The system is at a temperature of 295 K. The negative electrode, b, contains + V 2 and + V 3 (denoted V(II) and V(III), respectively, based on vanadium oxidation state) and the positive electrode,  ,, contains + VO 2 and + VO 2 (denoted V(IV) and V(V) respectively). Vanadium in any oxidation state is denoted ( ) V x . As the cell discharges, the reactions are ( ) ( )  V II V III at the negative electrode and ( ) ( )  V V V IV at the positive electrode. Sulfate species associate with vanadium ions and with protons to form various ion-paired products. 41 Table I details the electrochemical reactions and their half-cell potentials, the chemical reactions that occur when vanadium crosses over the membrane separating the electrodes (phase d), and the ion association equilibrium reactions.

Cell Performance
Cell performance is typically characterized by the power density on discharge, Y d , and the round-trip energy efficiency, e e . 42 Y d is the product of the discharge current density i d and the cell potential on The energy efficiency is the ratio of the integrated cell power density on discharge to charge and is typically treated as the product of the voltage efficiency, e v , and the Coulombic efficiency, e q , 6 where t is time, the subscript c or d denotes charge or discharge, and the cell voltage is approximated as constant during charge or discharge. e v characterizes efficiency losses due to cell overpotentials (e.g. kinetic and ohmic) and decreases with increasing current density. 3 e q quantifies efficiency losses due to vanadium passing through the membrane between electrode solutions, mixing, and reacting chemically rather than electrochemically. e q tends to increase with rising current density because the total charge extracted from the cell generally grows compared to the amount of charge lost to chemical reactions. 3 A typical operating goal is maximizing power density at a target energy efficiency. 42 The next two sections present the governing equations for cell potential that are required to calculate power density and voltage efficiency, and the state-of-charge loss during cycling that dictates Coulombic efficiency.
Cell Potential.-The cell potential is the difference in the electrochemical potential of the electrons me in the metal at the negative electrode a¢ and the positive electrode a; these are related to the electrochemical potential of species i, m f i , participating in the oxidation/reduction reactions in electrode f as outlined in Table I where F is Faraday's constant and h f is ionic potential drop between the metal and the electrolyte in electrode f due to mass-transport resistances, surface overpotentials, or other non-membrane cell resistances.
To put Eq. 1 in terms of the electrochemical potentials of species in the membrane, chemical equilibrium requires that the electrochemical potential of charged species (and the chemical potential of neutral species) must be equal at the interface of electrode phase f and membrane phase d 44 where m i is a function of chemical variables-temperature, pressure, composition-and, for charged species, the ionic potential. It is convenient to group electrostatic dependences into a single variable with the remaining variables depending only on composition, temperature, and pressure. We reference the electrochemical potential of species i to the electrochemical potential of the membrane M to define which is an electroneutral pairing of species and is independent of electrostatic potential. 43,45 To track the electrical state through the membrane, we set the ionic potential F equal to the electrochemical potential of the membrane, m F = d z F M M . This assignment is equivalent to measuring the ionic potential using a hypothetical reference electrode selective to the ionic groups on the membrane at the same conditions. 43 Using Eq. 4 and m d i,M , the cell potential is then groups non-membrane overpotentials in the cell. U is the equilibrium cell potential [( where d¢ and d are the membrane phase at the negative and positive electrode interfaces, respectively. DF = F -F d d  ¢ is the potential drop across the membrane due to ohmic and diffusion overpotentials. During discharge (charge), h and DF become more negative (positive) with increasing current density, causing the cell voltage efficiency to decrease. 3 We set U to the measured open-circuit potential (=1.42 V) 3 by neglecting the diffusion potential, which tends to be small. 43 η is modeled as an ohmic potential drop with resistance W R where h is negative for discharge and positive for charge.
State-of-charge loss.-Vanadium crossover causes the cell to self-discharge. The cell state of charge, SOC, characterizes the extent the cell is in the charged state; it is the ratio of the sum of the moles of vanadium in the charged state (namely, ( ) V II in the negative electrode and ( ) V V in the positive electrode) to the total moles of vanadium in the system where f n i is the moles of i in electrode f. The rate at which the state of charge changes during discharging or charging is a function of the cell current density i and the crossover current density i x where A is the active area of the membrane and n t is the total moles of vanadium in the system (i.e. ( ) ( ) ( ) ( ) . i x characterizes how vanadium crossover decreases the state of charge of the cell and is a function of the flux of vanadium species i through the membrane, N i , according to [1][2][3] [ ] The membrane/negative electrode interface is set at = z 0 and the membrane/positive electrode interface is at = z l M so that N i is positive for species that move towards the positive electrode (e.g. ( ) V II and ( ) V III ) and negative for species that move toward the negative electrode (e.g. ( ) V IV and ( ) V V ). With this convention for the direction of flux, the ionic current density i through the membrane is positive on discharge and is also the electronic current leaving the positive electrode normalized by A (i.e. positive on discharge). 1 A subsequent section shows that N i is constant across the membrane.
The first term on the right side of Eq. 10 is the consumption/ production of species due to electrochemical charge/discharge, the second term is due to the loss of ( ) V II and ( ) V V from crossover to the opposing electrode and from reactions with vanadium species that transport from the opposing electrode (Table I specifies the stoichiometry incorporated in Eq. 11). Based on the low concentration of absorbed vanadium in the membrane, the chemical reactions in Table I occur in the electrode solutions and not in the membrane. The supporting electrolyte provides an excess of protons, and the high concentration of vanadium in the electrode solution rapidly reacts with any vanadium that crosses the membrane into the electrode. 22 The factor of 1/2 on the left side of Eq. 10 results from considering the combined state of charge of the positive and negative electrodes. i x is a function of current density and cell state of charge (see Eq. 11). 13,36 As a result of vanadium crossover, the amount of charge extracted from the cell during discharge, ò i Adt Because i x is a function cell state of charge, the condition for Eq. 12 of constant i x is achieved by a small change to SOC over charge/discharge.

Membrane Phenomena
The following two sections discuss how the electrolyte partitions into and how it transports across the membrane. These expressions provide for calculation of e q , e v , and DF. Although the following expressions are general for any flow-battery membrane, the specific microscale theories we use to calculate thermodynamic and transport coefficient are specific to phase-separated polymer cation-exchange membrane that Parts I and II present. 32,33 The most common class of these materials are perfluorinated sulfonic-acid (PFSA) ionomers, such as Nafion. 46 These materials consist of nanoscale, interconnected hydrophilic domains filled with aqueous electrolyte solution and side-chains terminated with anionic moieties that are tethered to the PFSA polymer. 46 Surrounding these domains are hydrophobic polytetrafluoroethylene backbone that provides structural support. 46 Membrane thermodynamics.-Requisite for species to transport across the membrane, they must sorb from the electrode electrolyte into the membrane as Eq. 4 specifies. Chemical thermodynamics quantifies how these potentials are related to measurable quantities. The chemical potential of a neutral species is 44 [ ] where R is the universal gas constant, T is temperature, and x i is the mole fraction of i ( / =å n n j i j where n j is the moles of species j). The mole fraction of charged polymer groups in the membrane is x M and the charge for most cation-exchange membranes is =z 1 M . The reference chemical potential of i, m q i , is for a hypothetical ideal solution at unit mole fraction x i and is a function of temperature and pressure. 44 The second term on the right side of Eq. 13 is the ideal solution contribution. 44 In this paper, the ideal-solution contribution is specified a for fully-dissociated electrolyte (i.e. the Experimental Construct discussed in Part I). 32 The excess chemical potential accounts for ion/ ion, ion/membrane, solvent/membrane, and solvent/ion interactions including ion-pair formation (see Table I), electrostatic interactions, ion solvation, steric confinement, and membrane swelling. Part I details calculation of the excess chemical potential m .  For a charged species, the electrochemical potential has the same form as Eq. 13 but further depends on the electrical state of the phase and charge number z . i 43 However, expressing the electrochemical chemical potential using only neutral pairings of two charged of species (i.e. m m -z z i i j j ) eliminates the dependence on electrical potential. 43 Substitution of these neutral pairings of Eq. 13 into Eq. 4 relates species mole fraction in the membrane and electrode electrolyte 47 For neutral species (e.g. water) = z 0 i and, therefore, the jth species does not need to be specified in Eq. 14. Electroneutrality inside the membrane adds a constraint that fully specifies ion partitioning, 44 [ ] If oppositely charged ions are used for i and j (such as the anion and cation of a salt) and G ij is unity, the oppositely charged pairs i and j partition into the membrane according to ideal Donnan equilibrium. 47 If G > 1 ij , i and j favorably partition into the membrane whereas if partitioning is unfavorable, G < 1.
ij 47 For context, Part I shows that for a Nafion cation-exchange membrane in an aqueous hydrobromicacid solution at concentration of 5 mol kg −1 -solvent, the term on the right side of Eq. 14 is ∼0.75 for HBr and 0.95 for water. 32

Membrane Transport
Upon sorption into the membrane, species transport according to their molar flux N i across the membrane (z-direction in the 1-D model). The molar flux obeys species conservation at steady-state 48 Although the concentration of species in the electrodes changes during cell operation, transport in the membrane is assumed to be at steady-state because the concentration and pressure gradients in the membrane develop much faster (i.e. pseudo steady-state approximation). 48 Absent temperature and hydraulic-pressure gradients, which are typically small, transport is driven by electrochemical potential gradients. Friction between species introduces drag that reduces total flux. Stefan-Maxwell-Onsager theory describes the balance between driving force and drag on i not equal to M as 43,49,50 and for the membrane as where c i is the molar concentration of species i and c T is the total molar concentration in the membrane. (For convenience in this section, we do not superscript variables for quantities in the membrane phase d). As discussed in the preceding section, m i is a function of temperature, pressure, composition, and ionic potential F. D ij is the binary diffusion coefficient between species i and j. D jM is related to the friction coefficient between species i and the membrane, K iM , that is described in Part II, according to The membrane is stationary (i.e. is the reference velocity) and the membrane concentration c M is set as the molar concentration of charged groups on the polymer. Onsager reciprocal relations dictate that the diffusion coefficients are symmetric so that D D = ij ji . Consequently, for a VRFB with eight species (water, sulfate ions, protons, the membrane, and four vanadium species), there are 28 transport coefficients. In this paper, Eqs. 17 and 18 consider ionic species as fully dissociated (i.e. uses the Experimental Construct discussed in Part II). 33 Transport coefficients D ij include the effect of ion-pair formation to ensure that this model is consistent with the various ion-paired species that exist. Part II outlines calculation of D ij and K iM as a function of membrane water content and ion concentration. 33 Equations 17 and 18 are rigorous but inconvenient because they involve gradients in electrochemical potential that are not readily characterized and frame the driving force in terms of species fluxes, whereas experiments measure fluxes under applied forces. Appendix A shows that for constant total molar concentration, c , T eliminating the electrostatic dependence of the driving forces (i.e. making the substituting m m m M as shown in Eq. 5), 45 and expanding chemical potential gradients in terms of composition variables gives the flux of species i as the sum of migration and diffusion terms: where n is a reference species with charge number different from M, t i is the transference number of species i, and D ij are multicomponent diffusion coefficient between species i and j based on concentration driving forces. We assign protons as reference species n. Although t 0 and z 0 are zero for water, the ratio / t z 0 0 is definite and equal to the electroosmotic coefficient x. ). This formalism generalizes the treatment of concentration gradients in binary electrolytes 43 to an arbitrary number of species; the concentration gradient of j, M reduces to the concentration gradient of the salt for a binary electrolyte with | | -z cations and | | + z anions. Appendix A shows how D ij is related to the binary diffusion coefficients D ij 's and chemical potential gradients. Specifically, D ij 's are components of the -N 2 by -N 2 diffusion matrix D that is the product of the matrix inverse of B, an -N 2 by -N 2 matrix containing transport coefficients, and c, an -N 2 by -N 2 matrix containing multicomponent thermodynamic factors reduces to a well-established form for the salt diffusion coefficients in the case of a binary electrolyte. 43 The -N 2 vector t of charge number-normalized transference ij Since Eq. 19 is applicable to ¹ i n, we specify the flux of n using the current density and fluxes of species ). For negligible mass-transport resistance between the electrode electrolyte and membrane, chemical equilibria (i.e. Eqs. 14 and 15) specifies c i,M at the interface of the membrane with the two electrodes. Solving Eqs. 16 and 19 using these concentrations as boundary conditions gives the concentration profiles and fluxes of species across the membrane for a set current density.
With specified fluxes, Eq. 18 calculates the ionic potential gradient because it is proportional to the electrochemical potential of the membrane, Integration of Eq. 27 across the membrane shows that DF is a function of applied current density (ohmic overpotential) and flux of species at the open-circuit cell potential (diffusion overpotential). In the absence of concentration gradients, substitution of Eq. 19 into Eq. 27 identifies the membrane conductivity, k, as As k increases for a fixed current density and concentration gradient, the potential drop across the membrane decreases.

Parameters
Parts I and II discuss calculation of m i ex and D , ij respectively, as functions of concentration and water content as well as membrane equivalent weight and the modulus of the hydrophobic matrix of a dry membrane; 32,33 the required system parameters for these calculations are discussed therein. The properties specific to VRFBs are the infinite dilution diffusion coefficients in water and viscous volume of vanadium ions and vanadium-bisulfate ion pairs in solution, D ¥ i,0 andṼ , i respectively, the vanadium-sulfate binary diffusion coefficient D i,SO 4 , the vanadium-sulfate and vanadiummembrane interaction parameters b i,SO 4 and b i,M , ( ) V III -and ( ) V IV -bisulfate association constants -K i HSO 4 , and solvent/vanadium binding constant k .
i Unfortunately, there are relatively few measurements of vanadium thermodynamic and transport properties in all its oxidation states at well-defined conditions. As described in Supplemental Material (SM), reported experimental values are used when available, and otherwise, parameters of cations of similar charge number approximate those of the vanadium ions. For ( ) V IV and ( ) V V we fit b i,M to measured vanadium, sulfate, and water uptake in the membrane as a function of sulfuric acid and vanadium sulfate concentration detailed in Part I. 32 Similarly, for ( ) 4 and k i to isotherm measurements. These parameters are plausible values for this system; however, they are estimates. Consequently, there is a strong need for fundamental thermodynamic studies of vanadium-ion properties in their different oxidation states (i.e. activity and osmotic coefficients as a function of concentration).
Throughout this paper, we will consider a reference VRFB system containing electrodes at a fixed composition and a membrane with fixed properties, which are provided in Table II. Calculations are made for this system unless stated. The electrolyte concentrations match those in experimental cells as described in Refs. 1, 3, 51 unless otherwise stated. In those studies, an initial solution of 1.5 mol dm −3 vanadium IV sulfate in 2.6 mol dm −3 sulfuric acid was twice charged with the positive electrolyte replaced after each charge. As outlined in SM, after this charging process and at a state of charge of 50%, the negative electrolyte has 0.75 mol dm −3 of both V(II) and V(III) and a total sulfate concentration of 4.66 mol dm −3 . The positive electrolyte has the same concentrations of V(IV) and V (V) and a total sulfate concentration of 4.11 mol dm −3 . The electrolyte solution mass density was 1.19 g cm −3 . 51 The separator in the reference system is a Nafion 212 membrane with a drymembrane thickness, l M 0 , of 51 μm. The membrane equivalent weight, EW (mass-polymer per mole-sulfonate group) is 1100.
Other membrane properties considered are hydrophobic matrix modulus of the dry membrane, E , b 0 Archie's tortuosity scaling parameter, c, hydrophilic domain geometric transport factor, G, and spacing between hydrophilic domain in dry membrane, d 0 . EW sets the intrinsic concentration of ions in the membrane absent co-ions, E b 0 limits the extent the membrane can swell with water from the surrounding environment, G affects friction between the membrane and absorbed water and aqueous ions, and d 0 and c dictate the hydrophilic domain size and network tortuosity, respectively, at a given water volume fraction. Unless otherwise stated, membrane properties are those of the Nafion PFSA separator and are detailed in Parts I and II. 32,33 Numerical Implementation Figure S1 (available online at stacks.iop.org/JES/167/013549/ mmedia) outlines calculation of membrane properties, fluxes, and performance metrics. For a specified electrolyte composition in the negative and positive electrodes, Eqs. 14 and 15 specify x i at d¢ and d. Part I outlines calculation of m .   To solve Eq. 19 for each species, we extend to electrolyte systems a method outlined by Krishna et al. 52,53 for uncharged systems. Using this approach, Appendix B shows that for a specified current density, the species fluxes are specified by iterating over a set of six transcendental algebraic equations and neglecting changes across the membrane of: the thermodynamic factor c ij because the solution composition from the most abundant species (sulfate, the membrane, water, and protons) changes relatively little across the membrane; c T because the molar density of the membrane change only slightly with composition; D ij because the dominant factors its depends upon-ionic strength and water content inside the membrane-are relatively constant; and / ¶ ¶ x z ln M because the membrane concentration changes little between the electrodes. We evaluate D ij , c ij , and λ using the mean of the composition of the membrane at the interface with the two electrodes¯d , calculating the fraction of vanadium and protons that form ion pairs with sulfate at this composition. SM shows that the error introduced by these assumption is small. The derivative of x ln M is approximated as In Results and Discussion, we show this assumption is also reasonable. Upon specifying species fluxes at a given charge and discharge current density, Eq. 11 gives i x and, using Eq. 12, e . q Integration of Eq. 27 gives F. F is a linear function of z for the assumption of constant D .
ij Upon calculating DF, Eqs. 6, 2, and 1 give V , e v , e e , and Y . d To find the maximum discharge cell power density Y d at e = 80% e (our adopted design criteria), we iteratively solve these equations for varying i until finding the maximum. 54 As previously stated, the electrode composition and N i are assumed constant during charge or discharge. This condition is achieved for incremental changes in SOC during cycling.
We performed sample-based sensitivity analysis for the impact of membrane properties EW , E b 0 , c, G, and d 0 on transport properties and cell performance. SM details the sampling technique and range for these properties. The sensitivity analysis uses a cell current density of 200 mA cm −2 and the electrolyte composition is that of the reference conditions (see Table II).

Results and Discussion
Representative VRFB Membrane Transport.-This section discusses transport in the reference system (see Table II) consisting of the quintessential VRFB membrane, Nafion, under typical electrolyte concentrations. For Nafion in a VRFB, Fig. 2a plots the mean calculated diffusion coefficients D ij on a heat map colored from dark blue (large positive coefficients) through white (coefficient equal to zero) to dark red (large negative coefficient) and the thermodynamic factor c ij is given in parenthesis. The calculations show that transport in this system is far more complex than diffusion in an ideal, dilute solution. Although some off-diagonal terms of the diffusion matrix are small, they are not categorically negligible and in some cases are larger than the on-diagonal terms. For example, a concentration gradient of ( ) V IV M 2 causes more than twice the flux of ( ) V V compared to an equal concentration gradient in ( ) V V M. A few of the diffusion coefficients are negative, indicating that the flux of one species will take place up the concentration gradient of another, holding all else constant. In particular, for all vanadium species i, D i SO 4 is negative and large in magnitude. Thermodynamic nonidealities contribute to transport coupling in particular. c ij for vanadium species i and j are close to ideal because positively charged ions do not significantly change the excess chemical potential of other positively charged species. c ij for vanadium/water and /sulfate are mostly large because of the significant effect interactions between these species have on the excess free energy of the system. Figure 2b shows a heat map of charge number-normalized transference numbers. Absent concentration gradients, protons will carry most of the current because = Transport at open circuit.-The concentration difference between electrodes creates concentration gradients in the membrane. Figure 3 shows the mole-fraction profile of each species in the membrane between the interfaces with the negative electrode ( = z 0) and the positive electrode ( = z l M ) calculated using the procedure outlined in Appendix B with the no current density (i.e. = i 0). For transport that is not coupled (i.e. dilute-solution  Table II). approximation), the mole-fractions of the species will decrease linearly between the two electrodes. This behavior is not present in this system, indicating that transport coupling plays an important role.
The vanadium species mole fractions in the negative electrode (V(II) and V(III)) goes to zero on approaching the positive electrodeinterface, and the vanadium species mole fractions in the positive electrode (V(IV) and V(V)) decrease towards the negative electrode/ membrane interface. The mole fractions of protons, sulfate, and membrane fixed-charged groups are relatively constant across the membrane with the latter two showing a slight decrease/increase towards the negative electrode and having a maximum/minimum away from the electrode interfaces. The dotted line in Fig. 3 shows and illustrates that the assumption equating the two is reasonable.
The water mole fraction decreases from the negative electrode to the positive electrode/membrane interface. Although the water mole fraction of the two electrode electrolytes are nearly the same, G 0 at the negative electrode is larger than at the positive electrode causing a concentration gradient inside the membrane.
The concentration gradients (see Fig. 3) drive fluxes according to the matrix of diffusion coefficients (see Fig. 2). Using ( ) V II as an example, Fig. 4 plots the contribution to ( ) N V II (see Eq. 19) from the concentration gradient of species j referenced to the membrane (i.e.  ). As the large magnitude off-diagonal elements in Fig. 2 show, concentration gradients of species i causes transport of species j. This finding is further evidence that dilute-solution theories in which a species transport is driven by only its own concentration gradient (i.e. is not sufficient to calculate transport accurately. To quantify this error, the term in parenthesis on the left side of Fig. 5 shows that the dilute-solution approximation (i.e. that the only friction on a species is due to the membrane, D  ¥ ¹ ij M ) predicts fluxes that, on average (i.e. mean), deviate from the concentratedsolution model by 772%.
For all species, the principle transport driving force is the water concentration gradient. Water drives positively charged species to the positive electrode, resulting in a curvature of the vanadium species concentration gradients towards the positive electrode. To maintain electroneutrality, the water concentration gradient pushes sulfate to the negative electrode.
Mass-transport of ions and water is coupled because both the frictional interactions between species are nonzero and the excess chemical potential of species depends on the concentration of other species. When describing the solution as thermodynamically ideal (circles in Fig. 5), the behavior of the coupling between species is observed, since the thermodynamic nonidealities lead to quantitative differences from the ideal-solution approximation. i 0 for the reference system (see Table II). Note the scale change in the mole-fraction axis. In practice, an operating cell will undergo numerous cycles and electrolyte rebalancing such that the composition in the electrodes varies, changing the concentration gradients. Consequently, the specific direction and contributions to fluxes in Fig. 5 do not apply throughout VRFB operation. Incorporating concentrated-solution theory into a fully transient model of VRFB operation would give insight into how these driving forces change during cell operation.
Transport under applied potentials.-Under an applied potential, charged species migrate with the current. Figure 6 shows the flux of each vanadium species i (multiplied by z F i ) as a function of current density i. When current is positive, the concentration gradients acting on V(II) and V(III) are aligned and increasing current linearly increases vanadium flux. At negative currents, migration and diffusion are opposed causing the fluxes of V(II) and V(III) to go to zero. V(V) and V(IV) follow opposite scenarios than V(II) and V(III). This description is qualitatively consistent with previous work using a dilute-solution theory framework. For these electrolyte concentrations, the diffusional contributions causes a net flux of vanadium towards the positive electrode when the current is zero as shown in Fig. 5. Figure 7 shows the mean contributions outlined in Eq. 19 to the V (II) flux ( ) N V II from each diffusion term, Fig. 5) and Figure 6. Vanadium flux N i multiplied by charge number z i and Faraday's constant F as a function of current density i for the reference system (see Table II). V II averaged across the membrane as a function of current density for the reference system (see Table II). ) for the reference system (see Table II V II as a function of current density. As the current density increases, the fraction of flux due to migration increases. The contributions to flux from diffusion change slightly as the current density increases because migration alters the concentration gradients across the membrane. Above 250 mA cm −2 , the majority of the vanadium V(II) flux is due to migration. Figures 4, 5, and 7 show that coupling plays an important role in transport though VRFB membranes. The contributions from water concentration gradients are particularly large. 44 The calculations show that in cases that require only semi-quantitative predictions, neglecting the terms ¶ D c ij j,M for ¹ j i, 0 will cause ∼25% error to the calculated fluxes at OCP with the error decreasing with increasing current density. As a result, mathematical models and experimental analysis that rely on extended versions of dilutesolution theories that incorporate the effect of water on vanadium fluxes will be qualitatively correct in many cases.
Cell performance.- Figure 8 shows that the calculated (lines) and measured 3 (symbols) voltage and Coulombic efficiencies, e v and e , q respectively, as a function of current density. The calculated and measured efficiencies correspond to slightly different scenarios. The model calculates performance at a constant 50% SOC (i.e. transport coefficients and concentration gradients are constant), whereas in the experiments charge and discharge between two voltage windows corresponding to different SOCs. Despite these slightly different conditions, the model and experiment are in good agreement. As the current density increases, the Coulombic efficiency generally increases as the effects of crossover are reduced in comparison to the energy extracted from the cell. However, rising current density decreases voltage efficiency because the membrane and cell incur larger ohmic losses. As a consequence, the energy density increases initially, reaches a maximum, and then decreases with increasing current density. In the absence of resistance from the membrane = R 0 M (i.e. DF = 0), e v is higher, but non-membrane resistances limit the benefits of ultra-high conductive membranes. Figure S2 shows voltage and Coulombic efficiency calculated using dilute-solution theory that only accounts for friction with the membrane (i.e. D  ¥ ¹ ij M ) and an extended dilute-solution theory that also accounts for frictions with water (i.e. D  ¥ ¹ ij 0,M ). The dilute-solution theory severely miscalculates cell efficiencies. For Coulombic efficiency, the extended-dilute-solution theory closely agrees with concentrated-solution theory but slightly overestimates voltage efficiency. Based on this analysis, researchers should account for the effect of coupled-transport modes and, in particular, the role of ion-water interactions on cell performance.

Impact of System Properties
The molecular transport and thermodynamic models presented in Parts I and II calculate properties as a function of external electrolyte composition and membrane chemistry. 32,33 The following sections use these models to determine how changes away from the reference system (see Table II) to operating conditions and membrane properties affect performance.
Impact of electrolyte concentration.-In this section, we calculate the effect of changing the composition of the electrolyte in the electrodes. Specifically, varying the vanadium concentration (i.e. ) between 3.95-4.82 mol dm −3 . This analysis neglects vanadium solubility limits in sulfuric acid that may be exceeded under certain conditions at high sulfate concentrations. 41 Furthermore, we assume fast mass transfer from the bulk electrolyte to the carbon electrodes of the VRFB; at low vanadium concentrations mass-transport limitations may be important depending on the cell design. 35,38,51,55 Transport Properties Concentration of species in the external electrolyte solution strongly impact membrane transport properties. Figure 9 shows membrane conductivity k (a), V(II)0 diffusion coefficient and V(II) transference number ( ) t V II (c) calculated at membrane composition¯d x i as a function of the total vanadium molarity in the electrodes, and the lines are varying arithmetic average sulfate concentration in the electrodes. Figure S3 shows the membrane water content (a), proton molality (b), and V(II) molality (c) at the same conditions. For brevity, we plot the properties of V(II), which are representative of the transport properties of the other vanadium species.
As the concentration of vanadium in the external electrolyte rises, the proton content in the membrane decreases, as Fig. S3b shows. The loss of protons, which have high mobility, decreases membrane conductivity. As Fig. S3c shows, the vanadium content in the membrane also increases with rising vanadium concentration in the electrolyte leading to larger vanadium transport properties ( ) D V II 0 and ( ) t . V II Adding sulfuric acid increases the membrane proton concentration. Higher proton concentrations displace vanadium from the membrane, as Figs. S4b and c show. Consequently, at high vanadium concentrations in the electrode, rising sulfate concentration increases conductivity and decreases vanadium transport properties, as Fig. 9 shows. However, adding sulfuric acid in the external electrolyte also dehydrates the membrane, as Fig. S3a shows and Part I discusses. 32 At low vanadium concentration in the electrodes, the dehydrating effects of additional sulfate dominates, causing membrane conductivity to decrease. Figure 10 shows the calculated Coulombic and voltage efficiencies, e q and e v respectively, (a) and discharge power density Y d (b) at the current density that maximizes Y d at e e =80% (termed optimal power density) as a function of the vanadium concentration and for different mean sulfate concentrations in the electrodes. To show how sensitive the optimal power density is to the target energy efficiency, shaded regions in Fig. 10b show power density at ±0.2% of the target energy density. As Fig. 9 shows, rising vanadium concentration in the electrode electrolyte increases the vanadium transference numbers and diffusion coefficients, decreasing Coulombic efficiency. To maintain an 80% energy efficiency, the current density decreases to lower ohmic losses and increase voltage efficiency. The Coulombic efficiency e q (circles), voltage efficiency e v (squares), and energy efficiency e e as a function of current density for the reference system (see Table II). Dot-dashed line shows e v for DF = 0 (i.e. without ohmic losses in the membrane).

Cell Performance
cumulative effect lowers the optimal power density with increasing vanadium concentration. Moreover, the lower Coulombic efficiency at high vanadium concentration makes the optimal power density less sensitive to the required energy efficiency (i.e. narrower shaded region in Fig. 9).
Increasing sulfuric acid concentration in the electrolyte (shown in Fig. 10 as moving from dot-dashed to dashed to solid lines) exchanges vanadium in the membrane for protons, increasing membrane conductivity and decreasing vanadium transport coefficients, as Fig. 9 shows. The larger resulting Coulombic efficiency allows the system to incur more ohmic losses by operating at higher current densities while maintaining the target energy efficiency. As a result, increasing sulfuric-acid concentration in the electrolyte increases the optimal power density. However, at low vanadium concentrations, increasing sulfate concentration lowers membrane conductivity, as Fig. 9 shows. The addition of sulfuric acid will therefore be less effective at low active-species concentrations. Moreover, the improved optimal power density at low vanadium and high sulfuric-acid concentrations may not result in higher return on investment for VRFB systems if balance-of-plant costs increase due to lower energy density or higher pumping costs. 42 Impact of membrane properties.-The chemistry of the membrane impacts its transport properties. Here we focus on two common changes to the membrane: equivalent weight (EW, mass dry polymer per mole of charged group in the membrane, i.e., inverse of the ion-exchange capacity) and membrane modulus. EW is tuned by changing chemistry. 6,7,9 For PFSA membranes, the EW is typically in the range 700 to 1500 g mol −1 . 46 PFSA membrane pretreatment changes its modulus by varying the number of physical crosslinks in the materials or its crystallinity. 2,8,56 Annealing tends to increase and boiling tends to decreases the membrane modulus. 56 In practice, EW and modulus are not independent because lower EW membranes tend to have lower moduli because the additional ionic groups on the polymer disrupt crystalization and crosslink formation. 46 Transport Properties Increasing EW decreases membrane water content, as Fig. S4a shows, thereby decreasing hydophilic water volume fraction, increasing tortuosity of the hydophilic domains, and decreasing the pore size. As a result, membrane conductivity k and V(II)/water diffusion coefficient ( ) D V II 0 decrease, as Figs. 11a and 11b show. Furthermore, increasing EW decreases the concentration of vanadium and protons in the membrane because there are fewer oppositely charged sulfonate groups for these ions to interact with. As EW increases, the vanadium concentration decreases more  rapidly than proton concentration in the membrane, as Fig. S4b and c show, decreasing the vanadium transference number, as Fig. 11c shows. Figure 11 shows (moving from dot-dashed to solid lines) that increasing the membrane modulus impacts transport properties similary to increasing EW. Raising the membrane modulus decreases water content, as Fig. S4a shows. The lower water content of the higher modulus membranes results in lower membrane conductivity and vanadium-water diffusion coefficients. As with EW, increasing the membrane modulus decreases the vanadium content in the membrane more rapdily than the proton content, decreasing the vanadium transference number. Figure 12 shows the calculated Coulombic e q and voltage e v efficiencies (a) and discharge power density Y d (b) at the current density that maximizes Y d at e e =80% as a function of equivalent weight EW and dry-membrane modulus E . b 0 As Fig. 11 shows, increasing EW decreases the vanadium transference numbers and diffusion coefficients, increasing Coulombic efficiency. However, conductivity also falls with increasing EW, lowering voltage efficiency and cell potential. The benefits of higher Coulombic efficiency do not outweigh the costs to voltage efficiency and cell potential. The net effect of increasing EW is a lower optimal power density.

Cell Performance
Decreasing the membrane modulus (shown in Fig. 12 as moving from dot-dashed, to dashed, to solid lines) increases membrane conductivity, vanadium diffusion coefficients, and vanadium transference numbers, as Fig. 11 indicates. At high EWs, the resulting increase in cell potential and voltage efficiency more than compensates for lower Coulombic efficiency, leading to a higher optimal power density. These findings agree with experimental studies of cell performance with different EW membranes and varying mechanical reinforcement. 6,9,30 At low EWs (<900 g-polymer/mole SO 4 ), ameliorating the poor Coulombic efficiency by increasing membrane modulus offsets the voltage-efficiency losses and results in highest optimal power density.

Structure-Property-Performance Relationships
The model results show that transport properties in VRFB membranes (e.g. conductivity, diffusivities, etc) result from a complex interplay of multiple phenomena that are mediated by structural characteristics of the membrane (e.g. membrane modulus and EW). These transport properties, in turn, dictate VRFB performance (e.g. Coulombic efficiency). To provide general structure-property-performance relationships, we perform sample-based sensitivity analysis. This analysis (see Numerical Implementation) calculates properties and the performance of membranes over a range of equivalent weights, moduli, tortuosities, and hydrophilic domain shapes and sizes. This method is analogous to synthesizing a thousand membranes with attributes spanning this parameter space and measuring their properties and in-cell performance. Figures S5, S6, and S7 show pair plots of structural attributes (EW , E , b 0 c, G, and d 0 ), exemplar transport and uptake properties  (k, t , II D , 0II D , IIII m , II m , H and l), and cell performance at a current density of 200 mA cm −2 (Y , d e , q and e v ) and reference conditions (see Table II). Figure S8 gives the least-squares linear regression coefficients between the natural log of each variable. Based on the sample-based sensitivity analysis, Fig. 13 shows a chord diagram of correlation between (a) structural attributes of the membrane and transport and uptake properties and (b) membrane properties and cell performance metrics scaling the chord size by the Spearman rank correlation coefficient r s squared (Fig. S9 provides all values of the Spearman coefficients between variables). Spearman coefficient is a measure of correlated monotonicity between variables. Figure 13a illustrates that dry-membrane domain spacing, d , 0 and membrane modulus, E , b 0 are highly correlated to the molality of ions in the membrane and water uptake, respectively, because they directly affect the energetics of membrane partitioning (see Part I). 32 The shape and tortuosity of the hydrophilic domains, quantified by G and c, respectively, correlate with transport properties. The domain shape affects the transference number of vanadium, whereas the tortuosity influences the conductivity and diffusion coefficients. EW correlates with all of these membrane properties, but typically to a lesser extent than the other membrane attributes. Figure 13a suggests which membrane design approaches are likely to influence a property of interest. However, the structure-function correlations are not high (i.e. at most 0.8), indicating that varying the structure of the membrane does not guarantee altered properties. Moreover, structural attributes impact multiple transport properties such that properties cannot be tuned independently. Figure 13b correlates the membrane properties with in-situ VRFB performance. Membrane conductivity directly governs voltage efficiency, e , v and discharge power density, Y , d at a current density of 200 mA cm −2 (see Eqs. 1, 6, and 28), giving = r 1 s between these variables. Because conductivity is highly correlated with vanadium diffusivity, vanadium diffusion coefficients are also correlated with e v and Y d but do not have a causal effect (see Fig. S7). The sensitivity analysis shows that Coulombic efficiency is highly correlated with vanadium transference number ( =r 0.98 s ) and water-vanadium diffusion coefficient ( =r 0.84 s ) and only weakly correlated with vanadium-vanadium diffusion coefficient ( =r 0.55 s ). These correlations are consistent with Fig. 7 in showing the crucial roles of migration and water-driven transport in cell operation. Water and ion uptake are weakly correlated to cell performance metrics. As such, screening membranes using ex-situ measurements of conductivity and vanadium transference number are most indicative of cell performance in these conditions.

Conclusions
This work develops a multicomponent, concentrated-solution model of transport in VRFB membranes that accounts for Stefan-Maxwell-Onsager transport couplings between species and thermodynamic nonidealities. The molecular-thermodynamics model outlined in Part I calculates ion and water partitioning into the membrane from the electrolyte in the electrodes and provides thermodynamic factors that influence species diffusion. The molecular-transport model outlined in Part II describes binary-diffusion coefficients between ions, water, and the membrane as a function of water content and composition. 32,33 The resulting multicomponent diffusion coefficient matrix has large off-diagonal elements that dilute-solution theory neglects and that impact the net transport of the species through the membrane. Transport is highly coupled because of both thermodynamic and frictional interactions between species. In particular, water concentration gradients play a dominant role in vanadium crossover. Under an applied current, contributions to net flux from these diffusional modes and migration are aligned or opposed, depending on the direction of transport. Migration dominates at high current densities (>250 mA cm −2 ). The magnitude of coupling suggests that predictions made using dilute-solution theories for VRFB will be misleading. Extended-dilute-solution theories that account for the impact of water chemical potential gradients on vanadium transport are more reliable. Furthermore, membrane-permeation experiments that measure flux under illdefined concentration gradients will not provide diffusion coefficients that are predictive of transport under other conditions. Even in the case of carefully designed experiments, the magnitude of measured diffusion coefficients should be interpreted in the context of the collective interactions between all species present in the system. Sensitivity analysis indicate that membrane conductivity and vanadium transference number are the best predictors of cell performance metrics.
The cell performance in terms of power density and energy efficiency depends on the complex interplay of the thermodynamic and frictional interactions of all species and their gradients. We studied the effect of sulfate and vanadium concentration in the electrode electrolyte and membrane equivalent weight and modulus on membrane properties, species transport, and cell performance. For a metric of maximum power density at 80% energy efficiency, low-vanadium and high-sulfuric acid concentration electrolytes perform best. However, the solubility of vanadium in sulfuric acid limits the success of this strategy. Membranes with low equivalent weights (i.e. high ion-exchange capacity) perform better. Specifically, decreasing EW from 1500 to 900 improves performance by 20% by increasing membrane conductivity. At high equivalent weights, a low modulus is superior, although perhaps not readily synthesizable. where I is the identity matrix, the superscript -1 denotes the matrix inverse, and the exp is the matrix exponential. The concentration gradient vector is then We calculate the matrix exponential by diagonalizing the matrix argument and computing the exponential element-wise of the matrix eigenvalues. Eqs. B1 and B2 give the species concentration and concentration gradients as a function of species fluxes.
Equation B2 In conjunction with Eq. 19 calculate flux. The equations are an implicit set of algebraic equations that calculate flux. Krishna et al. 52,53 showed that the following numerical method specifies flux: (1) calculate Y for = N 0; (2) calculate x' from Eq. B2 at = z 0; (3) use calculated x' at = z 0 to calculate N using Eq. 19; (4) recalculate Y using the calculated N; (5) repeat steps 2-4 until achieving convergence (here, set to 1 × 10 -3 relative change in flux over an iterations).