The wall-jet electrode

The hydrodynamics of the wall-jet electrode (WJE) are considerably more complex than in the case of the channel or rotating disc electrodes. Fortunately most experiments are conducted at one of two limits. The first is formally referred to as a 'wall-jet' where the jet is very small compared to the electrode radius. Analytical expressions (given in section 8.1.1) have been derived for the velocity profile near the electrode surface in this limit. The second is referred to as a 'wall-tube' or 'micro-jet' where the electrode radius is much smaller than the jet radius. Analytical expressions (given in section 8.2), based on two empirically determined parameters, have been derived for the velocity profile close to the electrode surface. For an intermediate case, no analytical solutions appear to have been attempted, though such a system could be simulated numerically using computational fluid dynamics (usually finite element) software to solve the Navier-Stokes equations.

As will be shown in the conclusions chapter, both of the limiting geometries offer experimental advantages - the former spans a very wide range of timescales and the latter can reach very high mass-transport rates. Unfortunately, the lack of theory available has limited their experimental adoption. This chapter attempts to redress this situation, acknowledging that some modelling has already been done for the wall-jet using the BI method, ignoring radial diffusion (see below). The importance of radial diffusion is therefore assessed and efficient modelling methods are introduced utilising conformal mappings. Finally, working curves and surfaces are computed for common electrochemical mechanisms allowing mechanistic analysis of steady-state voltammetric data.

Previous attempts^{1} at a two-dimensional WJE simulation incorporating radial diffusion, using a rectangular mesh in cylindrical polar co-ordinates, failed to give acceptable convergence due to the highly non-uniform concentration distribution in the vicinity of the electrode.

In order to overcome this problem, the WJE simulations by Compton et al.^{2},
relied on the frontal nature of the BI method (section 2.3)
to incorporate an expanding grid via an interpolation technique. The frontal
nature of the simulation, however, meant that radial diffusion could not be
included. Given the decrease in radial velocity towards the disc edge, one might
expect this to be significant and this was postulated as a possible explanation
for discrepancies between simulations and experimental data at low flow rates
(a simulation study was also conducted using the BI method using a displaced
radial diffusion term)^{3}.

Figure 8.1: Expanding grid used by Compton et al., based on interpolation.

Figure 8.2: Wall-jet velocity profile.

A more flexible and accurate method to account for the non-uniformity than
interpolating between frontal solution vectors is to perform the simulations
in a conformal space based on the hydrodynamics. This was done by Laevers et
al.^{4}, using a frontal application of the Crank-Nicolson
method (see section 2.1.3) to simulate the
chronoamperometric response without radial diffusion. The transformation is
generalised in this work and applied to the full mass transport equation (including
radial diffusion). The following hydrodynamic equations form the basis of the
transformation.

Glauert^{5} studied the velocity distribution in a wall-jet, shown in Figure 8.2, and found an exact solution of the boundary layer equation. The radial component of the convection is given by^{2,5}:

(8.1) |

and the normal component by:

(8.2) |

where the curvilinear co-ordinate h, is related to r and z by:

where , and and n=5/4. | (8.3) |

The functions f'(h) and h(h) are related though the following equations:

(8.4) |

(8.5) |

(8.6) |

in terms of the quantity g, a function of h, defined implicitly by:

(8.7) |

Figure 8.3: Functional dependence of h on g. The 1

This is plotted in Figure 8.3 together with one and three term series expansions of g:

(8.8) |

where a = 3.086x10^{-3} and b = 4.8991x10^{-5}. Thus 2-term expressions for v_{r} and v_{z} may be found:

and | (8.9) |

which, as h(r)0, collapse to the 1-term expressions:

and | (8.10) |

The one and two term expressions are plotted in Figure 8.4. Compton et al.^{2} found that the diffusion layer was within h=0.4, which would imply that a 1^{st}-order approximation of the convection would be adequate. It was found in the simulations reported here that the difference in the current, using 1^{st} or 2^{nd} order approximations of the convection, was indistinguishable (to 6 or 7 s.f.) over the range of Peclet numbers studied.

