## Abstract

We present closed-form solutions for high Schmidt number mass transfer in a hydrodynamically fully developed turbulent flow. Governing equations for the near- and far-field are developed for a large class of boundary conditions (BCs) for which the mass flux is a function of the concentration at the wall. We show that for this class of BCs, which includes nonlinear wall reactions, the mass transfer coefficient is independent of the BC and the Sherwood correlation is therefore universal. For Dirichlet, Neumann and Robin BCs, the far-field solutions are in good correspondence with the method of separating variables and near-field solutions are in good agreement with numerical simulations. However, in contrast with the far-field solutions, the Sherwood correlation in the near-field depends on the specific BC. As an example of nonlinear BCs, solutions for a second-order wall reaction are derived which are compared with numerical simulations and found to be in excellent agreement.

## 1. Introduction

The exchange of mass and/or heat of a turbulent flow with its bounding surface occurs in several areas of engineering. One classical example is the corrosion of pipe walls owing to the presence of a corrosive species in water, such as carbon dioxide (Sydberger & Lotz 1982, Sarin *et al.* 2004) or chlorine (Rossman *et al.* 1994; Clark & Haught 2005; Al-Jasser 2007). Another classical example is the heat transfer through conduit walls, as is relevant, for example heat exchangers (Bergmann & Fiebig 1999) and district heating (Webb & Kim 2005). In this paper, we will concentrate on turbulent mass transfer problems, although the results are equally applicable to heat transport problems where buoyancy effects are negligible.

A typical quantity of interest is the realized mass (heat) flux at the wall, represented by the mass (heat) transfer coefficient and in dimensionless form by the Sherwood number *Sh* (Nusselt number). They will depend on the Reynolds number *Re* and the Schmidt number *Sc* (Prandtl number), the former representing the ratio between inertial and viscous forces, and the latter the ratio between kinematic viscosity and molecular (thermal) diffusivity. These relations often are presented as power-laws *Sh*=*b*_{1}*Re*^{b2}*Sc*^{b3}, where *b*_{1},*b*_{2} and *b*_{3} are coefficients. Determination of these coefficients has been the objective of many experimental, numerical and theoretical investigations and a selection of *Sh* correlations for high *Sc* is presented in the electronic supplementary material. The power-law exponents depend critically on the near-wall behaviour of the eddy diffusivity and several theoretical models based on theory of turbulent boundary layers have been developed over the years (Kader & Yaglom 1972; Aravinth 2000).

The mass transport equation (see §3) is a linear partial differential equation with variable coefficients. Many analytical methods therefore make use of the method of separating variables (Sleicher *et al.* 1970; Notter & Sleicher 1971, 1972; Biswas *et al*. 1993; Weigand *et al.* 2001; Weigand 2004), and with great success: the classical power-law relationship *Sh*=0.016*Re*^{0.88}*Sc*^{0.33} proposed by Notter & Sleicher (1972) remains widely used today (Rossman *et al.* 1994). However, there are some drawbacks and limitations to this method. First, the determination of each of the eigenfunctions and eigenvalues has to be done numerically because of the variable coefficients in the problem. Second, the method of separating variables is applicable to linear boundary conditions (BCs) only; obtaining solutions for nonlinear BCs, as of importance, for example, biofilm growth (Munavalli & Mohan Kumar 2004; Noguera & Morgenorth 2004) is generally not possible.

In this paper, we will focus on the development of asymptotic solutions for turbulent mass transport at high *Sc* and *Re* for linear (notably Dirichlet, Neumann and Robin) and nonlinear BCs. Closed-form solutions for the near- and far-field will be presented, where the near-field is the region where the concentration boundary layer is developing and the mass transfer coefficient will be dependent on the streamwise direction. The far-field is defined as the region where the mass transfer coefficient has become constant.

There is a long history of asymptotic solutions for developing concentration boundary layers. Analytical solutions for a thermal boundary layer in a fully developed turbulent flow were developed by Linton & Sherwood (1950). The derivation involved neglecting the streamwise advective, diffusive transport and all turbulent transport, thereby essentially reducing the problem to that solved by Lévêque (1928) in the context of heat transfer in a laminar boundary layer. By assuming self-similarity, closed-form solutions for the boundary layer growth and concentration profile and heat transfer could be provided. Electrochemical experiments of high *Sc* and *Re* mass transfer in the entrance region (Shaw *et al.* 1963; Berger & Hau 1977) showed good agreement between the predictions from the asymptotic solutions and the measurements. Kestin & Persen (1962) developed asymptotic solutions for the heat transfer across a developing boundary layer over a flat plate, including viscous entrance effects. Several solutions were presented, including a transition from a laminar to a turbulent boundary layer and stepchanges in the wall temperature.

