We provide and analyse a model for the growth of bacterial biofilms based on the concept of an extracellular polymeric substance as a polymer solution, whose viscoelastic rheology is described by the classical Flory–Huggins theory. We show that one-dimensional solutions exist, which take the form at large times of travelling waves, and we characterize their form and speed in terms of the describing parameters of the problem. Numerical solutions of the time-dependent problem converge to the travelling wave solutions.
Micro-organisms, in their various forms, constitute most of the living matter on the Earth, and they are found virtually everywhere: in the atmosphere, oceans, soils, rocks, ice. Popular culture is enthused by the concept that they may reside in deep Antarctic lakes, or on the frozen wastes of Mars.
They are also of extreme practical importance, both in the environment and in industrial processes. In soil, bacteria play a primary role in the nitrogen cycle, and their metabolic processes serve to clean groundwater of contaminants, both natural and man-made (Chapelle 2001). In wastewater and sewage treatment, bacterial action causes the breakdown of complex organic contaminants (Henze et al. 2002). However, bacterial growth is also a cause of infections in the medical and pharmaceutical industries, as well being a primary cause of disease in the human body—tooth decay, cystic fibrosis and urinary infections being some relevant examples (Costerton et al. 1999; Bryers 2008). It is therefore of importance to characterize the way in which bacteria grow and affect the environment.
Bacteria often grow on surfaces, forming a layer which is known as a biofilm. Depending on the nutrient conditions, a biofilm can be of typical thickness of order 100 μm. In the usual situation the bacteria are supplied with nutrient by a fluid flow across the free surface of the biofilm, whether this be a slow groundwater flow or a fast flow in an industrial pipe. The bacterial film grows as it is supplied with nutrient from the fluid flow. Nutrient uptake causes a nutrient concentration gradient within the film. As the biofilm becomes thicker, growth is localized near the biofilm surface. Bacteria at the base of the film reliant on the same nutrient become progressively starved and may die and lose their adhesion to the wall; in this case portions of the film may ‘slough off’, thus providing a sink term for the biofilm growth. Erosion at the biofilm surface owing to flow can also provide a sink term. In other situations, the growth of the biofilm may restrict flow, causing ‘pore clogging’ in the small pore space of a soil and ‘biofouling’ in industrial systems; such growth may be self-limiting owing to the consequent restriction in nutrient supply.
Because of their widespread prevalence and importance, it is desirable to be able to characterize the growth and structure of biofilms in terms of the natural parameters of the system, and this is the intention of the present paper. Existing models of biofilm growth based on physical mechanisms for spreading tend to be computationally very intensive. To the extent that rational mathematical simplifications can be made, such models will be more easily applied to address issues of direct practical importance.
(a) Microbial metabolism
Like all living things, microbes survive by generating energy from nutrients through a variety of metabolic reactions. This process involves a network of redox reactions, and involves the overall exchange of electrons between two distinct chemical fuels, which are consumed in the reactions; the metabolic process is in this case called respiration. While there may be a number of such fuels, there is a hierarchy in their use. Dissolved oxygen is commonly the terminal electron acceptor (as the externally sourced oxidant is typically referred to), while an organic carbon compound is the electron donor. When these preferred substrates are absent or depleted, other compounds can be used instead. When the same organic compound is used as both donor and acceptor, the metabolic process in this case is called fermentation. Many bacteria are able to use several reaction pathways independently, giving them a degree of flexibility to differing conditions. This capability is very species dependent, and competition ensures that the species best adapted to local conditions dominate in natural environments.
Microbial growth depends heavily on energy metabolism but also requires the uptake of other substrates needed to generate new biomass. Growth rate is generally limited by the supply of one or more substrates, but saturates to a maximum growth rate in conditions of ample supply. The dependence of bacterial growth rate is commonly taken, by analogy with simple enzyme kinetic uptake rates, as proportional to c/(K+c), where c is the relevant nutrient concentration and K is a constant. Such kinetics are called Monod kinetics (Monod 1949); when two nutrients control growth, as in respiration, it is usual to take the growth rate as proportional to the product of two Monod factors (Bader 1978). The rate of nutrient uptake is commonly modelled as being proportional to growth rate, with the constant of proportionality (1/Y) expressed in terms of a yield Y for that substrate. A variety of enhancements to this simple model have also been proposed to account for maintenance (i.e. non-growth) nutrient consumption, inactivation of cells in adverse conditions and other observed effects (Beeftink et al. 1990; Wanner et al. 2006).
(b) Biofilm constitution and rheology
Biofilms typically grow by the initial attachment of bacteria to a substratum, such as a vessel wall. The adhered bacteria then proliferate, and, as they do so, they produce a polysaccharide-rich matrix, within which the bacteria continue to grow. Thus biofilm is characterized by a matrix of extracellular polymeric substance (EPS), whose nature depends on the micro-organisms as well as on the physical and chemical environment (Sutherland 2001a). Biofilm composition is highly variable, but the cellular content of mature biofilm is typically only 2–5% (Sutherland 2001b). The extracellular component is a viscous gel-like slime formed by perhaps 1–5% EPS (representing 50–90% of the total organic matter in the biofilm) in water. The EPS is normally dominated by long-chain polysaccharides (40–95%) secreted by the cells and can also have significant protein, lipid and nucleic acid content. A standard biofilm model organism is Pseudomonas aeruginosa, whose primary EPS constituent is alginate polysaccharides with molecular weight 1–2 MDa, so that individual molecules are as long as 5 μm (Flemming & Wingender 2001). These polysaccharides absorb and retain a large volume of water and are the main source of a biofilm’s physical properties: the rheology can be characterized as that of a viscoelastic fluid, possibly with a yield stress, and with a relaxation time of the order of a few minutes (Klapper et al. 2002; Stoodley et al. 2002; Shaw et al. 2004). On the time scales associated with biofilm growth (days), it is reasonable to treat the biofilm as a viscous material.
In our study of biofilms, we will consider the EPS to have the properties of a polymer solution, in which long-chain polymeric molecules sit within a fluid solvent. The structure of our model is, therefore, that of a two-phase material. A simple description of this material uses Flory–Huggins theory (Flory 1953; Tanaka & Fillmore 1979; Graessley 2004) to describe the free energy associated with the presence of a volume fraction ϕ of long-chain polymer molecules in solution based on a lattice model. This free energy of mixing has entropic and enthalpic contributions. The entropic contribution favours imbibing more solvent at all concentrations, giving a swelling (or dilution) tendency. The enthalpic contribution results from the net effect of monomer–monomer, solvent–solvent and monomer–solvent interactions, and can give either a swelling or a contraction tendency. Biofilms do not swell indefinitely—they have a characteristic slimy but cohesive texture—so we are interested in the ‘poor solvent’ regime for which the dominance of like-interactions over unlike-interactions provides a free energy minimum in the range ϕ>0 corresponding to a swelling equilibrium. Alternative mechanisms such as transient chain entanglements or effective cross-linking between polyanionic chains via cations in solution may provide plausible explanations for an equilibrium composition (Körstgens et al. 2001), but current experimental evidence is inconclusive as to the true basis so we initially adopt our model for simplicity. The Flory–Huggins free energy per unit volume is manifested in the dynamic equations as an osmotic pressure, which is a measure of the additional pressure that must be applied to equilibrate the polymer solution with pure solvent across a semipermeable membrane. In the Flory–Huggins theory, the osmotic pressure of a polymer solution is given by 1.1where EL>0 is the lattice energy density and χ is the Flory interaction parameter. Owing to the high polymerization of many bacterial polysaccharides, the ratio n of polymer and solvent molar volumes is numerically large. The poor solvent regime describes the situation in which there exists a swelling equilibrium Ψ=0 at some ϕ>0, rather than a tendency for the polymer simply to dissolve. In the long-chain limit this corresponds to values of χ>1/2. In the poor solvent regime the solution can be subject to a phase separation instability at low polymer volume fractions such that Ψ′<0, in which case an additional free energy term in |∇ϕ| due to compositional inhomogeneity is significant and regularizes the instability (Cahn & Hilliard 1958). In this paper, we follow Cogan & Keener (2004) in considering pure growth without external stresses (due to flow interaction, for example), which results in EPS volume fractions well above the phase separation regime. As a result, we do not address phase separation issues here and we neglect the inhomogeneous contribution, as do Cogan & Keener.
(c) Previous work
There have been numerous efforts to model the growth of biofilms. Mostly these start from the concept of a nutrient which diffuses into the growing biofilm, being absorbed as it does so. Early models were one-dimensional (e.g. Rittmann & McCarty 1980; Wanner & Gujer 1986), but the application of confocal microscopy techniques and the consequent imaging of biofilm growth led to the revelation that biofilm growth often occurs through the heterogeneous growth of three-dimensional, tower- and mushroom-like structures. Thus, more recently, three-dimensional models have been put forward. Eberl et al. (2001) develop a continuum model based on diffusion and uptake of nutrient, with a biomass spreading term described by an ad hoc nonlinear diffusion coefficient. Dockery & Klapper (2001) and Klapper (2004) provide a similarly motivated model, in which the biofilm is modelled as a fluid medium which is porous in the sense that its motion is governed by Darcy’s Law, although a separate phase is not identified.
A more sophisticated type of model is that of Cogan & Keener (2004, 2005), who specifically consider the biofilm as a biological gel consisting of EPS and water, and whose conceptual model forms the basis for our study.
They provide a two-fluid continuum model based on polymer solution theory, in which they make an implicit assumption of an appropriate asymptotic limit which provides a dominant balance between viscous stress and the osmotic pressure term. They find a one-dimensional steady state with constant interface growth velocity and examine the instability of this solution to spatial perturbations. Our aims are essentially the same as those of Cogan & Keener, insofar as we seek to provide a description first of a uniformly growing biofilm, and subsequently study its stability to spatially heterogeneous modes. Our approach deviates markedly from theirs in that we use an explicit non-dimensionalization, choosing ‘natural’ scales where possible, to identify a different dominant balance, as we elaborate below.
Klapper & Dockery (2006) use a polymer solution model to investigate the effect of phase separation on steady states of biofilm spatial configuration in the plane of the solid surface. Since no parameters are estimated, the model serves as a theoretical analysis rather than an assertion that this mechanism is at play in real biofilms.
More recently, Zhang et al. (2008a,b) proposed a mixture model in which a single momentum equation describes the motion of the fluid mixture, and the relative velocity of the polymer and solvent components is given by an osmotic pressure term which includes a contribution due to solution inhomogeneity. The use of a mixture model is motivated by difficulties with boundary conditions in two-fluid models. They provide a framework which allows for a range of constitutive laws for the deviatoric stresses on the mixture components, but ignore the effect of these stresses on the relative velocity between components. They investigate one-dimensional solutions analytically and numerically, and simulate the interaction with fluid flow in the bulk solution in two dimensions.
In parallel to this development of continuum models has been an increase in discrete biomass models described by cellular automata, and more recently an exuberance of individual-based models, which have provided a framework for simulating a vast array of scenarios including heterogeneity, multi-species consortia, detachment and disinfection (e.g. Kreft et al. 2001; Picioreanu et al. 2004; Xavier et al. 2004, 2005; Alpkvist et al. 2006; Chambless et al. 2006), and arguably these models represent the current mainstream for theoretical approaches to biofilm growth.
The structure of the present paper is the following. In §2, we present the polymer–solvent model, following Cogan & Keener (2004). It is then non-dimensionalized, and in §3 we reduce the model to a simpler form. It is at this point that we diverge dramatically from Cogan and Keener’s work. We then study the reduced model in one dimension, and show in §4 that at large times a travelling wavefront solution exists. A useful simplifying approximation is identified in a parameter regime corresponding to denser biofilms. In §5 we confirm numerically that the travelling wave solution does develop from solutions to the time-dependent one-dimensional model.
2. The polymer–solvent model
We consider a biofilm in 0<z<s, as shown in figure 1. The biofilm consists of a matrix of EPS surrounded by nutrient-rich water, behaving as a gel, and with an osmotic pressure described by the Flory–Huggins theory. We also take both the polymer and the interstitial water to be viscous, and they interact via an interfacial drag term. For simplicity, we suppose that the biomass growth is rate-limited by one nutrient, for example the dissolved oxygen in water acting as the electron acceptor. Biomass normally constitutes a small volume fraction; in any case, we lump the biomass in with the EPS in the model and in so doing ignore the distinction between cellular growth and EPS production.
Let ϕ be the volume fraction of polymer, and 1−ϕ that of water. We then have two conservation of mass equations for polymer and solvent, which take the form 2.1where v and w denote the phase-averaged velocities of the polymer and solvent; g is a growth term whose form is discussed below. The concentration c is that of the assumed rate-limiting nutrient, and satisfies a conservation law of the form 2.2where D is a diffusion coefficient and r represents nutrient uptake by biomass.
The growth and uptake terms are related by a yield coefficient, and we assume Monod kinetics. A maintenance term is often included to account for nutrient uptake which does not lead to growth (e.g. Picioreanu et al. 1998), and this term can be significant in conditions of low nutrient concentration. For simplicity and comparison with the models of Cogan & Keener (2004) and Zhang et al. (2008a) we neglect this term in the current work.1 Specifically, we take 2.3where the constants G, R and K take typical values which are given in table 1. We take the momentum equations for slow two-phase flow in the form 2.4
In these equations, p is the fluid pressure, Ψ is the osmotic pressure given by equation (1.1) μ is the long-time viscosity of the polymer matrix, μw is the viscosity of water and the term in f is an interfacial drag term due to microscale viscous resistance. The interfacial drag refers to the (equal and opposite) forces exerted by each phase on the other owing to shear traction at the interface. In the normal way of things, the solvent viscous term is negligible, in which case equation (2.4)2 is just Darcy’s Law, and f is given by 2.5where k is the permeability. This relation simply defines what is conventionally meant by permeability. Values of the various parameters are given in table 1. In particular, we will assume that f is constant.
The equations (2.4) are the same as those given by Cogan & Keener (2004). At this stage, Cogan & Keener (2004) argue that, since ϕ is small, the interfacial drag term can be ignored. This leads them to an approximate model which is quite different from the one we derive here, where we find that in fact the interfacial drag term is dominant.
We seek scales for the eight variables ϕ, v, w, t, x, p, Ψ and c, and these are chosen by suitable balances in the equations. Specifically, we balance all the terms in the polymer mass conservation equation, scale w with v, diffusion with uptake of nutrient, interfacial drag with pressure gradient and osmotic pressure gradient. A comment on the osmotic pressure scale is appropriate. We can write equation (1.1) in the form 2.6Apart from the fact that the derivation of this expression assumes numerically small ϕ, it is in any case generally found that ϕ≪1 in most biofilms (Körstgens et al. 2001; Sutherland 2001b). Taking, in addition, the long-chain limit implies a solution of sparse, long molecules. We follow Cogan & Keener (2004) in adopting this limit and note that it assumes . Therefore, we may expect that, as an approximation, we may take 2.7Further, if Ψ′>0 for all ϕ (a good solvent), there is a tendency for the gel to swell indefinitely in contact with pure solvent, whereas if χ>1/2, then a stable gel fraction can exist where Ψ=0. This suggests that we balance 2.8expecting that this will be appropriate if the poorly constrained Flory parameter χ is close to (and greater than) one-half. Adopting these balances, we choose the scales 2.9
The length scale is based on a balance of the uptake and diffusion terms in the nutrient equation, and represents the depth of the active layer of biofilm. The time and velocity scales are chosen to balance the growth rate. The pressure scale reflects a balance between the interfacial drag and pressure terms in equation (2.4)2. The scale ϕ0 is chosen by balancing the osmotic and interfacial drag terms in equation (2.4)1. Based on the values in table 1, we can derive the values of the scales in table 2.
The parameters are defined by 2.11and typical values based on the values in tables 1 and 2 are given in table 3. α is the ratio of growth to uptake time scales. β gives the relative importance of polymer viscous stress to osmotic stress, and γ gives the relative importance of water viscous stress and pressure. With the exception of κ and, presumably, λ, the parameters are all small, and a reduced model is derived in the following section. From the equilibrium volume fraction ϕeq derived from equation (2.7), we can write λ=ϕeq/6ε. Typical EPS volume fractions of 1–5% then suggest values for λ in the range 0.5–2.3. The assumption that β is small is in contradiction to the assumption of Cogan and Keener, whose model is based on the implicit assumption that β is large (and Ψ and p are rescaled with β; Cogan & Keener did not explicitly non-dimensionalize their model, however, assuming essentially that the interfacial drag term would be small for sufficiently small ϕ0). From the definition of β, which can be written as β=μRϕ0/fDc0, we see that in fact β is small when ϕ0 is small. Biofilm viscosity is highly variable, and there is corresponding variability in the value of β. Nevertheless, our estimate of β≈1.1×10−4 using the parameter values of table 1 suggests that polymer viscosity would have to be dramatically increased before this term could become significant.
Boundary conditions for the model are discussed in §3.
3. Reduction of the model
The parameters α, β, γ and ε are all small. If we put them to zero, then the model is reduced to 3.1The approximations that β=0 and γ=0 are, at least potentially, singular limits, and they raise the possibility that boundary layers may occur in the model, depending on the nature of the boundary conditions which are imposed. This is an issue whose detailed consideration is postponed to a future three-dimensional study. In the present one-dimensional analysis, it is redundant. The equation for p uncouples, and we are left with three equations for ϕ, c and w in the form 3.2We note that for values ϕ<4λ equation (3.2)1 shows negative diffusion and is ill-posed as a result. This corresponds to the phase separation regime which we do not address here. The instability is regularized by inclusion of an inhomogeneous free energy term which provides a fourth-order spatial derivative of ϕ in this equation (Klapper & Dockery 2006).
Suitable boundary conditions for these equations are as follows. We suppose the biofilm is attached to the wall at z=0, where there is no normal flow or flux of nutrient. From these we find 3.3At the free surface, we have a kinematic condition that the interface moves with the polymer, a prescribed pressure and nutrient concentration,2 and we suppose that the polymer fraction in contact with free solvent is at the stable swelling equilibrium Ψ=0. These lead to the conditions 3.4The issue of the prescription of ϕ at the interface is based on the idea that the interface between biofilm and pure liquid is at local thermodynamic equilibrium, so that the free energy is continuous. This implies either ϕ=6λ or ϕ=0, and we presume the latter is precluded because it is unstable if . The kinematic condition provides the extra condition to determine s, and we see that we have two conditions for ϕ and c at z=0 and z=s, as necessary for equations (3.2).
There is an issue in the simplification afforded by equation (3.1), and that lies in the determination of the solvent velocity. w is only described by the equation (3.1)2 for its divergence, and while this is sufficient to determine w (=0) in one dimension, it is insufficient in three. The question thus arises as to what extra condition has gone missing.
The answer to this appears to lie in one of the vagaries of two-phase flow models. If we return to the two momentum equations in (2.10), keep γ=0 (corresponding to Darcy’s Law) but retain the viscous term in β, the equations can be written equivalently as 3.5from which it still follows that Ψ≈−p, and thus 3.6but also, by taking the curl of equation (3.5)1, we have the exact result 3.7This provides the extra relation that we seek to determine w. Its derivation is a little like the use of the Prandtl–Batchelor theorem in showing that closed streamlines contain fluid with constant vorticity in high Reynolds number flow.
If we adopt equation (3.5), then we need the extra appropriate boundary conditions for a viscous flow; these are the no-slip condition at the wall, and a no-stress condition at the free interface, 3.8where n is the normal to z=s, and ti, i=1,2 are the two tangent vectors there. The issue of how to carry through the limit is not addressed here, since, for the present, the one-dimensional solutions which we give below do not involve the use of equation (3.7), which is then satisfied identically.
4. Travelling waves
In one dimension the reduced model can be written in the form 4.1where 4.2subject to 4.3
It follows that 4.4and thus that ϕ satisfies the nonlinear diffusion equation 4.5where we write 4.6(since ϕP′=Ψ′). We therefore require ϕ>4λ (no phase separation) for the current problem to be well posed.
It is natural to expect that the solution of the one-dimensional model will tend towards a steady state at large times, in which ϕ and c are functions of the travelling wave variable 4.7and the appropriate conditions in the far field are 4.8both U and must be determined as part of the solution. With primes denoting differentiation with respect to η, the equations become 4.9
A first integral of these equations implies that 4.10while the kinematic condition implies that 4.11where ϕ′0=ϕ′(0). Together with equation (4.10), this implies that 4.12It follows that equation (4.10) can be written in the form 4.13
Our numerical strategy for solving the equations (4.13) for ϕ and (4.9)2 for c subject to the boundary conditions, from equations (4.3) and (4.8), 4.14is then to use finite differences on a truncated domain η∈[0,H] solving the iterative scheme 4.15with H chosen sufficiently large such that c′ and ϕ′ are small at the truncated boundary. We use the value ϕ(n)(H) to approximate , and in lieu of equation (4.14)2 we apply an approximate boundary condition in terms of a ghost point 4.16
Figure 2 shows the numerically computed travelling wave profiles of ϕ and c for λ=1 and λ=2 and constant κ=0.35.
Since , we see that λ may be large for dense biofilms. It thus makes sense to seek an approximate solution when λ is large, so long as is still small. This is easy to do; by rescaling 4.17we find that Φ is almost constant, assuming 36λ3≫1, and, in terms of the unscaled variables, the solution is approximately given by 4.18
This approximation is shown in figure 2 as dashed lines, and is reasonably accurate even for only modestly large values of λ: the maximum relative error between the numerical solution for ϕ and the large λ approximation decreases from 2×10−3 at λ=1 to 3×10−5 at λ=2.
5. Time-dependent solution
For numerical solution of the time-dependent one-dimensional problem we transform to a fixed spatial domain ξ∈[0,1] by setting ξ=z/s(t), and we define σ=s2. The system becomes 5.1with boundary conditions ϕ=6λ, c=1 at ξ=1, and ϕξ=cξ=0 at ξ=0. Equation (5.1)2 results from the kinematic condition. A second-order finite difference discretization of the spatial variable yields a tridiagonal system of stiff nonlinear ODEs for equations (5.1)1 and (5.1)2 coupled to a quadrature for c from equation (5.1)3. We integrate the system with Matlab’s standard stiff solver ode15S, choosing suitable initial conditions and ϕ(0,ξ)=ϕ0(ξ).
Numerical solutions of the time-dependent problem converge to the travelling wave solution, as typified in figure 3. This convergence is robust to varying initial conditions provided ϕ(0,ξ)>4λ for all ξ, and transients due to initial conditions decay very rapidly. Figures 4 and 5 show the result of varying λ and κ, respectively.
The form of the travelling wave velocity (4.12) suggests the quantity |cz|z=1/ϕz=0 as a plausible estimator of the time-dependent free surface velocity st. This approximation, shown with crosses ‘×’ in figures 3c, 4a and 5a, is remarkably good in our simulations and only differs appreciably from the calculated velocity when λ is small (figure 4a). This can be explained by applying the scaling of equation (4.17) to the time-dependent problem, which reveals that the time derivative term in equation (4.1)1 is of the order of (36λ3)−1. In the 36λ3≫1 regime the biofilm, therefore, evolves quasi-statically, and a relation equivalent to the travelling wave velocity (4.12) still applies. Figure 4a shows that the approximation is still relatively good for values as low as λ=0.5.
The proposed model describes biofilm growing on an impermeable substratum based on the material properties of a polymer solution. The biofilm cellular volume fraction is assumed to be small. Estimates of typical scales lead to a different mathematical description from earlier models.
Solutions of the model converge to a travelling wave solution strongly dependent on the dimensionless parameter EPS content of the biofilm is predicted to increase away from the biofilm surface, though for realistic values of λ this increase is likely to be too small to be experimentally verifiable. For EPS-rich biofilms, a formal assumption of λ≫1 leads to a significant simplification of the model while still providing a reasonable approximation for values of λ only modestly greater than 1. This may have useful relevance to many real biofilms.
The travelling wave solution provides an expression 6.1for the growth velocity of the biofilm surface in dimensional terms, where the nutrient gradient is evaluated at the biofilm surface and the EPS volume fraction is taken at the biofilm base. This expression also provides a good approximation in the time-dependent case of relatively thin biofilms for which nutrient penetrates the full biofilm depth, and thus provides reasonable scope for comparison with experiment.
We acknowledge the support of the Mathematics Applications Consortium for Science and Industry (www.macsi.ul.ie) funded by the Science Foundation Ireland mathematics initiative grant 06/MI/005. This publication has emanated from research conducted with the financial support of Science Foundation Ireland under grant no. SFI/09/IN.1/I2645. H.W. acknowledges the support of the Life Sciences Interface DTC (www.lsi.ox.ac.uk) funded by EPSRC grant GR/R96149/01.
↵1 We note, however, that it may lead to phase separation at the base of thick biofilms and be a mechanism for sloughing.
↵2 A diffusive boundary layer assumption is more appropriate in many situations. For fast or well-mixed flow a diffusive boundary layer Robin condition effectively reduces to a Dirichlet boundary condition.
- Received June 22, 2010.
- Accepted November 2, 2010.
- This journal is © 2010 The Royal Society