In order to formulate a simulation method (which incorporates radial diffusion effects) in a conformal space similar to that of Laevers et al., we seek a curvilinear transformation:

a = a(r,z) (r) a = a(r,y) where y = y(r,z) | (8.11) |

As shown in Appendix 3, under this transformation, the mass-transport equation becomes:

(8.12) |

(symbols are defined in the Appendix). The mixed second derivatives present an added complication.

Figure 8.5: Stencils on which the finite difference representation of mixed second derivatives may be formulated. On a 9=point stencil the four corners may be used. On a 7 point stencil, the average may be taken between two sets of four points.

These may be represented as finite differences at the corners of a 9-point stencil:

(8.13) |

or the average of the top-left and bottom-right squares of a 7-point stencil:

(8.14) |

Note that although the two partial derivatives are not necessarily equal, they reduce to the same finite difference formula. Introducing the coefficients (for brevity):

(8.15) |

The finite difference form of the mass-transport equation is therefore:

(8.16) |

This may be represented on a 7-point stencil (such as that for MGD1, introduced in section 4.1.1) as:

A^{1}_{j,k} | (k_{yd} - k_{ys} + k_{x} - k_{vy}) u_{j-1,k} |

A^{2}_{j,k} | (-k_{x}) u_{j-1,k+1} |

A^{3}_{j,k} | (k_{rd }- k_{rs}) u_{j,k-1} |

A^{4}_{j,k} | (-2k_{yd} - 2k_{rd} - 2k_{x}) u_{j,k} |

A^{5}_{j,k} | (k_{rd }+ k_{rs}) u_{j,k+1} |

A^{6}_{j,k} | (-k_{x}) u_{j+1,k-1} |

A^{7}_{j,k} | (k_{yd} + k_{ys} + k_{x} + k_{vy}) u_{j+1,k} |

Since the curvilinear transformation physically corresponds to a stretching of the z-co-ordinate moving along the electrode, the simulation space (depicted in Figure 8.6) and the boundary conditions (shown in Table 8.1) remain as they would be in cylindrical polar co-ordinates.

Figure 8.6: Simulation space for the WJE.

Table 8.1: Boundary conditions for WJE simulations in conformal space.

Zone | Boundary condition |

r = 0, all y | (symmetry) |

r = r_{max}, all y | (waste fluid) |

y = 0, r < r_{e} | a = 0 (electrode) |

y = 0, r > r_{e} | (insulator) |

y = y_{max}, all r | a = 1 (bulk solution) |

The implementation of the boundary conditions, shown in Table 8.2, is complicated by the cross-derivatives as the 'extra' boundary nodes (2 and 6) must reinforce the boundary conditions. Note that at r=0 there is a singularity in the Laplacian operator (since 1/r = infinity), but since the simulation space stops one node short of the boundary, this problem is avoided (if the boundary were to be simulated, a Maclaurin expansion could be used to avoid the problem).

Disc centre | k=1, j=1...N_{}GY | A^{}>4_{j,k} = A^{}>4_{j},k + A^{3}_{j,k}; A^{3}_{j,k} = 0 |

k=1, j=1...N_{GY} | A^{}>7_{j,k} = A^{}>7_{j},k + A^{6}_{j,k}; A^{6}_{j,k} = 0 | |

Waste | k=N_{GX}, j=1...N_{}GY | A^{}>4_{j,k} = A^{4}_{j,k} + A^{}>5_{j},k; A^{5}_{j,k} = 0 |

k=N_{GX}, j=2...N_{}GY | A^{1}_{j,k} = A^{1}_{j,k} + A^{2}_{j},k; A^{2}_{j,k} = 0 | |

Bulk solution | j=N_{GY}, k=1...N_{}GX | A^{7}_{j,k} = 1 |

j=N_{GY}, k=2...N_{}GX | A^{}>6_{j,k} = 1 | |

Electrode | j=1, k=1...k_{E} | A^{1}_{j,k} = 0 |

j=1, k=2...k_{E-1} | A^{}>2_{j,k} = 0 | |

