## Abstract

A model system for theory and experiment which is relevant to foam fractionation consists of a column of foam moving through an inverted U-tube between two pools of surfactant solution. The foam drainage equation is used for a detailed theoretical analysis of this process. In a previous paper, we focused on the case where the lengths of the two legs are large. In this work, we examine the approach to the limiting case (i.e. the effects of finite leg lengths) and how it affects the performance of the fractionation column. We also briefly discuss some alternative set-ups that are of interest in industry and experiment, with numerical and analytical results to support them. Our analysis is shown to be generally applicable to a range of fractionation columns.

## 1. Introduction

The process of foam fractionation uses a flowing foam to separate out or concentrate surface-active components of a solution. In the most basic form of fractionation, the liquid solution is foamed through a vertical column (with the surface-active molecules preferentially attaching to the foam films), and the overflow at the top is collected (and collapsed to liquid). This approach has proved difficult to model and simulate numerically owing to the uncertainty in accurately describing the overflow, where drainage, bubble rupture and rheology may all contribute (although careful analysis may allow such complexity to be reduced to a more tractable system [1], which can then be applied to the modelling of fractionation). This difficulty has hindered development of full analytical treatments, despite the widespread use of foam fractionation in several industries such as chemical engineering and biotechnology [2].

In a previous paper [3], we presented an analytical model of a foam fractionation apparatus consisting of an inverted U-shaped tube that connects two reservoirs (shown schematically in figure 1). Gas is bubbled into the feed reservoir, causing foam to flow through the tube at some velocity *V* . This leads to the transfer of enriched liquid from the feed to the output reservoir with a liquid flux *J*(*V*) (volumetric flow rate per unit area). This arrangement enables us to use well-defined boundary conditions. The analysis of the foam and liquid flow in the U-tube is based on a modified form of the foam drainage equation [5]. The previous paper [3] also contained experimental results from U-tubes of the type shown in figure 1 which qualitatively and semi-quantitatively agreed with the numerical and analytical predictions. Foam drainage theory was previously used as the basis for analysis of fractionation by Neethling *et al.* [6], and others as cited in [3].

In this study, we examine in more detail the effects of finite system size, and specifically the effect on the performance or efficiency of the system.

## 2. Background of the inverted U-tube

### (a) Main analytical results

When the gas flow rate (i.e. the gas velocity divided by the cross-sectional area of the tube) is constant, there is a steady flow of both gas and liquid from left to right. The key question in this context is: *what is the liquid flow rate that is delivered?*

The main result of the previous work [3] provides the liquid flux *J* as a function of gas velocity *V* ,
2.1where *J*_{0} is the limiting value for the liquid flux *J*(*L*) as the length of both legs [3]. The parameter *c*_{1} and other relevant physical quantities are discussed in §2*b*. For this work, we look at set-ups where the leg lengths *L*_{l} (left, inflow) and *L*_{r} (right, outflow) are finite. It is shown in §3*a* that the left or inflow leg is key to the flow analysis. We thus write
2.2where as (leading to *J*_{0} as expected).

The *J*=*J*_{0} condition has been identified in various forms recently by other authors, as noted in reference [3]. There is also at least one much older statement to the same effect, by Desai & Kumar [7]. Lacking from all of these treatments, however, is an analysis of the effect of finite tube length, which is the primary objective here. As equation (2.1) indicates, there will generally be a compromise between low *J* (which results in a higher output concentration) and high *V* (resulting in a higher rate of delivery).

The design of fractionation columns may be considered to be one of *optimization*. The surface-active molecules associated with the foam surfaces are delivered at a rate proportional to the gas velocity *V*, whereas the bulk dilute liquid is delivered at a rate proportional to *J*. We wish to increase the delivery of surface-active molecules while keeping the bulk liquid delivery to a minimum.

The choice of leg length enters such considerations, through equation (2.2). However, practical considerations (the available space for the column, for example) constrain the values of *L* and *V* , and there may exist threshold gas velocities below which the U-tube will not operate as drainage and coarsening will destroy any foam before it overflows the column (or flows through the bend). A further factor that must be considered in real applications is the adsorption time, i.e. the time it takes for the surface-active chemicals to attach to the foam film surfaces [8,9].