Far-field asymptotic solutions for high *Sc* solutes, however, seem to have escaped attention until recently. Garcia-Ybarra & Pinelli (2006) used matched asymptotic expansions to derive a closed-form solution for the concentration profile for a fixed concentration (Dirichlet) BC at high *Sc*. Sookhak Lari *et al.* (2010) independently arrived at the same closed-form solution by observing that the scalar flux was approximately constant across the concentration boundary layer in a study on the decay of residual chlorine in pipes as a result of a first-order reaction with the wall (Robin BC). Despite the different BCs, both studies report the same *Sh* correlation.

The aim of this paper is to generalize the work of Garcia-Ybarra & Pinelli (2006) and Sookhak Lari *et al.* (2010) to much more general BCs and to provide simple closed-form solutions for the mass-transfer, decay coefficients and concentration profiles. After a discussion on the appropriate velocity and turbulent diffusivity profiles (§2), the governing equations for the far-field are derived (§3). It will become apparent that the *Sh* correlation reported in Garcia-Ybarra & Pinelli (2006) and Sookhak Lari *et al.* (2010) is in fact representative for a very large class of BCs (§4). Closed-form solutions for Dirichlet, Neumann and Robin BCs are presented in §5*a*. We then use the Von Karman–Pohlhausen method to pose the governing equation for the near-field and present solutions in §6. As the solution method is not limited to linear BCs, closed-form solutions for a second-order wall reaction are presented (§§5*b* and 6). Concluding remarks are made in §7.

## 2. Near-wall profiles of velocity and eddy diffusivity