Insulator | j=1, k=k_{E}...N_{}GX | A^{4}_{j,k} = A^{4}_{j,k} + A^{1}_{j,k}; A^{1}_{j,k} = 0 |

j=1, k=k_{E}-1...N_{GX} | A^{5}_{j,k} = A^{5}_{j,k} + A^{2}_{j,k}; A^{2}_{j,k} = 0 |

Figure 8.7: Curvilinear 'boundary layer' co-ordinate

Figure 8.8: Convergence plot comparing trapezoidal and rectangular integration.

So far, all the formulation has been for a general function, y. This approach has the advantage that functions can readily be interchanged and thus optimised without reworking the algebra. We now introduce a specific function for y which allows the full 2-D steady-state problem to be solved efficiently.

As Leavers et al. suggested, the hydrodynamic co-ordinate, h, (plotted in Figure 8.7) may be used as the basis of the transformation:

(8.17) |

(where n=5/4). The partial derivatives of this variable are

| (8.18) |

These may be substituted directly into the finite-difference equation (8.16). The current may be evaluated from:

(8.19) |

or as a first-order finite difference form using rectangular quadrature:

(8.20) |

As in the case of the channel electrode (see section 3.9), significant error cancellation occurs if rectangles are used to perform the integration rather than trapezia. This is illustrated in Figure 8.8.

Figure 8.9: Slice through Figure 8.11, at one node above the electrode surface. |
Figure 8.10: Slice through Figure 8.11 in the normal co-ordinate at the disc centre. |

Figure 8.11: Simulated concentration distribution in (r,h) space at log_{10 }P_{s} = 5.45. Parameters are: N_{GX}=121, N_{GY}=121, l_{d}=1, h_{max}=1, n=0.0089cm^{2}s^{-1}, D=2x10^{-5}cm^{2}s^{-1}, [A]_{bulk}=1x10^{-6 }molcm^{-3}, v_{f}=0.05cm^{3}s^{-1}, r_{e}=0.2cm, r_{jet}=0.0345cm, 2^{nd} order convection, upwind differences in r, central in h.

At high Peclet numbers, the concentration distribution in the conformal space resembles a rotating disc electrode (RDE). The flux is uniform across the electrode surface, shown in Figure 8.9, until the electrode edge where it rises sharply. The flux in the normal co-ordinate, shown in Figure 8.10, is linear near the electrode surface, then curving towards its bulk value away from the electrode. At lower Peclet numbers, radial diffusion becomes more significant and the concentration distribution tends to that of a microdisc. This suggests that secondary transformations may be applied based on those used to improve the efficiency of microdisc or RDE simulations.

One possibility would be to use a second conformal transformation such as Amatore and Fosset's^{6} to transform the (r,h) space into a closed form. There are a number of disadvantages with this, however:

- At high Peclet numbers the concentration distribution more closely resembles the RDE than a microdisc.
- The conditioning of the linear system in (r,h)
space is already rather poor - stability of the MGD1 solver was a constant
problem throughout the simulations presented in this chapter. It is highly
likely that the equation system would be too ill-conditioned even for the
Krylov subspace solvers (introduced in section 4.2)
in the doubly transformed space, as was found when attempting to simulate
kinetics at a microdisc electrode in Verbrugge and Baker's transformed space
^{7}(see section 9.2.2). - The mathematics would be extremely tedious since the equations defining Amatore and Fosset's transformation cannot be inverted algebraically - the partial derivatives must be solved for simultaneously.

Therefore, based on the work of Taylor et al.^{8}, expanding grids were used to refine the mesh in the diffusion layer and around the boundary discontinuity at the electrode edge.

Figure 8.12: Sigmoidal co-ordinate transformation to concentrate finite difference nodes at the edge of a disc. Several values of g are shown.