In this paper, analytical results for the length dependence are presented, and an explicit expression for *ϵ* in equation (2.2) has been obtained. Our analytical results are derived from a steady-state version of the channel-dominated foam drainage equation [5]. In a previous paper [3], we made some useful generalizations to the equation, but we do not discuss them in detail here. Fits of numerical results are compared with this prediction, with good agreement observed.

### (b) Physical parameters of the theory

(Here, we set out the units to be used; the reader who is mainly interested in the analytic theory may wish to proceed to §2*c*.) The theory developed in [3] involved two physical constants, *c*_{1} and *c*_{2}. *c*_{1} has dimension of velocity and *c*_{2}/*c*_{1} dimension of length,
2.3and
2.4

These constants contain the following physical parameters: the density difference of gas and liquid *ρ*, surface tension *γ*, effective viscosity *η**≃150*η* (with liquid viscosity *η*), gravitational acceleration *g* and average bubble volume *V*_{b}, together with numerical constants related to foam structure and geometry: a constant related to the cross section of Plateau borders, , and the total length of Plateau borders per unit volume of foam .

*c*_{2}/*c*_{1} is (to within a constant of order unity) equal to the height of that section of wet foam that exists owing to capillarity when a foam is in contact with underlying liquid under gravity [10]. The parameter *c*_{1} is (again to within a constant of order unity) equal to an important constant in elementary drainage theory [10]. In steady uniform drainage, the flow velocity is proportional to liquid fraction with the constant of proportionality given by *c*_{1}.

In this study, we present results scaled by the constants *c*_{1} and *c*_{2}. The reader may therefore apply our results to a system of interest simply by calculating these constants and rescaling appropriately. To provide some feeling for the typical practical values, we may use the following: *η*=0.001 Pa s, *γ*=0.05 N m^{−1}, *ρ*=1000 kg m^{−3}, *V*_{b}≈4×10^{−9} m^{3}, giving *c*_{1}≈3.2×10^{−2} m s^{−1} and *c*_{2}≈4.7×10^{−5} m^{2} s^{−1}. This gives a value for *c*_{2}/*c*_{1}≈1.47×10^{−3} m.

Using these values, the simulations presented here are for fractionation columns with leg lengths between 10 cm and 1 m, with gas velocities between 1 and 20 mm s^{−1}.

### (c) Drainage theory and boundary conditions

As in previous work [3], our analysis is based on a modified version of the channel-dominated foam drainage equation [5]. This is a partial differential equation for the spatial and temporal variation of the local liquid fraction *ϕ*(*x*,*t*) [11]. For completeness, a brief background of the theory is reproduced in the following text.

The steady-state form of the equation (i.e. ∂*ϕ*/∂*t*=0) relates the liquid flux (volume flow rate per unit cross-sectional area) *J* to liquid fraction *ϕ*. Hutzler *et al.* [3] extended the steady-state foam drainage equation by adding a gas velocity *V* , and a variable gravity term to allow for the orientation of the flow with respect to gravity, taking account of the flow around the bend. This gas velocity is treated as *constant*, neglecting its small variation owing to the variation of the liquid fraction along the column. The model also treats the U-tube as ‘one-dimensional’ (i.e. the liquid fraction across a cross section of the tube is taken to be homogeneous). In real-world columns inhomogeneities can arise in the bend owing to gravity, leading to a liquid boundary layer and potentially a reduction in the liquid fraction in the outflow leg (modelling any such non-uniformity is not carried out in this work). Note that, in what follows, *x* is the height above the surface of the left or right liquid reservoir.

For the *left-hand side* (liquid flux and gas velocity in positive *x*-direction), the equation takes the form
2.5whereas for the *right-hand side*, where both liquid flux and gas velocity are now in the negative *x*-direction, we have
2.6

Finally, the boundary conditions must be set. As both left and right tube exits are in contact with liquid reservoirs, we may set the liquid fraction to some critical value *ϕ*_{c} (taken here to be 0.36, widely used in the literature [12,13], noting that key results are insensitive to this precise value),
2.7This is referred to, in this work, as the bottom boundary condition.

The top boundary condition used here for the U-tube requires simply that the liquid fraction at the top of the left column be equal to the liquid fraction at the top of the right. For a more detailed description of how this condition may be set, see appendix A.

Other fractionation set-ups require different boundary conditions, which will be presented as needed.

## 3. Performance of a fractionation column