Consider a fully developed turbulent flow field. When *Sc*≫1, as is the case for many mass transfer problems (and heat transfer in e.g. heavy oils), the scalar diffuses much slower than momentum. The associated layer of thickness *δ*_{m} near the wall where molecular diffusion dominates over turbulent transport, hereafter referred to as the mass transfer boundary layer (MTBL) will then be entirely nested in the viscous sublayer (Pope 2000; Schlichting & Gersten 2000; Garcia-Ybarra & Pinelli 2006; Sookhak Lari *et al.* 2010). Outside the MTBL, the turbulence causes sufficient mixing to assume a uniform concentration. The two fundamental parameters governing the viscous wall region are the kinematic viscosity *ν* and the friction velocity , where *τ*_{w} is the wall shear stress and *ρ* is the fluid density. In the viscous sublayer, the turbulent stress is negligible, and therefore the average velocity profile is given by (Pope 2000; Schlichting & Gersten 2000; Bird *et al.* 2002)
2.1where *u*^{+}=*u*/*u*_{τ} and *y*^{+}=*y*/*δ*_{v}, *y* represents the distance from the wall and *δ*_{v}=*ν*/*u*_{τ} is the viscous length scale.

A second property of the viscous sublayer is that the turbulent momentum flux , and therefore the eddy-viscosity *ν*_{T}, has a cubic dependence on the wall distance
2.2which can be shown using Taylor expansions (Antonia & Kim 1991; Bird *et al.* 2002). The prefactor *b* has been approximated experimentally and numerically, and takes the value *b*≈9.5×10^{−4} (Bird *et al.* 2002). Using (2.2), the ratio of the turbulence diffusion coefficient *D*_{T} to the molecular diffusion coefficient *D* is given by
2.3where *Sc*=*ν*/*D* is the Schmidt number and *Sc*_{T}=*ν*_{T}/*D*_{T} the turbulence Schmidt number. The effective exponent for *D*_{T} is crucial for high *Sc* mass transfer as it influences the *Sh* correlation: *D*_{T}∼*y*^{m} implies that *Sh*∼*Sc*^{1/m}. Many laboratory experiments (Harriott & Hamilton 1965; Mizushina *et al.* 1971; Dawson & Trass 1972; Berger & Hau 1977; Zhao & Trass 1997) find that *Sh*∼*Sc*^{0.32−0.35}, thereby indirectly confirming (2.3). Shaw & Hanratty (1977) report a slightly lower *Sc* dependence *Sh*∼*Sc*^{0.29}, although it is not entirely clear what the cause is for the deviations between this and the other experiments. More information about *Sh* correlations for high *Sc* mass transfer, including the range of *Sc* and *Re* considered, can be found in the electronic supplementary material.

A Taylor expansion confirms that the turbulent scalar flux is indeed expected to vary as the cubic on wall distance for fixed concentration (Dirichlet) BCs (Antonia & Kim 1991; Bird *et al.* 2002; Garcia-Ybarra & Pinelli 2006). However, for flux (Neumann) BCs, a second-order dependence of the turbulent scalar flux (and therefore *D*_{T}) on the wall distance is obtained. As the *D*_{T} profile at high *Sc* has not been reported as yet for flux BCs, it is not known how dominant the second-order term is. The only available data are from simulations for the heat transfer across a fluid layer and a solid wall with finite thermal conductivity (Tiselj *et al.* 2004; Bergant & Tiselj 2007). The simulations show that *Sc*_{T} decreases very close to the wall, as expected for a quadratic *D*_{T} profile. However, Bergant & Tiselj (2007) report that the influence on the mean temperature profiles and the heat transfer coefficients is almost negligible.

For Dirichlet BCs, studies performed with Direct Numerical Simulation (DNS), show that *Sc*_{T} is indeed constant when *Sc* is of order unity (Antonia & Kim 1991; Kawamura *et al.* 1998; Schwertfirm & Manhart 2007). At *Sc*>10, an increase of *Sc*_{T} is observed very close to the wall (*y*^{+}<1), which becomes more pronounced for higher *Sc* (Na & Hanratty 2000; Crimaldi *et al.* 2006; Bergant & Tiselj 2007; Schwertfirm & Manhart 2007; Kozuka *et al.* 2009). Garcia-Ybarra (2009) used data from numerical simulations to show evidence that the fourth-order term overwhelms the cubic term a bit further away from the wall, which suggests that the effective exponent for *D*_{T} is larger than three. This is consistent with the experiments of Shaw & Hanratty (1977) but not with the other studies mentioned above. Further laboratory experiments and/or DNS at higher *Re*_{τ} and *Sc* will be required to settle this issue.

In what follows we will use the classical assumption (Kader 1981; Bird *et al.* 2002) that *D*_{T} is a cubic and that *Sc*_{T} is constant. Even though this excludes some of the phenomena described above, the net effect of *y*^{+} variation of *Sc*_{T} (the modification of the mass transfer coefficient, etc.) can be incorporated by tuning of the parameter *b*/*Sc*_{T} which is discussed in appendix A. The procedure maps the actual profile for *D*_{T} onto a cubic which has the same boundary layer thickness *δ*_{m} (defined below), thereby ensuring that integral quantities be predicted accurately. Note that if the profile for *D*_{T} differs significantly from a cubic, the parameter *b*/*Sc*_{T} will become dependent on *Sc* and *Re*, the consequences of which will be described in §7.

The properties of the viscous sublayer and MTBL for a high *Sc* solute are depicted in figure 1. Figure 1*a*,*b* shows the dependence of *u*^{+} and *D*_{T}/*D* on the wall distance, and figure 1*c* shows a typical concentration profile valid for e.g. a first-order reaction of chlorine with the wall (Robin BC) (Sookhak Lari *et al.* 2010). We follow Kader (1981) and define the typical thickness *δ*_{m} of the MTBL as the distance from the wall at which *D*=*D*_{T}:
2.4Equation (2.4) clearly demonstrates that the MTBL will be nested in the viscous sublayer for *Sc*≫1 because *δ*_{m}/*δ*_{v}∝*Sc*^{−1/3}.

Equations (2.1) and (2.3) will be used for the asymptotic solutions and are formally only valid very close to the conduit wall. To compare the predictions of the asymptotic solutions to the solutions to the full behaviour of the system, a model which accurately describes *u* and *D*_{T} throughout the entire conduit is required. Over the years, a multitude of models have been developed (Reynolds 1975; Weigand 2004); here a modified Van Driest mixing-length model has been selected. This turbulence model accurately reproduces the flow and turbulent diffusion in the inner layer, including the cubic dependence of *D*_{T} on *y* very close to the wall. For more details, see Hanna *et al.* (1981) and Sookhak Lari *et al.* (2010, 2011).

## 3. Derivation of far-field equations

Consider the transport of a high *Sc* solute through a conduit at high *Re* which exchanges mass with the conduit walls. For fully developed flow through a pipe with radius *R*, the governing equation is the axisymmetric Reynolds-averaged, steady-state mass transport equation (Bird *et al.* 2002)
3.1where *x* and *r* are the streamwise and radial directions, and *C*(*x*,*r*) is the (Reynolds-averaged) mass concentration. Streamwise diffusion has been neglected, which is permitted for a flow for which the Péclet number *Pe*=*ReSc* is very large. As stated before, we consider a hydrodynamically fully developed flow; *u* and *D*_{T} are therefore functions of *r* only. The axisymmetric coordinate system is used for convenience of presentation; the approach is equally valid for non-circular cross sections as long as the viscous wall region (∝*δ*_{v}) is much thinner than the local surface curvature.

Equation (3.1) is supplemented by a wall BC of the form
3.2where *C*_{w}=*C*(*x*,*R*), ∂*C*/∂*r*|_{w}=∂*C*/∂*r*(*x*,*R*) and *G*(*C*_{w}) is a generic function which depends on the wall concentration. Most BCs can be captured by equation (3.2), including the standard (linear) Dirichlet, Neumann and Robin BCs, but also higher order reactions and other BCs for which the wall-flux depends nonlinearly on the wall concentration. We will show that for BCs which satisfy (3.2), the Sherwood number *Sh*, which is the dimensionless mass flux, is universal and consistent with classical correlations of *Sh*(*Sc*,*Re*). The other two BCs are symmetry in the centre, ∂*C*/∂*r*|_{r=0}=0, and a constant concentration at the entrance: *C*(0,*r*)=*C*_{0}.