In the case of the microdisc electrode, the singularity at the electrode edge is the major impediment to convergence. The effect is less pronounced in the wall-jet as the radial convection increases with the distance from the jet centre - this effectively 'smears out' the concentration gradient at the electrode edge. Even so, the concentration profile is still steepest at the electrode edge and also (due to the 2pr integration factor) this is the zone of the electrode that contributes the largest fraction of the current, hence any error in this zone will be amplified. This suggests, particularly at low Peclet numbers, that the simulation accuracy could be improved by transforming the radial co-ordinate via a sigmoidal function to concentrate nodes at the electrode edge, as shown in Figure 8.12. Taylor et al.^{8} used a Fermi-Dirac function (sigmoid) for this purpose:

(8.21) |

As with the channel electrode (Chapter 6) a dimensionless expansion parameter, g_{d}:

g_{d} = r_{e}/g | (8.22) |

was used in the simulations, which is independent of the electrode radius.

In order to perform the transformation (which is presented in Appendix 3):

a(r,y) => a(r,y) where r=r(r) | (8.23) |

we require the first and second derivatives of the transformed variable with respect to r, since the mass-transport equation is now given by:

(8.24) |

The first derivative is:

(8.25) |

and the second derivative:

(8.26) |

A simulated concentration profile is shown in Figure 8.13, showing the contraction of the mesh at the electrode edge resulting in better definition of the rapidly changing concentration gradient in this region.

Figure 8.13: Simulated concentration profiles in (r,h) space and (r,h) space with an expansion parameter of 10. The concentrations at all h values are shown on the plot. Parameters are: N

Whilst conducting simulations in this space, it became apparent that at low Peclet numbers the large values of R (dimensionless) required to accommodate radial diffusion were resulting in values of r very close to 1 which were being rounded up due to finite machine precision. As was found with ChE simulations in section 6.3, this resulted in a singular system of equations.

If e is the machine precision, the point where the co-ordinate transformation breaks down is:

1-r_{max} < e | (8.27) |

where r_{max} is the value at the downstream boundary, given by:

(8.28) |

Thus the co-ordinate transformation may be used only when:

l_{d}.g_{d} < ln(e^{-1}-1) » ln(e^{-1}) | (8.29) |

Since l_{d} is fixed by the Peclet number, this effectively imposes an upper limit on the amount of grid expansion which decreases as the amount of radial diffusion increases:

g_{d} < ln(e^{-1})/l_{d} | (8.30) |

The value of ln(e^{-1})
was found to be approximately 86 using 8 byte double precision floating-point
variables on an SGI Indigo^{2 }(tm). This is inconvenient,
although it is partially offset by the concentration profile being more spread
out at low Peclet numbers, thus a lower grid expansion factor is required. It
should be noted that similar problems should occur with any other transformation
function which is used to expand the grid nodes 'exponentially' towards infinity,
but is fixed at a finite point (other than the origin) such as the electrode
edge. This is one argument for adopting the approach of Pastore et al.^{
9}, fixing the grid at infinity and allowing the electrode size
to vary with the mesh size. The main disadvantage of such a method is that a
transformation to an equivalent dimensionless problem is required to simulate
a given electrode length (for example in the ChE simulations in section 6.3
this was achieved by compensating the cell height). This is especially important
if one wishes to use 2^{nd} or higher-order convection
approximations in WJE simulations where the current is no longer a sole function
of the Peclet number.

A second disadvantage with this transformation (both for the WJE and microdisc) becomes apparent when one considers a simulation involving homogeneous kinetics, such as an ECE process. Although a sigmoidal transformation might improve the definition of species A, the concentration distribution of kinetically unstable species falls away rapidly along the radial co-ordinate, thus a transformation which shifts points from the centre to the edge of the disc would impede the definition of these species. Since Taylor et al.^{8} only considered heterogeneous kinetics, they did not encounter this problem. In Chapter 9, the conformal mapping of Amatore and Fosset^{6} is shown to be well-suited to the simulation of homogeneous processes at a microdisc electrode.