Performance in a fractionation column can be considered in terms of enrichment (the ratio of concentration between the solution at the end of the column and the feed) or recovery (the fraction of the desired product that is recovered from the outflowing foam).

As an illustrative example of how our analytical model may be applied to column performance, we look at the recovery performance of the U-tube. As our model does not include chemical concentrations, we use *V*/*J* as a proxy for the efficiency of recovery, as follows.

Let us aim to increase the concentration of the surface-active components in the outflow. The amount of surface-active molecules that are carried through the column is related to the available area of film surfaces, and therefore to the gas velocity *V* (as increasing the gas velocity causes more bubbles to move through the column, and hence more surface-attached molecules). However, the foam has a finite liquid fraction and thus liquid that does not contain much of the component to be concentrated is also carried through the column, reducing the outflow concentration. We can therefore construct a simple metric for the performance of the column, the ratio of *V* (representing the amount of surface-active components carried through the column) to *J* (representing the other components carried through),
3.1

In order to compute this metric, we need to derive the detailed dependence of flux on the inflow leg length, *ϵ*(*L*_{l}).

### (a) Dependence of liquid flux on length of left leg

Our first step in calculating the performance metric is to determine *ϵ*(*L*_{l}), that is, to quantify the effect of *finite* lengths for the two legs of the inverted U-tube discussed in reference [3]. In the interests of clarity, some essential details are dispensed with in what follows, but all are included in the appendices. As previously mentioned, it can be shown that the length of the right-hand (output) tube *L*_{r} has a relatively small effect on the boundary conditions of the U-tube above a certain threshold length. This is discussed in appendix A, and the nature of the threshold is given by equation (A7). Accordingly, *L*_{l} is the primary object of the leg length analysis, with *L*_{r} implicitly set to a sufficiently large value. The effect of finite length is to increase *J* from the minimum value *J*_{0}, as shown in equation (2.2).

We can define an integral for the leg length based on equation (2.5) (using the liquid flux given in equation (2.1)), 3.2

After setting the bounds of integration (the top and bottom values of *ϕ*) to appropriate values and solving (with some assumptions, given with other details in appendix B), we arrive at an asymptotic expression for *L*_{l},
3.3where

Rearranging the relation in terms of *ϵ*(*L*_{l}), we arrive at an equation describing the effect of left leg length on liquid flux
3.4neglecting terms of higher order in 1/*L*_{l}. We can see that *ϵ*(*L*_{l})∝1/*L*_{l} (in leading order, for finite *V*)—increasing *L*_{l} will decrease *ϵ*(*L*_{l}), with as . Equation (2.2) may be compared with numerical solutions. This is shown in figure 2, with good agreement over a wide range of *L*_{l}.

We may also explore the behaviour of *ϵ*(*L*_{l}) numerically. Fixing the gas velocity *V* , we set for some *ϵ*_{i} and integrate the flux equation (equation (2.5)). We then take the point at which the liquid fraction *ϕ*=*ϕ*_{3} to be the leg length *L*_{l}(*ϵ*_{i}), as this is the value of the liquid fraction at the top of the right leg (*ϕ*_{3} is more precisely the asymptotic value of the liquid fraction in the right leg as ; see appendix A for a more thorough discussion).

We can use this method to see how as , shown in figure 3. The figure shows (on a log–log plot) the shape of *L*(*ϵ*^{2}). We observe a slope of −1/2, corresponding to *L*∝*ϵ*^{−1}, as expected for finite *V* .

Figure 4 shows how the liquid fraction profiles for the left leg vary with *ϵ*. Note that decreasing *ϵ* brings the finite left leg solution closer to the infinite left leg result, *J*(*L*)=*J*_{0}.

### (b) Results for performance metrics

As discussed above, the performance of the U-tube fractionation column depends on both the gas velocity and the length of the inflow leg, *L*_{l}. Plotting equation (3.1) (substituting equations (2.1) and (3.4)) and plotting for the same physical parameters used throughout the paper gives us figure 5.

Figure 5 contains a lot of useful information about our metric of efficiency. First, reducing gas velocity *V* *increases* efficiency. Second, for any given *V* , reducing *L*_{l} can lead to large performance drops. This can be seen more clearly in figure 6.

These figures allow us to make the following recommendations for fractionation column operation based on applying our performance metric *V*/*J* to the model U-tube system. First, the gas velocity *V* should be reduced as much as possible (taking into account desired output rates and physical limitations of the foam). Second, the length of the inflow leg *L*_{l} should be chosen carefully to ensure that operation is (for example) in the 95% regime or better. There may of course exist physical limitations on gas velocity and column size in real fractionation columns. Our model does not take into account changes in the liquid fraction (and hence liquid flux) arising from either coalescence or coarsening, and thus the metric of efficiency derived above may not be applicable to foams where these effects play a significant role.