For the problem under consideration, the concentration is expected to be uniform, except within the MTBL. A dimensional analysis of the advective and diffusive terms of equation (3.1) in the MTBL results in
3.3where is a typical concentration and is the typical length scale for the streamwise variations. The typical velocity in the MTBL was estimated by evaluating equation (2.1) at *y*=*δ*_{m}.

The central premise of the approximation is that streamwise variations occur on much longer length scales than changes in the wall-normal direction, i.e. that (Notter & Sleicher 1972; Garcia-Ybarra & Pinelli 2006; Sookhak Lari *et al.* 2010). It follows from the estimates above that advection will be negligible relative to diffusion if the ratio satisfies
3.4where *Re*_{τ}=*u*_{τ}*R*/*ν* is the shear Reynolds number. The validity of this assumption will be established at the end of this section.

When (3.4) holds, equation (3.1) will no longer depend on *x*. Indeed, by neglecting the advective term, changing coordinates to mass transfer wall units *η*=(*R*−*r*)/*δ*_{m}, assuming that *δ*_{m}≪*R*, and using the cubic for *D*_{T} given in (2.3), the following linear partial differential equation is obtained
3.5where *c*=*C*/*C*_{0}. The equation above is equivalent to Garcia-Ybarra & Pinelli (2006, eqn 30) and Sookhak Lari *et al.* (2010, eqn 29). One BC is provided by equation (3.2) which in dimensionless form is given by
3.6where *g*=*δ*_{m}*G*/*C*_{0}. The second BC used is , which states that *c* tends to the bulk concentration *c*_{b} far away from the MTBL. Here, the dimensionless streamwise coordinate *ξ* is defined as . With these BCs, equation (3.5) admits the following closed-form solution (Garcia-Ybarra & Pinelli 2006; Sookhak Lari *et al.* 2010)
3.7where the wall concentration *c*_{w}(*ξ*) and *F*(*η*) are defined as
3.8and
3.9The function *F* increases monotonically from *F*(0)=−1 to , and .

To complete the approximation, an equation is required which governs the behaviour of *c*_{b}. Such an equation can be obtained by averaging equation (3.1) over the cross section
3.10where is the average streamwise mass flux. Because *c* is assumed constant throughout the cross section except in the MTBL, we can approximate 〈*uc*〉≈*Uc*_{b}, where *U* is the average velocity. This results in
3.11where *r*_{h}=*R*/2 is the hydraulic radius and *J*_{w} is the wall mass flux per unit area
3.12Substituting equation (3.12) into (3.11) results in
3.13and
3.14Here, we have used the standard pipe Reynolds number definition *Re*=2*UR*/*ν*. The typical length-scale can be interpreted as the distance for which the solute mass in the entire cross section is depleted by the flux through the wall. Equation (3.14) can be used to check the validity condition (3.4), which now takes the form of *b*^{2/3}*ReSc*^{2/3}≫1. Even for very low values of *Re*=2000 and *Sc*=100, *b*^{2/3}*ReSc*^{2/3}≈466 and therefore, condition (3.4) is satisfied for high *Sc* compounds in turbulent flows.

Equations (3.7), (3.9) and (3.13) comprise a set of coupled equations which can be used to construct asymptotic far-field solutions for the concentration profile and the mass transfer at the wall. Because the derivation does not rely on linear techniques such as separating variables, this includes nonlinear BCs. Examples will be discussed §5.

## 4. A universal Sherwood number equation

A universal expression for the Sherwood number *Sh* can be derived for BCs satisfying equation (3.6). The Sherwood number is defined by *Sh*=2*k*_{f}*R*/*D*, where *k*_{f} [*LT*^{−1}] is the mass transfer coefficient (Bird *et al.* 2002)
4.1Substituting equation (3.9) into equation (4.1), and using that *J*_{w}=*g*(*c*_{w})*DC*_{0}/*δ*_{m} results in
4.2i.e. *k*_{f} is *independent* of the type of BC. Using equation (4.2), the *Sh* equation for high *Sc* compounds is
4.3The universality is a direct consequence of the linear dependence of the wall concentration gradient and the concentration difference between wall and bulk, as is evident from (3.9). The underlying reason for this is the invariance of (3.5) to scaling because of its linearity.

The universality of equation (4.3) is confirmed by recent work with Robin BCs (Sookhak Lari *et al.* 2010) and Dirichlet BCs (Garcia-Ybarra & Pinelli 2006). It also compares favourably with experimental data. The Fanning friction factor *f* is defined as (Bird *et al.* 2002, eqn 6.1-4a). Substitution into equation (4.3) results in
4.4which corresponds well to the established correlation (Bird *et al.* 2002, eqn 14.2-5) upon substituting *b*=9.5×10^{−4} and *Sc*_{T}=1. By applying the Blasius formula *f*=0.0791*Re*^{−0.25} (Bird *et al.* 2002, eqn 6.2-12), we obtain *Sh*=0.016*Re*^{0.88}*Sc*^{1/3} which is in good agreement with the *Sh* correlations presented in the electronic supplementary material.

The Stanton number *St* is defined as the ratio of mass transfer coefficient to average velocity:
4.5Note that *St* is closely related to by . The value of *St* will therefore immediately give an indication of the appropriateness of neglecting the streamwise advection.

## 5. Far-field solutions

### (a) Linear boundary conditions

We will now provide closed-form far-field solutions for Dirichlet, Neumann and Robin BCs, by considering a linear BC of the form
5.1Here, *α*,*β* and *γ* are constants. Using equations (3.6), (3.9) and (5.1), the wall concentration and gradient are given by
5.2Substituting (5.2) into (3.13) and solving for *c*_{b} results in
5.3

The specific solutions for Dirichlet (*α*=1, *β*=0), Neumann (*α*=0 and *β*=1) and Robin (*α*=−*σ*, *β*=1 and *γ*=0) BCs are presented in table 1 as solution AS-D, AS-N and AS-R, respectively. The Neumann solution required performing a Taylor series expansion around *ξ*=0 and taking the limit of . Solution AS-R is equivalent to the solution derived in Sookhak Lari *et al.* (2010). AS-N is documented to a large extent in Bird *et al.* (2002, pp. 411–414), but that solution still contains an integral which needs to be approximated numerically. Garcia-Ybarra & Pinelli (2006) derived *F*(*η*) and (4.5) for Dirichlet BCs, but did not solve for *c*_{b} (solution AS-D).

The far-field asymptotic solution will be compared with solutions of equation (3.1) obtained with the method of separating variables. To allow comparison to the full solution, the method retains the cylindrical coordinate system and uses the realistic velocity and diffusivity profiles provided by the modified Van Driest Mixing length model. An expansion of the form
5.4is used, where *h*(*η*) is a function which maps (5.1) onto a homogeneous BC. The eigenvalues *k*_{n} and associated eigenfunctions *Y* _{n}(*η*) can be found by solving a Sturm–Liouville problem. The ordinary differential equation (ODE) governing *X*_{n} can be solved analytically and takes the form of a damped exponential if *k*_{n}>0 and a linear function if *k*_{n}=0. The eigenvalues and eigenfunctions do not have closed-form solutions and are determined numerically using a shooting method. Details of the method and implementation are discussed in the electronic supplementary material.

The solutions for Dirichlet, Neumann and Robin BCs using the method of separating variables are presented in table 1 and are denoted by SV-D, SV-N and SV-R, respectively. The parameter *Re*_{m}=*R*/*δ*_{m} represents the distance to the centre of the pipe in mass-transfer units.^{1} It is clear from table 1 that the structure of the two solution methods is very similar.