Returning to the WJE, the major obstacle to convergence of the (r,h) simulations was the accurate description of the concentration gradient near the electrode surface, hence rather high values of N_{GY} were required for satisfactory accuracy. An exponentially expanding grid in the normal co-ordinate mimics the wall-jet concentration profile in (r,h) space reasonably well and has been successfully applied to the RDE^{10}. This may be incorporated into the generic WJE transformation by transforming h through Feldberg's ln(1+ax) expanding grid function:

where and | (8.31) |

The first derivatives are:

(8.32) |

and

(8.33) |

The second derivatives are:

(8.34) |

and

(8.35) |

Again, these can be substituted directly into the finite difference form of the transformed mass transport equation. A simulated concentration profile is shown in Figure 8.14. The transformation packs more nodes close to the electrode surface, so that the curvature of the diffusion layer is reduced in the transformed co-ordinate.

Figure 8.14: Concentration profile in (r,y) space with an expansion parameter of 10. Parameters are: N

Figure 8.15: Convergence plot for (r,h) mapping and the (r,y) mapping with roughly optimal and greater than optimal expansion parameters. Of course as the expansion parameter tends to zero, the mapping tends to the (r,h) case.

Optimisation of the expansion parameter is made awkward by the fact that the
current continually drops with increasing grid expansion, thus unlike the case
with an exponentially expanding normal co-ordinate in the ChE (see section 6.1),
there is no turning point to signify the optimal value. For this example system,
by examining the rate of convergence in the r-co-ordinate, optimal values of
the expansion parameter were found to be around 2.5 (using h_{max}=1),
as shown in Figure 8.15, and since h_{max}
does not need to be varied, this value will suffice for most simulations. The
rate of convergence is considerably accelerated - the number of nodes required
is reduced by a factor of 6 in this case. The main advantage of this transformation,
however, is when fast kinetics are simulated, since it concentrates nodes in
the reaction layer. Of course the exponential function only roughly approximates
the shape of the concentration distribution (consider the linearity near the
electrode surface in Figure 8.10). It might be possible to construct a more
efficient (at least for an electron transfer) closed-space transformation for
the normal co-ordinate based on the Hale transformation for the RDE (see section
7.2).

It is shown in Appendix 4 that under a 1^{st} order convection approximation, the dimensionless current is a sole function of the shear-rate Peclet number, P_{s}. In this section WJE simulations in conformal (r,h) space are used to quantify this relationship. At high Peclet numbers, this should limit to Levich behaviour. The departure from Levich behaviour at lower Peclet numbers allows a quantitative assessment of the importance of radial diffusion. Simulations were run for the parameters shown in Table 8.3.

Kinematic viscosity | n | 8.9x10^{-3} cm^{2}s^{-1} |

Diffusion coefficient | D | 3.2x10^{-6} cm^{2}s^{-1} |

Electrode Radius | r_{e} | 4.025x10^{-1} cm |

Nozzle radius | r_{jet} | 3.45x10^{-2} cm |

Bulk concentration | [A]_{bulk} | 1x10^{-6} mol cm^{-3} |

Space downstream | l_{d} | 2 |

number of nodes in R | N_{GX} | 641 |

Number of nodes in h | N_{GY} | 641 |

Figure 8.16: Working curve for a simple electron transfer showing the
departure from Levich behaviour at low Peclet numbers due to radial diffusion.
The response of the channel electrode (transformed through an equivalent diffusion
layer thickness - see section 10.3) is also
shown for comparison.

The simulation data is summarised in Figure 8.16, a log-log working curve of dimensionless current vs. Peclet number. The linear Levich region is evident at high Peclet numbers as is the departure from this behaviour due to radial diffusion effects at low Peclet numbers. From the simulation data one may conclude that radial diffusion effects become significant (i.e. augment the current by >1%) below log_{10}Ps = 3.5. The experimentally accessible^{11} range of Peclet numbers span from ca. 1.4 (where radial diffusion augments the current by 11%) up to 9.0 (where the Levich equation is obeyed). Note that the amount of axial diffusion is greater at the equivalent channel microband than radial diffusion at the wall-jet.