## 4. Alternative fractionation columns

While we have had good success in our analysis of the U-tube, we have not yet generalized our approach to other types of fractionation columns. Any real industrial process is likely to diverge from our idealized model, and as such it is worthwhile to attempt to analyse a very different column in the same manner.

### (a) The skimmer

In certain fractionation applications, a further goal is to remove the surface-active components of the liquid phase as quickly as possible (e.g. removing proteins from aquaria). In those cases, ‘skimmers’ are often used. These devices collapse and remove foam from the top of a straight column, and with it liquid which is rich in surface-active molecules. A schematic showing a column with a skimmer can be seen in figure 7.

We now consider the case in which a skimmer defines the boundary condition at the top of a single vertical column. The skimmer removes foam, and with it liquid, at some rate *J*_{s}=*ϕ*_{t}*V* (per unit area), where *V* is gas velocity and *ϕ*_{t} is the liquid fraction at the top of the column. Conservation of mass then gives
4.1where *J*_{s} is the flux resulting from the action of the skimmer and *v*_{l} is the liquid velocity at the top of the column, or *v*_{l}=*V* at the top (similar boundary conditions have been applied to rising foam columns [15]).

The appropriate boundary condition is therefore 4.2

We proceed to calculate *J*_{s} in the same manner as before. We first determine the infinite leg solution, then add corrections for the effect of finite leg length.

In the infinite leg limit, *J*_{s}=*J*_{0} (see appendix C). The argument for *J*=*J*_{0} in the left leg of the full U-tube model holds again here. Then, considering a finite leg length *L*, we have a liquid flux of the form *J*_{s}(*L*)=*J*_{0}+*ϵ*_{s}(*L*)^{2} as before (the form of *ϵ*_{s}(*L*) is different from the U-tube case, and will be outlined below).

*ϕ*_{t} can, therefore, be written as
4.3

As in the full U-tube, the skimmer column solutions converge to an asymptotic value for the liquid fraction as . Setting *J*=*J*_{0} (and therefore *ϵ*_{s}=0) in equation (4.3) gives us the value of this asymptote, in this case *ϕ*_{t}=*V*/4*c*_{1}. For the U-tube, this value is *V*/2*c*_{1}, which implies that the liquid fraction at the skimmer *ϕ*_{t} is half the asymptotic value for the left leg of the U-tube (a relationship previously noted for overflowing single columns [1]). Figure 8 shows a numerically solved liquid fraction profile for a single column with skimmer.

We can derive a formula for *L*_{l}(*V*), following the same procedure as in §3*a* and appendix B. We arrive at an expression with a similar form to that in equation (3.3), as follows:
4.4where *μ*_{1} and *μ*_{2} are the same as in the U-tube, as only the top boundary condition has changed. The last term (*μ*_{3} in the U-tube) does depend on the top boundary condition, which changes from *ϕ*_{l}(*L*)=*ϕ*_{3} (for the full U-tube) to *ϕ*_{l}(*L*)=*ϕ*_{t}. It thus changes to *μ*_{4}, given by
4.5

A derivation for this quantity can be found in appendix D. From here, following the same procedure as used for the U-tube—solving equation (4.4) for *ϵ*_{s}—leading to equation (4.6) as follows:
4.6

Figure 9 shows numerical solutions for the skimmer column, along with theoretical predictions for *J*_{s}(*L*) (following from equation (4.4) as in the U-tube set-up). Good agreement can be seen over a wide range of leg lengths *L*.

### (b) Analysing fractionation modes using forced drainage

Throughout this work, we have looked only at the so-called simple mode of fractionation [16], with no liquid feed independent from the main solution reservoir. However, some fractionation columns incorporate a solution feed into the column, corresponding to what is called ‘forced drainage’ in accounts of the foam drainage equation. This can also be referred to as ‘deep foam washing’ [17] in a chemical engineering context.^{1} Such designs aim to ensure that all the surface-active components fully adsorb onto the bubble surfaces, increasing the enrichment performance of the column. An example is shown in schematic form in figure 10.