All results presented in this paper are for *Re*_{τ}=2000, *Sc*=1000 and *Sc*_{T}=1 unless stated otherwise. The prediction for AS-D is , which compares well with SV-D, which predicts *k*_{1}=0.8510. For AS-R with *σ*=2, , against SV-R which predicts *k*_{1}=0.5970. The small difference can be traced back to differences in the eddy diffusivity profile. The wall damping employed in the modified Van Driest mixing length model is purely empirical, and only satisfies cubic behaviour very close to the wall. Between 1<*y*^{+}<5, *D*_{T} is up to 30 per cent higher than a pure cubic. This will make the MTBL a bit thinner, which in return results in a slightly higher decay coefficient. This small difference could have been avoided altogether by adopting a modified value for *b*/*Sc*_{T} using the calculation method outlined in appendix A, but presenting the slight differences was deemed more instructive.

It is clear from (3.7) that (*c*−*c*_{w})/(*c*_{b}−*c*_{w})=1+*F*(*η*) regardless of the BC. For the method of separating variables (*c*−*c*_{w})/(*c*_{b}−*c*_{w}) in the far-field can be obtained by setting . The equations are presented in the third column of table 1 and are plotted in figure 2*a*. For SV-D and SV-R, (*c*−*c*_{w})/(*c*_{b}−*c*_{w}) is closely related to the first eigenfunction and the correspondence to the far-field asymptotic solution is excellent. The small differences originate again from the differences in turbulence model, and could be avoided using a modified value for *b*/*Sc*_{T} (appendix A). For SV-N, (*c*−*c*_{w})/(*c*_{b}−*c*_{w}) shows deviations very close to the wall which are caused by a truncation of the infinite sum. The convergence to the asymptotic solution can be seen to be quite slow as is evident from the profile for *N*=25 modes (plus symbols) and *N*=100 modes (squares).