Since radial diffusion effects are negligible over most of the experimentally-accessible
range of Peclet numbers, a working curve of the experimental observable vs.
dimensionless rate constant can be used for the analysis of most kinetic data.
This also allows a solver such as the BI or frontal PKS method (see section
4.5) to be used to conduct these simulations,
saving substantially in CPU time. One problem with using a frontal solver is
the choice of the initial concentration vector corresponding to the jet impinging
on the electrode surface. Both Laevers et al.^{4}
and Compton et al.^{2} used the bulk concentration
in their BI simulations, but Figure 8.10 shows that the concentration is uniform
in h across the disc surface at steady-state. This
means that in a steady-state simulation, if the bulk concentration is used for
the initial concentration vector, the current is artificially augmented and
a very high number of nodes in the x co-ordinate are required for accurate simulations.
This is illustrated in Figure 8.17. Since the simulations of Laevers et al.^{4}
were time-dependent (chronoamperometry), the error from the incoming concentration
boundary would increase as the current approached steady-state.

In such a simulation, the concentrations relax to the true value by the end of the electrode, and this should also, in fact, be the same concentration vector as at the disc centre. Therefore the incoming concentration vector may be 'pre-relaxed' before the simulation commences.

Figure 8.17: Electrode surface concentrations. Errors are evident in the BI simulations due to the jet boundary condition. The 'true' value was simulated using N

Figure 8.18: Convergence plot of frontal WJE simulations showing the importance of pre-relaxing the incoming concentration vector. 20 cycles of the Thomas Algorithm (on the mass-transport coefficients at k=1) were used to relax the initial concentration vector.

This can be done by sweeping the concentration vector of the incoming species until it reaches a steady value (since only 10-20 iterations are required, this is computationally inexpensive). The actual BI or Frontal PKS simulation may then be conducted (adding any kinetic terms) using this vector as the starting approximation. This dramatically improves the convergence of the simulations, as shown in Figure 8.18.

Working curves were simulated using the BI and frontal PKS methods for EC, EC_{2}, ECE, EC_{2}E, DISP1 and DISP2 mechanisms using the parameters shown in Table 8.4 (though since the working curves are dimensionless, many of these are arbitrary). An exponentially expanding grid in the y co-ordinate with an expansion parameter of 2.5 was used to improve accuracy. A uniform mesh spacing was used in the radial co-ordinate.

Number of nodes in r | N_{GX} | 1000^{†} |

Number of nodes in y | N_{GY} | 1000* |

Diffusion coefficient | D | 1x10^{-5} cm^{2}s^{-1} |

Kinematic viscosity | n | 0.0089 cm^{2}s^{-1} |

Electrode radius | r_{e} | 0.2 cm |

Nozzle radius | r_{jet} | 0.0345 cm |

Volume flow rate | v_{f} | 0.01cm^{3}s^{-1} |

Bulk concentration | [A]_{bulk} | 1x10^{-6} mol cm^{-3} |

(8.36) |

are summarised in Figure 8.19. At lower Peclet numbers, a working surface
(of N_{eff} vs. log K vs. log P_{s})
is necessary for the analysis of kinetic data. Surfaces were simulated for ECE
and DISP1 mechanisms using ILU(1) preconditioned, multigrid accelerated BiCGStab(4)
(introduced in section 4.3) on a mesh of 513x513
nodes.

The simulations were run in strips of fixed P_{s}, increasing the rate constant to traverse a working curve in log K. The concentration distribution at the previous rate constant was used as the starting approximation for the next, considerably reducing the number of V-cycles after the first rate constant was simulated. The tolerance had to be tightened to 1x10^{-14} to ensure satisfactory convergence (possibly due to the close starting approximation). Again, an exponentially expanding grid in the y co-ordinate with an expansion parameter of 2.5 was used to improve accuracy. Figure 8.20 shows how, moving across each working surface, the working curves change with Peclet number. Comparing the response at Ps = 10 with that in the limit of no radial diffusion, the working curve is distorted at most by 2.5%. Comparing this with the augmentation in current for a simple E process - approximately 18% at Ps = 10, shows that the radial diffusion effect is largely cancelled out in N_{eff}.

