The optimal electrode pore size and channel width in electrochemical flow cells

time and expedite the industrial design process. The optimal channel width and pore sizes we obtain, in the order of 100 microns and 1 micron respectively, are much smaller than those often used. This shows that there is a significant room for improvement of energy efficiency in flow cells that can sustain the resulting pressure drop.


Introduction
The intermittent nature of wind and solar based renewable energy sources demands the usage of energy storage devices that act as a reservoir or source of energy during times of high availability or demand of energy, respectively. These storage devices operate in small scales (kW/kWh) for individual usage to grid scales (GW/GWh). Electrochemical energy storage systems like redox flow batteries, electrolyzers, and fuel cells are attractive candidates for grid-based electrical energy storage since their storage and power capacities can be scaled up independently [1,2]. This results in versatile storage and power capabilities [3]. However, these electrochemical devices are subject to energy efficiency losses that need to be minimized to facilitate their wide-scale acceptance.
Electrochemical systems often utilize an ion exchange membrane to hinder cross-over of chemical species and avoid self-discharge. This membrane adds significantly to the cost and ohmic resistance of the battery [4,5]. Hence, 'membraneless' electrochemical systems were proposed to not only reduce costs but also to reduce the overall efficiency losses resulting from the ohmic potential drop across the membrane [3,6,7].
A common design used in membraneless systems relies on a laminar flow between the two electrodes to keep the electrolytes at the anode and cathode separated. The laminar flow ensures that only molecular diffusion is responsible for the transverse transport of species in the channel, while advection in the axial direction removes any reaction products in the channel. Different battery chemistries like-vanadium [6][7][8][9][10][11], hydrogen/bromine [12], alkaline hydrogen/oxygen [13] and formic acid/oxygen [14][15][16] systems have been investigated using such membraneless designs. Laminar flow channels have also been studied extensively for their use in micro fuel cells (MFC) or membraneless microfluidic fuel cells (MMFC) [17,18]. However, these designs are difficult to scale up as the width of the diffusion zone increases with the length of the channel. Hence, if the length of the channel is too long reaction products of the cathode will reach the anode and vice-versa. A nanoporous separator was added to further suppress the crossover of species in hydrogen/bromine [19], boronhydride/cerium [20] and methanol [21] systems. Flow channels are also used in electrochemical CO 2 -reduction to ensure hydration of the catalyst layer and to remove products [22][23][24][25].
3D porous electrodes are another common feature in electrochemical energy storage systems. In a porous electrode the total current density is an accumulation of the local current densities along the A. Bhadra and J.W. Haverkort current direction. A porous electrode with a substantially high surface area has a lower magnitude of local current density. Therefore, the overpotential for such an electrode is lower than for a smooth planar electrode where the local current density is equal to the total current density [26,27]. The working of porous electrodes in different electrochemical devices has been thoroughly explored in numerous previous works, see for example [26,[28][29][30][31][32][33][34][35][36]. Physical parameters such as porosity, internal pore diameter, and thickness characterize these electrodes and can be tailored to improve their performance [37][38][39][40]. Unlike in the 'flow-by' flow configuration, in 'flow-through' and 'interdigitated' designs the flow of reactants is through the porous electrode itself. Hence, significant pumping energy losses can occur in both these designs.
In addition to pumping losses, laminar channels and porous electrodes can also give rise to ohmic dissipation [41][42][43][44]. Additionally, electrodes give rise to activation losses. A proper design of their geometry can strongly reduce these losses. In the present work we present an explicit analytical expression for the optimal laminar flow channel width that minimizes the sum of pumping losses and ohmic dissipation. For porous electrodes, we present expressions for the volumetric surface area to minimize the combined pumping and activation losses in flowthrough and interdigitated flow configurations. We then verify these expressions using 2-dimensional numerical simulations in COMSOL Multiphysics.
Typical values of electrolyte channel widths used in experiments are of the order of 1 mm [10,11,19] and the pore size in the electrode is in the order of 10 microns [19,23,45]. We find that the optimum values of these parameters are usually an order of magnitude lower. This means that significantly less ohmic and activation losses can be achieved by further miniaturization. This does lead to higher pressure drops, requiring more careful engineering. It seems to be a general trend in many research papers in the literature that a higher efficiency is sacrificed for a more practical pressure drop. We hope that our work draws more attention to this underappreciated aspect of cell engineering and can further improve performance optimization.
The rest of the paper is organized as follows: in Section 2 we present the equations used for the numerical simulations and on the basis of which we develop analytical models; in Section 3 we put forward our analytical models along with the assumptions used in the optimization; in Section 4 we compare the obtained analytical expressions with the results from numerical simulations; in Section 5 we compare our results with various experimental results from the literature, and in Section 6 we summarize our work and present the final conclusions. Fig. 1 gives a pictorial description of the geometry considered to optimize the electrolyte channel width, ch and the volumetric surface area of the porous electrode, . Two-dimensional simulations are performed in the x-z plane. To optimize the specific surface area , we consider a porous electrode through which a mixture of the reactant, product, and electrolyte flows. Instead of a membrane, we have an electrolyte channel through which the ionic current can flow. The electrolyte flows between two nanoporous separators that inhibit advection and dispersion to and from the channel. Diffusion of reactants and products through the nanoporous separator into the electrolyte channel is facile due to the thinness of the layer, but flow through the nanoporous separator is strongly inhibited due to its high hydraulic resistance. An example of this design, showing an excellent power density of roughly 1 Wcm −2 and current density of 3 Acm −2 , can be found in Ref. [19].

Governing equations
The fluid flow in the open channel and the porous electrode are governed by the laminar steady Navier-Stokes equation, Eq. (1), and the Brinkman equation, Eq. (2), respectively. The incompressibility condition, Eq. (3), is assumed to be true in all regions of the flow.
Here, is the superficial flow velocity is the pressure, is the density, is the dynamic viscosity, is the porosity and is the permeability of the porous medium. The term = f √ with f = 1.75 is based on Ref. [46]. The inertial terms, ( ⋅ ∇) 1 2 and | | , due to the linear naure of the flow in the channel and small pores in the electrode, do not make significant contributions.
We use dilute solution theory with constant transport coefficients. The steady Nernst-Planck equation, ∇ ⋅ = is used for the transport of the th dilute species. The source term is written using Faraday's law as = − ∇⋅ , where is the number of electrons transferred per molecule and is the ionic current density vector. The Nernst-Planck flux reads Here, , , , are the flux, concentration, diffusivity, and charge number of the th species respectively; while and  denote the Faraday's constant and the gas constant respectively. For noncharged species, = 0 and a standard advection-diffusion-reaction equation results. Inside porous media, the effective diffusivities are related to the molecular diffusivities m, by = m, 1.5 corrected using Bruggeman's relation derived for transport through polydisperse spherical particles [47]. In the flow channel, the porosity = 1 and the effective diffusivity equals the molecular diffusivity m, .
The ionic current density perpendicular to the flow of electrolyte, , is given by For a binary electrolyte, due to quasi-neutrality, + = − = el so that advection does not contribute to the current.
Using Eqs. (4) and (5), in the absence of gradients in the electrolyte concentration, we have Here is the ionic conductivity and reads Here, + and − are the effective diffusivities of the positive and negative ions so that we can write = m 1.5 , where m is based on the molecular diffusivity. The electronic current density, , in the porous electrode is given by Ohm's law as Here, the effective medium electronic conductivity = m ( 1 − ) 1.5 is obtained from the material conductivity m using Bruggeman's relation for the solids fraction 1 − , and is the electronic potential. The superficial electronic current density and ionic current density are related by ∇ ⋅ = −∇ ⋅ . The ionic current entering the porous electrode is related to the local current density ⟂ according to where is the volumetric surface area of the porous electrode. We assume that the symmetric Butler-Volmer equation holds, which reads Here * is the exchange current density and red & ox are the concentrations of the reducing and oxidizing agent respectively, and red,in & ox,in are their inlet concentrations, and the Tafel slope = 2 . Taking = ,in , the concentration-dependent Butler-Volmer Eq. (10) can be solved for to give For low activation overpotentials , the relation becomes linear, while for high overpotentials the logarithmic Tafel equation applies.

Results: Analytical modelling
In this section, we develop an analytical model for the total relevant power loss [W], that can be optimized analytically. In order to do so, we make some simplifications. We assume that the reactant concentration in the porous electrode is uniform in the -direction. We also assume a uniform local current density ⊥ in the -direction.
The net power loss in the system is written as Here, act is the activation loss, fr is the pumping loss due to friction in viscous flow and res is power loss due to ohmic dissipation. These losses reflect the losses in just those parts of the cell that we aim to optimize, namely the flow channel and porous electrode. We have ignored losses due to concentration depletion. In Section 3.3.1 we have included mass transfer limitations to optimize the volumetric surface area of the porous electrode. We take the cross-sectional area of the channel/porous electrode as ch/pe = ch/pe and the cross-sectional area of the electrode as = ℎ. The current density magnitude ( ) may vary in the flow direction, due to reactant depletion, so we define an average value ⟨ ⟩ = 1 ℎ ∫ ℎ 0 . Assuming the average of the squared current density approximately equals the square of its average, ⟨ 2 ⟩ ≈ ⟨ ⟩ 2 , which is true if ⟨ ⟩ is uniformly distributed along the -direction, we show in Appendix A that the ohmic power can be approximated by We use the conductivity at the inlet in our model. This is valid if the electrolyte concentration does not vary significantly along the flow direction. The dissipation per unit area due to the average activation overpotential ⟨ ⟩ = ∫ pe 0 in case of a local current density ⊥ that is constant in the -direction reads This notation for the dissipation per unit area is understood to be valid also locally, in case varies with . The pressure gradient, ∕ℎ for a unidirectional laminar flow in a channel of length ℎ and average velocity ⟨ ⟩ is given by  [48,49]. The hydraulic pore diameter is related as h = 4 ∕ = 1− f so that for ≈ 0.25 the inverse 1∕ can be seen as a rough measure of a typical pore diameter.
In Appendix A we show that the frictional dissipation in the channel (ch) or porous electrode (pe) reads The average velocity ⟨ ⟩ is often chosen on the basis of the 'conversion' , which denotes the fraction of reactants entering the cell that is consumed.
We define as is the velocity-averaged or 'cup-mixing average' reactant concentration. The weighting with velocity accounts for the fact that the concentration of the faster flowing liquid in the centre of the channel contributes more than the slowing moving liquid near the walls. We hope that the only slightly different font used here for this quantity does not lead to confusion.
Since for every electrons a reactant molecule is used, we can relate ⟨ ⟩ and as follows Inserting Eq. (19) into Eq. (17) and Eq. (16), we have, A pump energy efficiency can be easily included by dividing the viscosity with it.

Optimal electrolyte channel width at given
In this section, we consider the optimization of a flow channel between the two electrodes, two separators/membranes, or between an electrode and a separator/membrane. This is relevant for membraneless and microfluidic electrochemical cells. We assume that a reactant in the channel achieves a particular degree of conversion, , but this can also be replaced by a certain outlet concentration of the reactant or product. We note that with very close to 1, concentration overpotentials may no longer be negligible. In such cases it remains to be seen how well our analysis will hold.
The electrolyte channel width, ch , is one of the factors influencing pumping and ionic resistance losses. The power loss impacted by the electrolyte channel width is defined as Keeping the flow-rate and the average current density constant and increasing the channel width increases the power loss due to ohmic dissipation but reduces the power loss due to friction. Hence, an optimum width exists that minimizes ch .
Using Eqs. (13) and (20), the overall power loss in the electrolyte channel reads, The optimal gap width, ch,opt is obtained for that ch which minimizes ch . Solving ch ch = 0 for ch gives,  [50] the flow-rate instead of the conversion was taken as an independent variable that does not change during the optimization, giving a similar weak dependence on and . Because the relation between flow-rate and velocity depends on the parameter ch/pe that is optimized, this relation cannot a posteriori be inserted as is done in Ref. [50]. Our analysis with as the independent variable is only valid for cases where the electrolyte is utilized in the reaction. Otherwise ⟨ ⟩ should be used as an input parameter for the model. A similar result to Eq. (23) was also obtained in Ref. [51] for a flow-through cell with flow parallel to the current.
We have not included a membrane in our channel, however our model would hold true even in the presence of a membrane.

Optimal electrolyte channel width avoiding boundary layer overlap
In case a flow channel is used to avoid cross-over of reactants or products it should be sufficiently wide and the velocity sufficiently high to avoid the boundary layers arising on opposite sides of the channel to overlap. In Appendix C we derive the conditions that minimize the sum of ohmic and pumping losses. We assumed that the species in both boundary layers have equal diffusivity and allow at the end of the channel 1% of the concentration tail to overlap. The optimal velocity required to ensure this reads and the associated optimal channel width, which reads ch,opt = 6.3 Note that the viscosity and conductivity enter Eq. (25) very weakly with a power of merely 1∕6. With = 1 mPas, = 1 S/cm, and = 10 −9 m 2 /s Eq. (25) gives ch,opt ≈ 43 μm ( and, from Eq. (24), These optimal values will typically lead to very small microfluidic gaps and relatively high velocities. Some example values for lab-scale cells are show in Table 3 and will be discussed later in Section 5.

Volumetric surface area
Next, we optimize the volumetric surface area of porous electrodes. A flow-through electrode is considered in which the electrolyte flows through the porous electrode, parallel to a membrane or a microporous separator. In Section 4.2.2 we apply this model to an interdigitated flow field. Increasing the volumetric surface area, , reduces the activation losses but increases the pumping losses. Thus, an optimum exists that minimizes the sum of these losses. We will neglect the flow through any additional surface area that may be present in the form of, roughness, a coating, or inside a fibre bundle of a carbon cloth. This internal surface area will effectively increase the exchange current density * of the external surface area, considered here. This allows the use of the same symbol and value for the surface area in Eq. (9) and the hydraulic surface area in Eq. (16). First, the concentration-independent case will be considered and the next section deals with the effect of mass transfer limitations.
The power loss influenced by is given by The ohmic drop is not influenced by the volumetric surface area and as such is not present in the expression for power loss influenced by . Note that Eq. (19), relating current density and velocity, does not involve . Therefore, we can keep the velocity as an independent variable and later inserting Eq. (19), without influencing the optimization. As a result, the optimization can be performed for a particular position . Therefore, we use the local current density instead of the channelaverage ⟨ ⟩ that we used before, and we obtain an optimum that will depend on through .
For the symmetric Butler-Volmer kinetics of Eq. (11), Eq. (26) becomes Eq. (A.7) which can be written for the dissipation per unit area as Herē= 2 pe * ,̄= act 2 0 asinh 1 = 2 0 , and Here * is the local exchange current density, so pe * is the effective superficial exchange current density [27]. The optimal value of the dimensionless surface areācan be obtained by solvinḡ= 0, which gives, This gives a cubic equation in̄2 that can be solved analytically bȳ . The limiting solutions for high and low values of̄are given bȳ A. Bhadra and J.W. Haverkort  Table B.6. The relative error between the approximate Eq. (29) and (33) is shown in the inset for different values. All powers between = 3 and = 10 show a reasonably low maximum relative error ≲ 10%, with a minimum of 0.6% error around ≈ 4.9.
which in dimensional form reads Note the absence of * in the Tafel regime, so that the optimum depends on the electrode kinetics only through the Tafel slope, . The optimal surface area in both cases increases with increasing porosity, as this increases the permeability, allowing smaller pores for the same pressure drop. A higher current density increases activation losses so favours more surface area. A higher velocity increases pumping losses, so favours less surface area.
Depending on the effective electrical conductivity and , beyond a certain thickness the overpotential will no longer be constant and the full electrode will no longer be effectively used. Ref. [27] derives the associated optimal electrode thickness. In Appendix D we show the resulting expressions for the optimal volumetric surface area in case that the electrode is also chosen to be of optimal thickness.
We can add the two limits in Eq. (32) with powers − < 0 to give most weight to the smallest one and obtain an expression that reduces to the proper limits and remains approximately valid in between: Fig. 2 shows that the relative error with the exact solution of Eq. (29) is minimized to 0.6% for = 4.9.

Effects of mass transfer in the Tafel regime
The mass transfer within a pore of the electrode is governed by the concentration profile near the pore wall. Especially, for high current densities, often in the Tafel regime, the effects of mass transfer limitations can become prominent. In order to include the effects of intra-porous mass transfer in the optimal volumetric surface area, we consider the concentration-dependent Tafel kinetics regime, obtained from Eq. (10) as = ln . Here, w is the concentration at the pore wall, see Fig. 3, and the reference concentration is taken to be the concentration at the inlet, in .
With m the mass transfer coefficient based on the concentration difference c − w , we obtain for the local molar flux, ⟂ Hence, the limiting current lim that can maximally be obtained when w = 0, is given by Combining Eqs. (36) and (37) we obtain, Combining Eqs. (14), (20), (34) and (38), we extend Eq. (27) in the Tafel regime into, where we define lim as Hence, ∕ lim = lim ∕ . To optimize Eq. (39) with respect tōwe will assume that lim is independent of̄. This introduces an error since typically the mass transfer coefficient m does weakly depend on the pore size.
Taking the derivativē= 0 gives which gives with̄o pt,Taf =̄1 ∕2 For̄l im = 0, we obtain the previous solution̄=̄o pt,Taf and if lim ≫̄1 ∕2 , we find̄o pt =̄l im . The final expression in Eq. (42) is a rough approximation with a maximum relative error of 25%. It shows that the optimal surface area is approximately the sum of the optimal surface area without mass transfer limitations and lim , so that it is mostly determined by the largest of these two contributions. In terms of pore size,̄o pt,Taf determines the optimum until they become too large and the maximum A. Bhadra and J.W. Haverkort   Fig. 4. Boundary conditions used in the simulations for optimizing (left) the channel width and (right) the volumetric surface area in a flow through configuration. For optimizing the volumetric surface area in the linear regime in the flow through configuration, the electrolyte channel and upper nanoporous separator were excluded from the geometry. This was done as the pressure-driven crossover, or leakage, of reactants into the electrolyte channel would significantly violate Eq. (19). For the Tafel regime, the leakage flux of the reactants less large compared to the higher current densities.

Table 1
The channel width ch and flow velocity ⟨ ⟩ are varied over the indicated ranges to make sure that ⟨ ⟩ and obtained within the simulations of each case are a constant. Correspondingly, the electric potential across the cell is also changed. Remaining parameters used in the simulations are shown in Table B For constant conversion the velocity and current density are proportional, so that the first term in Eq. (43) decreases with increasing while the second term, associated with concentration overpotential, increases. If the variation in and c with is known or can be estimated, Eq. (42) or Eq. (43) can be used to create a variable volumetric surface area that is optimal for all .

Results: Numerical modelling
We solve the tertiary current density distribution using the Nernst-Planck module available in COMSOL Multiphysics v5.5. In order to solve for the fluid flow, we use the Free and Porous Media Flow module. The boundary conditions are displayed in Fig. 4. The parameters used in the simulations are listed in the tables in Appendix B. These models allow us to take full account of concentration effects and mass-transfer limitations, and hence the predictions are more comprehensive than those of our simplified analytical model.

Optimal electrolyte channel width
For the optimal electrolyte channel width, numerical simulations are performed for the two cases shown in Table 1: For the 'low' case of low current density and conversion the assumptions of the analytical model are reasonably well satisfied. For the 'high' case, of high current density and conversion, reactant depletion leads to a highly inhomogeneous current and reactant distribution in the electrolyte channel, see Figs. 5(d) and 5(f). These conditions violate the assumptions made in the derivation of the ohmic power loss in Eq. (13), so will be a good test of the robustness and generality of the analytical results.
The power loss, ch = fr + res , for the numerical simulations is calculated using Eqs. It is seen that for both the 'high' and the 'low' cases, our analytical model underpredicts ch for higher channel widths (Figs. 5(a) and 5(b)). The effect is more pronounced for the 'high' case. This discrepancy is due to the violation of the assumption of constant current density. The ohmic losses described by the first term in Eq. (22) contain ⟩ in the exact result of Eq. (A.1). It can be shown that ⟨ ⟩ 2 ≤ ⟨ 2 ⟩ so our analytical expression underestimates the losses. In the 'high' case, as can be seen from Fig. 5(f), a limiting current is attained throughout most of the channel due to reactant depletion. This gives rise to a current density that, as shown in Fig. 5(d), drops dramatically over a very short distance near the entrance.
However, the calculated ch,opt of 7.83 × 10 −5 m and 3.64 × 10 −5 m from the analytical model equation (23) for the 'low' and 'high' cases respectively, are in the same order as those seen in numerical simulations; 7.91 × 10 −5 m and 3.01 × 10 −5 m respectively. It is seen in Fig. 5(a) that for the smaller channel widths, when the concentration remains close to the inlet concentration and the current distribution is relatively homogeneous, the results from the analytical and numerical calculations are similar.
From Figs. 5(a) and 5(b) we see that the power loss in the cell increases rapidly with decreasing channel width, ch below the optimum value. This is because the pumping losses scale inversely proportional to the cube of ch , as seen from Eq. (22). However, the power loss increase is more gradual when the channel width is higher than the optimum value due to the weaker linear dependence of ohmic losses on the channel width, see Eq. (22). Hence, the value of ch,opt should perhaps, as a design criterion, rather be used as a lower limit below which the power loss drastically increases.

Volumetric surface area
Next, we turn our attention to the determination of the optimal pore size or volumetric surface area of a porous electrode. We first consider the flow-through geometry of Fig. 4 (right) and next the interdigited flow-field geometry shown in Fig. 7. We assume ≪ lim so we can neglect the internal pore mass transfer limitations considered in Section 3.3.1. To test both limits of Eqs. (31) and (32) we consider the two cases shown in Table 2: a case of relatively low current density and high exchange current density for which linear kinetics will hold for all used values of and a case of high current density and a lower exchange current density so that Tafel kinetics will hold.
A complication is that for the small pores used in the nanoporous separator for the Tafel case the pressure drop in the electrode becomes much higher than that in the channel, resulting in a flow from the electrode to the channel. This poses a serious problem to membraneless cells with flow through the electrode normal to the current. To avoid this problem, in our simulations we chose an extremely small permeability of ns = 9 ⋅ 10 −24 m 2 for the nanoporous separator. We for a low average total current density of ⟨ ⟩ = 655 A/m 2 and high average total current density of ⟨ ⟩ = 3032 A/m 2 respectively . It is seen that in the case with a higher current density the power loss deviates from the analytical results for higher channel widths. This is attributed to current inhomogeneity, violating the assumption ⟨ 2 ⟩ ≈ ⟨ ⟩ 2 underlying Eq. (13). In (c) and (e) the normalized electrolyte current density and normalized reactant concentration are shown for ch = 2.34 ⋅ 10 −4 m for the low average total current density. In figures (d) and (f) these quantities are shown for the high average total current density at 1.4 ⋅ 10 −4 m. Note that in the 'high' case a limiting current is approached due to reactant depletion. This leads to a strongly inhomogeneous current distribution. for both the Tafel and linear kinetic regimes is compared to the analytically obtained results. In the simulation for the Tafel case, the leakage of the reactant through the nanoporous separator into the electrolyte channel adds to the deviation from the analytically obtained result. In the linear case the electrolyte channel was not a part of the simulation. In this case the deviation of our approximation for act , Eq. note that this is a much smaller value than commercial nanoporous separators have and actually unrealistically small using the formula of Eq. (16). Therefore, such a low value will likely require pores sizes of the order of the molecules. This predicts that porous separators are not effective in avoiding advection-driven cross-over. A membrane would be much better in this case as it can sustain larger differential pressures without as much liquid permeation. Therefore, we argue that usually a membraneless system in combination with flow-through electrodes is not possible without large advective cross-over, unless accurate pressure balancing is ensured between the channel and the electrode.
Another option to avoid strong advective cross-over is to making sure that the pressure in the electrode and channel are similar. Using Eqs. (15) and (16)

Table 2
The channel average current density and exchange current density used in the 'linear' and 'Tafel' cases studied to investigate numerically the optimal volumetric porous electrode surface area in a flow-through and interdigitated geometry. With pe = 0.75 mm, the superficial exchange current density pe * can be seen to be larger than ⟨ ⟩ for all values for studied in the linear case and smaller than ⟨ ⟩ in the Tafel case. In the linear case no flow channel was present. Remaining parameters used in the simulations are shown in Table B way that the pressures balance. Since in general the channel width is much larger than the electrode pore size, this will require a proportionally larger channel velocity. Flow-rate controllers activated by pressure difference sensors may in practice ensure that leakage can be largely avoided. In the interdigitated design the pressure drop is not linearly distributed over the channel so that this pressure balancing is impossible. Since the pressure drop over the electrode is usually much smaller in this design, this is however less of a problem.

Flow-through flow field
We verify Eqs. (27) and (31) using numerical simulations of the flow-through porous electrode configuration of Fig. 4 (right).
The power loss in the porous electrodes is calculated using Eq. (26). The analytical results are obtained by obtaining the overpotential, , for the Tafel regime and linear regime using Eq. (11). The average total current density, ⟨ ⟩, we keep constant. For the Tafel regime, we modify Eq. (19) to take into account the flux of reactants leaked into the electrolyte channel using Here, leak is the flux of reactants leaked into the electrolyte channel. This flux has been obtained from the numerical simulations. It can also be estimated using Eq. (15) by calculating the pressure drop between the channel and the electrode and using the permeability of the nanoporous separator. For the linear regime, the electrolyte channel was not included in the simulation geometry because the leakage flux was comparable to the flux of reactants producing electric current. This undesirable operating regime would cause significant deviation from the analytical model.
In Fig. 6 we see that the analytical prediction of act agrees reasonably well with the results from the numerical simulation. A minor discrepancy comes primarily from the difference between the exact and approximate expressions for act . The exact expression for the dissipation, Eq. (A.3), is an integration of the product of and ⟂ over the porous electrode area, while Eq. (14) multiplies the average of with an integration of ⟂ over the porous electrode area. For both the linear and the Tafel cases, the analytically obtained opt values A. Bhadra and J.W. Haverkort Table 3 A comparison of the electrolyte channel width and velocity used in experiments and the optimal values calculated using Eq. (25), assuming ( ) 1∕6 = 10 −2 N 1∕6 . This corresponds to for example = 10 −3 Pa s and = 10 −9 m 2 /s. The maximum experimentally used current densities are indicated and were used in the calculations. Clearly the predicted optimal gap thickness is much smaller, and the optimal velocity is much higher, than actually used in the experiments. of 3.1 × 10 4 m −1 and 1.7 × 10 5 m −1 respectively, are within 5% of the numerically obtained optima (2.9 × 10 4 m −1 and 1.6 × 10 5 m −1 ). The predictions can be improved by including the variation of and ⟂ in the model.

Interdigitated flow field
The interdigitated flow field is commonly used in electrochemical devices. Liquid flows in through a channel with a dead end, so it has to flow through the porous electrode to a nearby outlet channel. Fig. 7(a) illustrates this schematically. An advantage of this flow field is that the fluid travels over only a small fraction of the entire porous electrode length, resulting in a much smaller pressure drop compared to a flowthrough electrode where the entire electrode length is traversed. We will here illustrate how our analytical results for the optimal pore size can also be used for an interdigitated flow field.
We have to make a few modifications to use Eq. (27) for the power losses, we will do so by redefining some of the geometrical variables. Fig. 7(a) introduces the channel width and the distance r between rib centres. First, the distance over which the pressure drop arises changes to roughly ℎ = ( r + pe − ). This expression approximates with straight lines the average length of a path that a typical fluid parcel will traverse, in case the flow distributes well over the electrode. This trajectory describes roughly the middle of the dashed lines shown in Fig. 7(a) or the middle series of arrows in Fig. 7(d). Since for a well-designed interdigitated flow field the pressure drop over the electrode is much larger than that in the channel, the fluid will indeed usually distribute itself relatively homogeneously. In the dimensionless coefficient of Eq. (27),̄= ⟨ ⟩ r 2 0 , we take = r with r or the electrode area of just one of the repeating units shown in Fig. 7(a) between the blue dashed lines. Correspondingly, the frictional pressure drop fr will also be taken to correspond to only the flow indicated for this half-channel-pair so we take pe = ∕2. These replacements change Eq. (19) to ⟨ ⟩ = ⟨ ⟩ in r ∕2 and Eq. (17) becomes fr = . Therefore, 0 = fr ∕̄2, with̄= Fig. 7 shows the two-dimensional simulation results for an interdigitated geometry. The boundary conditions for the simulation are shown in Fig. 7(a). Only a part of the cell along the -direction is simulated. The symmetry in the flow field of the porous electrode is utilized to truncate the simulation geometry. Similar simulations were performed in [52] and the authors concluded that 2D simulations captured all the important features of the cell. A constant conversion of ≈ 0.048 was used and a current density of ⟨ ⟩ = 2 A cm −2 .
We see in Fig. 7(b) that our analytical model accurately predicts the power loss in the porous electrode. The analytical optimum value of 5.9 × 10 6 m −1 is in the same order of magnitude as the numerically obtained optimum of 5.1 × 10 6 m −1 .
Despite the complex flow configuration considered here, the analytical predictions turn out to be surprisingly good.

Comparison with literature experiments
The analytical expressions derived in Section 3 are used here to derive the optimal volumetric surface area for the porous electrodes used in some experiments from previously published literature.

Channel thickness
In Table 3 we consider three papers with an electrolyte channel. Refs. [23,45] are from the field of CO 2 electrolysis, while Ref. [19] considers a membraneless hydrogen-bromine flow battery. In all these publications it is desirable to avoid overlapping boundary layers growing from two sides of the channel. For this we derived the optimal velocity in Eq. (24) and the optimal channel thickness in Eq. (25). We immediately see that the used velocities are much lower than what is predicted to be optimal and the gap thickness much wider than the optimum. This shows that there is still a lot of room to decrease the ohmic losses while keeping the frictional losses acceptable. Due to the wide gaps and low velocities, the present pressure drops over the channel can be estimated to be always below 1 Pa. With the optimal parameters this will increase to ten to tens of millibars, which may be easily achievable. Therefore, our recommendations provide a simple and straightforward way to further improve these systems. In the case of Ref. [19] the pressure drop will be more than 1 bar, in which case it will be very hard to keep advective cross-over under control, even with the microporous separators present. Very delicate pressure balancing between the channel and the porous electrode will be required in this case, which may be the understandable reason why suboptimal conditions were chosen in this case.

Volumetric surface area
In Table 4 we show the parameters used for three papers using flowthrough porous electrodes. Unlike in our analytical model, Refs. [10,11] have flow inside the porous electrode parallel to the electric field instead of normal to it. Fortunately, the analytical model developed in our work can be applied to these situation as well. As seen in Eq. (17), fr is proportional to volume of the electrode through which the electrolyte flows so that it does not matter in which direction the flow is.
We see that, in the chosen examples, the volumetric surface area is several times, to more than an order of magnitude, smaller than is optimal. This means that the activation losses could have been substantially lower, while keeping the pumping losses acceptable. These more optimal conditions do come with a strongly increased pressure drop, of the order of 0.1, 0.3, and 1.6 bar for the three references in the table, respectively. Such significant pressure drops require stronger pumps and puts higher demands on seals, to ensure leak-tightness. In the case of Ref. [19] these pressure drops would additionally lead to excessive flows from the porous electrode into the channel. That is, unless the pressure at any position is exactly balanced by an equal pressure in the channel.
A. Bhadra and J.W. Haverkort Table 4 A comparison of the optimum volumetric surface areas and the parameter used in experiments and the optimal value calculated from Eq. (32) assuming Tafel kinetics with = 56.5 mV and KC = 5. The viscosities and the porosities that were not reported were estimated. The used volumetric surface areas for the used commercially available porous electrodes were obtained from [53].

Conclusion
We analysed the geometric parameters of an electrochemical flow cell and their effects on the power loss. Simple explicit analytical relations are provided to obtain the optimal values of these parameters that maximize the energy efficiency.
We studied the dependence of pumping power losses and ohmic dissipation on the electrolyte channel width. In the case of a constant conversion of reactants, we found that the optimum electrolyte channel width from the analytical model matches the results from the 2D binary electrolyte numerical simulations, even when the assumptions of a constant current density and electrolyte concentration are strongly violated, close to the limiting current density. It deserves recommendation to use this optimum as a lower limit in the design, since decreasing the channel width below the optimum rapidly increases the pumping losses.
We also obtained the optimum volumetric surface area, roughly the inverse of a typical pore size, considering the pumping and the activation overpotential losses. Both Tafel and linear kinetics regimes are considered. As with the electrolyte channel width, this model of activation losses is compared with the 2D model from simulations, which combine the Navier-Stokes equation, including the Brinkman term, with a Nernst-Planck tertiary current distribution electrode model using Butler-Volmer kinetics. The analytical model accurately predicts the optimum value of the volumetric surface area for both flow-through and interdigitated flow fields. A model taking into account the mass transfer inside the pores of the electrode is also considered, to extend the obtained expression for the optimum volumetric surface area. This model will be useful for cases in which the current density is close to the limiting value. The expressions obtained in this work give direct insight into what parameters are of most influence and allow designing flow cells that are optimal for the intended operating conditions. In comparing with several experiments from the literature we find that typically the channel thicknesses and pores sizes in experiments are chosen an order of magnitude too large, compared to what is optimal from an energy perspective. Therefore, we argue that more attention should be paid to operating cells with higher pressure drops, in order to lower overall energy use.

Declaration of competing interest
The authors declare the following financial interests/personal relationships which may be considered as potential competing interests: Abhiroop Bhadra reports financial support was provided by European Union.

Data availability
No data was used for the research described in the article.

Acknowledgement
The authors would like to acknowledge the funding support from the European Union's Horizon 2020 research and innovation programme under grant agreement No 875524.

Appendix A. Power losses
The power loss due to ohmic dissipation across the electrolyte channel is given by (A.1) Assuming uniform ionic conductivity and a one-dimensional flow of current in the -direction and ⟨ 2 ⟩ ≈ ⟨ ⟩ 2 = ⟨ ⟩ 2 , we obtain Eq. (13). Note that we use the symbol here for the -dependent magnitude of the current density in the current collector, not for the local magnitude of the vector which depends also on . Correspondingly, the averaging here is over .
If the electrode effectiveness factor [27] is close to 1, the local current density ⟂ can be assumed to be constant in the -direction. Integrating Eq. (9) in this case gives (A. 2) The power dissipated due to activation losses in the porous electrode is given by Assuming ⟂ does not depend on , using Eqs. (9) and (A.2), gives Eq. (14).
The power dissipated due to frictional losses in the electrolyte channel/porous electrode is given by .
Using the incompressibility condition, Eq. (3), and the fundamental theorem of calculus Here, is the superficial velocity along the length of the channel or porous electrode. If the pressure is considered constant along , its difference between the inlet and outlet of a channel/porous electrode is denoted by , and the average velocity in the channel/porous electrode as ⟨ ⟩ = 1 ch/pe ∫ ch/pe 0 ; the pumping loss fr can be written as Combining Eqs. (15) and (A.6), we obtain Eq. (17). Taking the symmetric case for the definition of the overpotential from Eq. (11), using Eq. (A.2), Eq. (A.2), and previously defined quantities in Eqs. (14) and (20) we find = asinh A. Bhadra and J.W. Haverkort which gives also the local dissipation per unit area at a given . Taking as the local concentration of the reactant entering either through the channel or the porous electrode, we define c as the cupmixing average concentration of the reactant in the channel/porous electrode of width ch/pe : By using Faraday's law and equating the difference between the rate at which reactants are entering and leaving the cell to the rate of charge leaving the cell, we have Here is the number of electrons transferred per reactant molecule. Using Eq. (18) we obtain from this expression Eq. (19).

Appendix B. Parameters
See Tables B.5 and B.6.

Appendix C. Optimal electrolyte channel width for a given boundary layer thickness
In membraneless systems, a function of the laminar flow channel is to avoid product cross-over. As products are formed and enter the electrolyte channel they continue diffusing in the transverse direction while also being advected in the streamwise direction. The result is a boundary layer thickness, defined by the distance at which the concentration has decreased by 99% compared to its wall value, given by where is the product species diffusion coefficient. This expression assumes a linearized flow profile ≈ 6⟨ ⟩ ∕ ch near the wall of a Poiseuille flow (Lévêque approximation) and constant wall concentration (Dirichlet) boundary conditions. See e.g. Ref. [54] for the analytical solution. Only the pre-factor changes slightly in case of constant flux (Neumann) boundary conditions.
In case two such boundary layers arise from products with similar diffusion coefficient, from the top and bottom of the channel, and we require that their 99% 'tails' do not overlap at the end of the channel, we require ch ≥ 3.2 When this is marginally satisfied it gives Parameters to find the optimal volumetric surface area of a porous electrode in both the flow-through and interdigitated simulations. Here ns , ns , and ns are the width, porosity, and pore size of the nanoporous separator in between the porous electrode and channel, where the density and viscosity are taken to be that of the electrolyte.̇is the rate of mass flow in both the electrolyte channel and porous electrode.  (25), where √ 3.2 3 ∕0.82 ≈ 6.3. The optimal velocity and gap thickness both increase with channel length ℎ. A higher current density increases the optimal velocity, but decreases the optimal gap thickness.
Note that the energy efficiency of a pump can be easily included in the definition of the frictional power losses by dividing the viscosity with the pump energy efficiency. This increases the power losses associated with friction and will result in a slightly higher optimal gap thickness.

Appendix D. Optimal electrode thickness
In Ref. [27] an expression is derived for the thickness of a porous electrode that minimizes the activation overpotential. It was assumed that the ionic current enters the electrode from the direction of the counter-electrode and sees a constant effective conductivity while the electronic current enters from the opposite side, facing the current collector, and sees a constant effective conductivity . The optimal electrode thickness is given by For the high conductivity typical of most metals, this optimum can become very large. This happens because the electronic ohmic drop is very small so that the electrode can be made very thick in order to benefit from a small decrease in the activation overpotential. Since cost and practical arguments usually also play a role, a more pragmatic optimum is that beyond which the decrease in activation overpotential with increasing electrode thickness becomes small. This leads to the expression given by [27] pe,opt ≈ In this expression, contrary to Eq. (D.1), the electronic conductivity no longer plays a role when it is much larger than the ionic conductivity. Note that either way, the optimal porous electrode thickness is inversely proportional to the current density so that the higher current densities call for thinner porous electrodes. Inserting these expressions into Eq. In case of linear kinetics the optimal electrode thickness reads [27] pe, opt = Here, opt ≈ 2 when and are of similar magnitude, but slightly higher if they differ by orders of magnitude [27]. In this case the optimal thickness depends on the volumetric surface area, so that we cannot simply insert it after optimizing . Since now the relation between current, velocity, and conversion, Eq. (19), depends on we have to choose which of these to maintain constant during the optimization. For linear kinetics asinh ( where in case of a constant current density ∕⟨ ⟩ = 1.