### (b) Nonlinear boundary conditions

Approximate analytical solutions can be constructed for nonlinear BCs. As a demonstration, a second-order wall reaction
5.5is studied, where *a* is a dimensionless reaction coefficient. Substituting (5.5) into (3.9) gives
5.6Substitution of (5.6) and (5.5) into (3.13) results in
5.7The change of variables simplifies equation (5.7) to
5.8

which has a solution
5.9where , and *W*_{0} is the Lambert *W* function (Corless *et al.* 1996).

The far-field asymptotic solution in equation (5.9) was compared with the solution of a finite-volume approximation of the full three-dimensional axisymmetric partial differential equation (3.1) and mass transfer BC (5.5). As for the method of separating variables, the modified Van Driest mixing-length model was used to determine velocity and turbulent diffusivity profiles. The advective term was discretized using a first-order upwind scheme, which allows for explicit marching in the *x*-direction. The *r*-direction is discretized using second-order central scheme which is solved using direct matrix inversion. The nonlinearity of the BC (5.5) is incorporated using a simple iterative method.

The problem was solved for *a*=10^{−1}, 10^{0} and 10^{1}. Grid convergence was observed at *N*_{x}=1200 and *N*_{r}=600, although logarithmic spacing was required because of strong variations very close to the wall and near the entrance. The cell-sizes vary up to eight orders of magnitude. A conservative method such as the finite volume method is crucial for such extreme stretching (Mathias & van Reeuwijk 2009). The Grid Convergence Index (GCI; Roache 1994) for these simulations is GCI<1.2% in the far-field based on . Note that in the calculation of the GCI, we assumed that the method is entirely first order; the reported value for the GCI is therefore a conservative estimate.