Homann^{12} and then Frössling^{13} derived an asymptotic solution to Navier-Stokes equations, valid close to the electrode surface as a series expansion:

(8.37) |

and

(8.38) |

where a is an undetermined flow rate-dependent parameter with units of s^{-1},

(8.39) |

and the coefficients (according to Frössling) are:

a_{2} = 0.656, a_{3} = -0.16667, a_{6} = 0.0036444 and a_{7} = 0.00039682 | (8.40) |

Albery and Bruckenstein^{14} proposed use of the just the first term in the series expansion:

and | (8.41) |

where

(8.42) |

Chin and Tsang^{15} found empirically that for 0.1 < r_{jet}/z_{jet} < 2.5:

(8.43) |

and this was also confirmed experimentally by Macpherson et al.^{16,17}. Note however that Chin and Tsang's determination was based on only four experimental data points. In order to assess the validity of Frössling's boundary layer approximation and Chin & Tsang's empirical result for the hydrodynamics independently, Coles^{18} is currently simulating the MJE numerically using the commercial fluid dynamics package, FIDAP^{19}. For now, it will be assumed that the velocity field predicted by the above equations is a reasonable approximation within the diffusion layer and use this as a basis to investigate over what range of Peclet numbers the MJE is uniformly accessible.

Since the wall-tube geometry is at least approximately uniformly accessible, there is little to be gained from using a hydrodynamic conformal transformation as was applied to the wall-jet. Therefore simulations were conducted on a simple rectangular mesh in cylindrical polar co-ordinates. The simulation space and boundary conditions are exactly as in Figure 8.6 and Table 8.1 for the WJE except that the normal co-ordinate is now z rather than y. The central finite difference form of the steady-state mass transport equation is:

(8.44) |

which has no cross derivatives hence the finite difference implementation boundary conditions, shown in Table 8.5, is simpler than in the WJE. Since the mass transport equation may be represented by a 5-point stencil in this case, a wider range of solvers could be applied. It was found in practice that MGD1 was unstable at both low and high Peclet numbers even if diagonal dominance of the coefficient matrix was maintained by upwind convection terms, so Stone's SIP (see section 2.6.2.4) was used as the solver for these simulations.

A simulated concentration profile is shown in Figure 8.21, clearly conforming to the hypothesis of uniform accessibility - above the electrode surface the concentration contours are linear.

Disc centre | k=1, j=1...N_{GY} | A^{4}_{j,k} = A^{4}_{j,k} + A^{3}_{j,k}; A^{3}_{j,k }= 0 |

Waste | k=N_{GX}, j=1...N_{GY} | A^{4}_{j,k} = A^{4}_{j,k} + A^{5}_{j,k}; A^{5}_{j,k }= 0 |

Bulk solution | j=N_{GY}, k=1...N_{GX} | A^{7}_{j,k} = 1 |

Electrode | j=1, k=1...k_{E} | A^{1}_{j,k} = 0 |

Insulator | j=1, k=K_{E}...N_{GX} | A^{4}_{j,k} = A^{4}_{j,k} + A^{1}_{j,k}; A^{1}_{j,k }= 0 |

The current may be evaluated from simple rectangular quadrature:

(8.45) |

Given the uniformity of the concentration distribution over the electrode, there is little to be gained from higher order integration formulae.

Figure 8.21: Simulated concentration profile at a MJE including radial diffusion effects. Parameters are as in Table 8.6 with v

It is shown in Appendix 5 that, as for the WJE, RDE and ChE, the dimensionless current is a sole function of the shear-rate Peclet number under a 1^{st} order convection approximation. This means that it is possible to quantify the effects of radial diffusion by examining the departure from Levich behaviour at low Peclet numbers. The simulation parameters are shown in Table 8.6.

Diffusion Coefficient | D | 1x10^{}-5 cm^{}>2s^{}-1 |

Kinematic viscosity | n | 8.9x10^{}-3_{}> cm^{2}s^{}-1 |

Electrode Radius | r_{e} | 25 mm |