This may be analysed straightforwardly, building on the preceding results for the simple U-tube, which apply above the point where additional solution is introduced. We assume *L*_{2} (the length of the leg segment above the liquid addition) is large enough to take *J*=*J*_{0}, as discussed in §3*a*.

We may model the left (inflow) leg of such a column using numerical solutions based on equation (2.5), with *J* replaced by *J*(*z*),
4.7where
This differential equation can be solved numerically, resulting in liquid fraction profiles for the left leg. The right (outflow) leg solution remains identical to that of the simple U-tube. An example numerical solution is shown in figure 11, with a solution for the same column leg without forced drainage for comparison.

The additional flux *J*_{in} can be chosen to increase the liquid fraction in the lower part of the leg (within limits of stability). Controlling the liquid fraction in this manner may allow the column adsorption efficiency to be improved. In order to gain the most from the additional liquid flow, the length *L*_{1} must be sufficiently long to ensure complete adsorption. Large increases in enrichment are possible, as reported for forced drainage experiments with reflux. Reflux in a fractionation context describes columns where the foam outflow (at a higher concentration than the original solution) is collected and added in the manner discussed above [18].

## 5. Conclusion

We have presented an analysis of the effects of varying leg length on the operational efficiency of foam fractionation. Our results describe how the *ϵ* term in the liquid flux equation *J*=*J*_{0}+*ϵ*^{2} varies with left leg length. Numerical simulations were carried out, and agree closely with our analytical results.

We also provide an alternative boundary condition in the form of the single column with skimmer. The analytic approach used to examine the U-tube is shown to be valid in the skimmer set-up, and numerical results again agree closely with predictions.

When building a model for processes such as fractionation as used in industry and chemical engineering, it is important to consider real-world uses. In this work, we have attempted to both keep the model and analysis as widely applicable as possible, while keeping in mind constraints that exist in industry.

We do not yet consider chemical processes that could affect the performance of a fractionation column in our analytical models. Of particular interest to column design may be the adsorption time *τ* (how long surface-active molecules take to fully bind to the foam films). Combined with the gas velocity *V* , we can write a minimum operational leg length *L*_{a}=*τ*×*V* . Below this length, the efficiency of the column is reduced, as not all surface-active molecules are carried through the foam, but instead drain out of the foam.

## Funding statement

Research is supported the European Space Agency (MAP grants no. AO-99-108:C14914/02/NL/SH and no. AO-99-075:C14308/00/NL/SH) and COST Action MP1106 Smart and Green Interfaces.

## Acknowledgements

We thank David A. Whyte for the Mathematica calculations used for the left leg integral, and Sara McMurry for discussions related to the contour integration method. We also thank Peter J. Martin for several helpful conversations, and for pointing us towards the work of Desai and Kumar.

## Appendix A. Variation of liquid fraction in the right leg

In order to fix an upper boundary condition for the top of the left leg, we first consider the right leg. (In the following, *x* is the upwards vertical coordinate, with *x*=0 at the liquid reservoir surface). Hutzler *et al.* [3] give the equation for the variation of liquid fraction as
A1

Writing *J*=*J*_{0}+*ϵ*^{2} (with *J*_{0}=*V* ^{2}/4*c*_{1}), we obtain
A2Here, *ϕ*_{3} and *ϕ*_{4} represent the roots of the expression in parentheses in equation (A1) in the limiting case *J*=*J*_{0}, and are
A3and
A4

In the limit of infinite leg lengths, as . For large enough leg length *L*_{r}, we may assume Δ*ϕ*=(*ϕ*_{r}−*ϕ*_{3}) is small as *ϕ*_{r} asymptotes to this value, and can therefore approximate equation (A2) by
A5neglecting terms of order *ϵ*^{2}. Then,
A6