The asymptotic solution for *c*_{b} (5.9) and the three-dimensional axisymmetric simulations are in excellent agreement (figure 3*a*), even though the asymptotic solution does not take into account entrance effects. In §6, we will show that this is the case because the entrance length is so small that it does not influence *c*_{b}.

## 6. Near-field solutions

In the near-field, streamwise advection will not be negligible. In dimensionless variables, (3.1) is given by
6.1Here, (2.1), (2.3) were used for the velocity and eddy diffusivity profiles, respectively. The small parameter *ϵ* is given by
6.2Hence, a suitable near-field coordinate is *ξ*_{e}=*ξ*/*ϵ*, for which (6.1), (3.13) become
6.3and
6.4The equations above are parameter-free which suggests universal behaviour, although it should be noted that the BCs may still introduce a parameter dependence. As *ϵ*≪1, (6.4) immediately results in *c*_{b}(*ξ*_{e})≈1. The entrance length-scale is given by
6.5Note that depends on flow properties only. This can be understood by realizing that is related to the time *T*_{D} it takes for the mass in the MTBL to deplete: . During this time, the boundary layer section moves at a typical velocity *δ*_{m}/*δ*_{v}*u*_{τ}, which is the velocity at the edge of the MTBL. The entrance length-scale can therefore be estimated by
6.6Note that the validity condition (3.4) for the approximation can be expressed using and as
6.7

In order to obtain closed form solutions for (6.3), the Von Karman–Pohlhausen integral method (Lighthill 1950; Spalding 1954; Schlichting & Gersten 2000) will be used. This method is not exact as it involves substituting the assumed concentration profile *F*(*η*). However, (6.3) does not admit self-similar solutions, because (i) the BC (3.6) is nonlinear and (ii) the total diffusion (1+*η*^{3}) does not allow power-law behaviour for the boundary layer thickness.

By introducing the concentration deficit
6.8and integrating from *η*=0 to *λδ*(*ξ*), (6.3) becomes
6.9here, *δ*(*ξ*) is the typical boundary layer thickness and *λ*>1 is a coefficient.

The concentration is assumed to be of the form *c*(*ξ*_{e},*ζ*)=*c*_{b}+(*c*_{b}−*c*_{w})*F*(*ζ*), where *ζ*=*η*/*δ*. In terms of the concentration deficit, the assumed profile is
6.10

By substituting (6.10) into (6.9) and changing variables to *ζ*, we obtain
6.11where
6.12Combining (3.6) and (6.10) results in
6.13and substituting the expression above into (6.11) finally yields
6.14

Equation (6.14) is a first-order nonlinear ODE. For the general BC (3.6), *g*′=*g*′(*c*_{w}) which introduces another dependence on *δ*. It may therefore be impossible to derive closed-form solutions for complicated BCs. For the linear BC (5.1) however, *g*′=−*α*/*β* evaluates to a constant. In table 2, the solutions to (6.14) are presented for Dirichlet (*α*=1, *β*=0), Neumann (*α*=0, *β*=1) and Robin (*α*=−*σ*, *β*=1, *γ*=0) BCs. The Dirichlet and Neumann BCs have the classic *x*^{1/3} dependence (Linton & Sherwood 1950; Kestin & Persen 1962; Shaw *et al.* 1963; Berger & Hau 1977). The Robin BC is more complex because *g*′ is finite, but essentially behaves like a Neumann BC when *σ*≪1 and like a Dirichlet BC when *σ*≫1.

Although the *Sh* correlation for the far-field is the same for all BCs satisfying (3.6), the near-field correlation *Sh*_{x}(*ξ*_{e}) is different. Indeed, the mass transfer coefficient is given by *k*_{f}=*DF*′(0)/(*δ*_{m}*δ*) and therefore
6.15

The coefficient *λ* controls, via *B*, the growth rate of the boundary layer. As the integral (6.12) is divergent for , *λ* has to be tuned to results from three-dimensional axisymmetric simulations of the full problem. The resolution of the simulations is *N*_{x}=1200 and *N*_{r}=600 and the GCI<1% based on . Good agreement is found for *λ*=1.5 and therefore *B*=0.4. The near-field and far-field solutions can be combined using . The correspondence of the equation for *Sh*_{x} with three-dimensional axisymmetric solutions for Dirichlet, Neumann and Robin BCs is good (figure 2*b*). Note that *Sh*_{x} for the Dirichlet and Neumann BCs are scaled by a factor 100 and 10, respectively. Also shown are the solutions obtained using the method of separating variables using the first 100 modes. There is good agreement with both the asymptotic solutions and the simulations, although more modes are required to describe the behaviour for *ξ*_{e}<10^{−2}.

Using the same substitutions as for the *Sh* correlation, the *Sh*_{x} correlation for Dirichlet BCs is given by
6.16This correlation is in good agreement with electrochemical mass transfer experiments which report *Sh*_{x}=0.184*Re*^{0.58}*Sc*^{1/3}(*x*/2*R*)^{−1/3} (Shaw *et al.* 1963; Berger & Hau 1977).^{2} For the asymptotic solution for Neumann BCs, the correlation is identical to (6.16) but the prefactor is 0.22.