Nozzle Radius | r_{jet} | 52 mm |

Electrode-nozzle separation | z_{jet} | 0.3 cm |

Bulk concentration | [A]_{bulk} | 1x10^{}-6 molcm^{}-3 |

Space downstream | l_{d} | 1.5 |

Number of nodes in r | N_{GX} | 600 |

Number of nodes in z | N_{GY} | 600 |

Figure 8.22: Levich plot for MJE.

The Levich plot shown in Figure 8.22, allows the radial diffusion effect to be assessed. The geometries and flow rates used by Macpherson and Unwin^{16} spanned Peclet numbers in the range 3.5 < log_{10}Ps < 7. At the lowest, radial diffusion augments the current by approximately 2.5%, and the augmentation drops below 1% at log_{10} Peclet numbers above 4.5.

It may be desirable to increase the rate of flow in order to achieve higher
rates of mass transport, which should be reasonably easy given the relatively
modest volume flow rates employed by Macpherson and Unwin (a maximum of 5x10^{-2}
cm^{3}s^{-1} was used
in ref. 16). The unimportance of radial diffusion, together with the assumption
that the dominating v_{z} convection term is quadratic
in z within the diffusion layer, means that kinetic analysis may be conducted
using theory generated for the RDE^{14}. This is
discussed in Chapter 10. Until the validity
of the hydrodynamic approximations has been confirmed, however, caution is advised.

^{1}J. Alden, Part II Thesis, Oxford University, (1995).^{2}R.G. Compton, C.R. Greaves, A.M. Waller,*J. Appl. Electrochem*.,**20**, (1990), 575.^{3}R.G. Compton, A.C. Fisher, M.H. Latham, R.G. Wellington, C.M.A. Brett, A.M.C.F. Oliveira Brett,*J. Appl. Electrochem.*,**23**, (1993), 98.^{4}P. Laevers, A. Hubin, H. Terryn, J. Vereecken,*J. Appl. Electrochem.*,**25**, (1995), 1023.^{5}M.B. Glauert,*J. Fluid. Mech.*,**1**, (1956), 625.^{6}C.A Amatore, B. Fosset,*J. Electroanal. Chem.***328**,^{7}M.W. Verbrugge, D.R. Baker,*J. Phys. Chem*.**96**, (1992), 4572.^{8}G. Taylor, H.H. Girault, J. McAleer,*J. Electroanal. Chem.*,**293**, (1990), 19.^{9}L. Pastore, F. Magno, C.A. Amatore,*J. Electroanal. Chem*.,**301**, (1991), 1.^{10}S.W. Feldberg, C.I. Goldstein and M. Rudolph,*J. Electroanal. Chem.*,**413**, (1996), 25.^{11}These values (which are listed in the conclusions chapter) are calculated from values given by C. Brett (e-mail brett@cygnus.ci.uc.pt) in a personal communication from the experimental parameters given in R.G. Compton, A.C. Fisher, G.P. Tyley,*J. Appl. Electrochem.*,**21**, (1991), 295.^{12}F. Homann,*ZAMN*,**16**, (1936), 153;*Forsschg. IngWes*,**7**, (1936), 1 [both in German]^{13}N. Frössling,*Lunds Univ. Avd.*,**2**, (1940), 35. [in German]^{14}W.J. Albery, S. Bruckenstein,*J. Electroanal. Chem.*,**144**, (1983), 105.^{15}D.T. Chin, C.H. Tsang,*J. Electrochem. Soc.*,**125,**(1978), 1461.^{16}J.V.Macpherson, S.Marcar, P.R.Unwin,*Anal.Chem.*,**66**, (1994), 2175^{17}J.V.Macpherson, M.A. Beaston, P.R.Unwin,*J. Chem. Soc. Faraday. Trans.*,**91**, (1995), 899.^{18}B.A. Coles, Oxford University, unpublished work. For further information, e-mail bac@physchem.ox.ac.uk.^{19}FIDAP is a commercial CFD package (based on a finite element method) made by Fluent Incorporated (http://www.fdi.com).