Accordingly, Δ*ϕ* decreases exponentially with height *x*, and we can write a decay length as
A7Therefore, by increasing leg length from *L*_{d} to 2*L*_{d} the liquid fraction is decreased by a factor 1/*e*. Similar decay lengths have been found for analysis of froth floatation [14]. The error in this length is of order *ϵ*^{2}. For *L*_{r}≫*L*_{d}, this implies that, at the top of the leg,
A8and that, for such a case, it would be a good approximation to take *ϕ*_{r}=*ϕ*_{3}, therefore fixing the boundary condition for the top of the left leg at the same value that was found for the limit of infinite legs. A comparison of a full numerical solution and the simple exponential may be seen in figure 12.

## Appendix B. Full solution of the finite left leg integral

We start from an integral for *L*_{l}, following from equation (2.5),
B1where *ϕ*_{t} and *ϕ*_{b} are the top and bottom values for the liquid fraction, respectively, and
B2We evaluate this integral by making a drastic approximation, which is then corrected: we set the lower boundary value for *ϕ*_{l} to infinity, and the upper bound to zero. Previous work on fractionation columns similar to our model's left leg have dealt with the case where *ϕ*_{1} and *ϕ*_{2} are real [14] (corresponding to ); we must here handle complex values for these quantities.

Substituting equation (2.1) into equation (B2) gives B3

Expressing these complex values of *ϕ* using Euler's notation gives
B4where *a* is the magnitude of the complex number, and *b* is the argument (or phase). These quantities are given by
B5and
B6

Substituting equation (B4) into equation (B1) results (after some straightforward manipulation) in B7where .

Integration using Mathematica yields B8where .

Alternatively, we may evaluate the integral quite neatly by contour integration (the method of residues), using the contour illustrated in figure 13. The cut along the real axes is such that has a positive real root on the upper side, and a negative real root on the lower side. The denominator of the integrand has two complex roots at B9

The residues associated with each are (taking into account the nature of square roots in a complex plane) B10and B11

Our contour integral (from Cauchy's residue theorem, with the above caveats taken into account) is therefore (noting that such a contour results in twice the correct integral) B12

Evaluating this expression leads to precisely the same result found by Mathematica shown in equation (B8). Gathering all the constants into a single *μ*_{1}, we obtain
B13

Note that, for given *V* , the limit of infinite *L*_{l} corresponds to *ϵ*=0, in accord with the previous theory.

We can now proceed to examine the leading approximation for the two end corrections that are required when the proper boundary conditions are reintroduced. We correct equation (B1) for the correct bottom boundary condition, *ϕ*_{l}(0)=*ϕ*_{c}, by adding
B14to the integral of equation (B1).

Neglecting terms of order *ϵ*^{2} and taking into account that and we may approximate the correction as
B15
B16

Finally, we add a further correction for the top boundary condition, *ϕ*_{l}=*ϕ*_{3}, as follows (again neglecting terms of order *ϵ*^{2}):
B17

In this case, in the range of integration, hence B18 B19 B20

Summarily, we have arrived at
B21where
This is a good approximation for the left leg length, provided that the right leg length *L*_{r}≫*L*_{d}, where *L*_{d} is given by equation (A7) (see appendix A), and hence that the boundary condition used for the top correction is valid.

## Appendix C. Generality of the limiting case

It is worth noting that the initial approximation (integrating from 0 to ) for *L*_{l} in the previous section results in an upper bound for *L*_{l}, given by . As the corrections for the top and bottom bounds of integration given by *μ*_{2} and *μ*_{3} are negative, . Therefore, if , irrespective of any boundary conditions chosen.

This implies that, for a wide range of boundary conditions (indeed, almost all), the steady-state solution must have *J*=*J*_{0} in the limit . Indeed, numerical simulations for the skimmer boundary conditions bear this out, recovering this relationship for sufficiently long columns. This limiting result was presented by Hutzler *et al.* [3], derived for the specific boundary conditions of the U-tube, but appears to be far more general (previous works have also observed such generality [14]). Outside of the limiting case of infinite leg lengths, the boundary conditions are not equivalent and greater care must be taken. Numerical simulations have shown that the limiting case is approached relatively quickly (see, for example, figure 4).

## Appendix D. Top boundary correction for skimmer

In the case of a skimmer, the top boundary condition changes from *ϕ*_{l}(*L*)=*ϕ*_{3} to *ϕ*_{l}(*L*)=*ϕ*_{t},
D1where . Therefore, over the range of integration. If we assume that , we can calculate an approximate correction for the top boundary condition,
D2
D3
D4
D5
D6

## Footnotes

- Received September 20, 2013.
- Accepted January 21, 2014.

- © 2014 The Author(s) Published by the Royal Society. All rights reserved.