For nonlinear BCs, *g*′ depends on *c*_{w}. For the second-order wall reaction (5.5), substitution into (6.13) results in a quadratic in Δ_{w}, of which the physically relevant root is given by
6.17Therefore, *g*′ is given by
6.18Substitution of the equation above into (6.14) and solving the ODE results in
6.19The good agreement of *Sh*_{x}/*Sh* based on the equation above with the numerical simulation is shown in figure 3*b* for three different values of *a*. Note that *B* was kept the same value as for the linear BCs.

## 7. Concluding remarks

In this paper, asymptotic solutions for high *Sc* scalars for both linear and nonlinear BCs were developed, which provide simple closed-form solutions and predict accurately the concentration profile and mass transfer in the near- and the far-field. It was shown that in the far-field, the mass transfer coefficient *k*_{f} and associated dimensionless Sherwood number *Sh* is *independent* of the specific wall BC. This result is valid for all BCs satisfying (3.6), which include virtually all known BCs including those of Dirichlet, Neumann and Robin type. This is a remarkable result which emphasizes that for high *Sc* turbulent flows, the mass transfer coefficient *k*_{f} (4.2) depends only on the molecular diffusivity *D* and the thickness of the MTBL *δ*_{m}.

The method presented in this paper assumed that the eddy diffusivity *D*_{T} is a cubic in the distance to the wall. As discussed in §2, the precise near-wall behaviour of *D*_{T} is subject to some uncertainty. However, the universality of the mass transfer coefficient is independent of the exact profile for *D*_{T}. Indeed, making use of a different profile for *D*_{T}, for example, the fourth-order polynomial suggested in Garcia-Ybarra (2009), would have led to the same conclusion.

The solutions reported in this paper can be used even when the actual *D*_{T} profile is not cubic by calculating the effective value for *b*/*Sc*_{T} using the procedure described in the appendix. This procedure matches *δ*_{m} from the actual *D*_{T} profile to that of the assumed cubic, thereby ensuring that the integral parameters are predicted accurately. Note that this would make the parameter *b*/*Sc*_{T} dependent on *Re* and *Sc*. As discussed in §2, the *D*_{T} profile may be different for specific BCs. In that case, *b*/*Sc*_{T} would have a different *Sc* and *Re* dependence for each BC and therefore the *Sh* correlation (4.3) would not be universal. It is therefore desirable that further research focuses on the near-wall profile of *D*_{T} at high *Sc* turbulent mass transfer, in particular the influence of different BCs.

## Appendix A. Inferring the value of *b*/*Sc*_{T} from turbulence models

The coefficient *b* has a strictly defined physical meaning stated in (2.2). However, it can also be treated as a free parameter representing a measure for the ‘conductivity’ of the MTBL, in which case the value can be determined from the turbulence model employed. We start from (3.1), neglect horizontal advection and change to plus-units using the change of variables *r*=*δ*_{v}(*Re*_{τ}−*y*^{+}). The result is
A1

Integrating twice and using (3.2) results in
A2and the mass transfer coefficient *k*_{f} is therefore given by
A3Many theoretical studies on mass transfer approximate the integral on the right-hand side to develop *Sh* correlations (Kader & Yaglom 1972; Aravinth 2000).

By equating the equation above to (4.2) and using (2.4), we obtain
A4which allows for the calculation of an effective conductivity parameter *b*/*Sc*_{T} for each given turbulence model.

## Footnotes

- Received September 6, 2011.
- Accepted January 23, 2012.

- This journal is © 2012 The Royal Society