Convective instability in sheared foam

Abstract This work provides evidence that anisotropic drainage in sheared foam is at the origin of convective instability in very long foam channels. Convective instability occurs in foam under forced drainage when a critical liquid fraction is exceeded. Liquid spontaneously accumulates at one side of the channel. The weight imbalance induces convection rolls in the foam. Experiments in a very long vertical foam channel demonstrate that the critical liquid fraction is smaller than in previous findings by a factor of five. The critical liquid fraction depends on both the channel length and the inhomogeneity of the liquid feed. Well below the critical liquid fraction, a static, elastic shear deformation of the foam structure occurs. At the critical liquid fraction, initial steady convection rolls are located at the lower region of the channel and expand as the liquid fraction further increases. Combining the drainage equation with both the elastic response of the foam and a model for anisotropic drainage, a critical liquid fraction for the growth of an initial liquid imbalance is derived analytically, which corresponds very well to experimental findings. Numerical simulations of the drainage equation and the elastic response of the foam reproduce these experimental and analytical findings.

of the foam, which in turn depends on the liquid fraction. Parameters η f , ρ f and g in (1.1) denote the fluid dynamic viscosity, density and gravitational acceleration, respectively. If a small flow of liquid is added homogeneously to the top cross-section, a static and homogeneous liquid fraction is established in the foam. However, when the flow of liquid exceeds a critical value, an inhomogeneous liquid fraction results and convection rolls set in. These convection rolls bear similarities to Rayleigh-Bénard cells occurring in thermally stratified liquids when the density gradient, respectively the Rayleigh number, exceeds the critical value (Chandrasekhar 2013). But in foam drainage, no vertical density gradient is present when the convection rolls occur. Consequently, the mechanism behind these rolls is still a matter of debate. When liquid is added to the top region of an aqueous foam column, it drains downward through the network of Plateau borders. This effect is described by the drainage equation (Cantat et al. 2013): It links the local liquid fraction φ l and the superficial liquid flow rate q with the properties of the liquid (density ρ f , viscosity η f , surface tension γ ) and those of the foam (permeability α, bubble radius R b ). The geometric constant δ b = 1.74 results from the relation between the Plateau border cross-section and the liquid fraction. The drainage equation (1.2) considers two effects: the vertical flow due to gravity g and a diffusion of the liquid fraction due to capillary pressure gradients. The diffusion term balances inhomogeneities in the liquid fraction. Experiments in steady drainage (Hutzler, Weaire & Crawford 1998;Vera, Saint-Jalmes & Durian 2000;Hutzler et al. 2007) have investigated whether the foam structure remains static. They have observed that if the liquid fraction exceeds a certain critical value φ crit , a stable convection roll is formed, known as convective instability (CI). This is linked to an inhomogeneous liquid distribution over the channel cross-section, which cannot be explained by the drainage equation. The formation of this convection roll biases experiments in foam at a high liquid fraction, because it disturbs the assumption of homogeneous liquid fraction.
Although detailed studies have not been carried out in industrial foam and froth applications, CI can be anticipated to have a severe impact. To illustrate this, we take deep-froth flotation as an example, which is an important step in the beneficiation of valuable mineral particles, e.g. copper-bearing ores. Frequently, wash water is added to the particle-laden foam (froth) to remove hydrophilic contaminants, such as undesired entrained gangue particles. It can be presumed that the CI reduces the efficiency of wash water addition. As a result, the gangue particles are not completely washed out. Thus, the grade of the particle concentrate in the froth is reduced which has a significant economic impact.
The critical liquid fraction φ crit for the onset of CI has been determined to depend on the bubble size (Hutzler et al. 1998) and can be reduced by tilting the channel (Cox et al. 2006). Convective instability results in increased effective drainage and can cause bubble sorting (Hutzler, Weaire & Shah 2000). However, the actual onset mechanism is still a matter of debate. Embley & Grassia (2006) assume that a sudden scarification of individual Plateau borders leads to a localization of the drainage flow. Weaire & Hutzler (2003) demonstrate that dilatancy combined with elastic deformation of the foam can cause a growing imbalance in the liquid fraction. Neethling (2006) predicts that anisotropic drainage in sheared foam causes an imbalance in the liquid fraction.
Using a combined experimental, analytical and numerical approach, we revisit the onset of CI in detail and provide a conclusive explanation of the mechanism underlying CI. To facilitate understanding, the stages of the instability which are passed in a typical experimental run are illustrated in figure 1. Starting from a dry state (I), the drainage flow rate is increased in small steps, leading to a stepwise increase in the liquid fraction, which causes the foam to expand by the volume of additional liquid content (II). At each drainage level, a static foam structure is reached. Above a critical liquid fraction of φ crit ≈ 0.65 % the deformation becomes inhomogeneous (III). However, the foam structure still reaches a static state for each drainage level. In this static state the shear rate is zero, but the shear angle with reference to the initial stage I is non-zero. The inhomogeneity and shear angle increase with increasing drainage flow rate. Above a second critical liquid fraction of φ CI ≈ 1.05 %, the critical yield strain is exceeded near the bottom of the column and a steady convection roll is formed in the lower region. In the upper region, the foam maintains its static inhomogeneous deformation (IV). With further increase in drainage flow rate, the convection roll covers more and more of the channel height but remains steady for a constant drainage flow rate.
We claim that the primary instability is the transition from homogeneous (II) to inhomogeneous (III) static deformations of the foam structure, which is linked to an inhomogeneous liquid distribution. In stage III, the inhomogeneity grows in the vertical direction. If the column is long enough and the growth rate of the inhomogeneity with respect to the vertical position is high enough, the yield stress is exceeded at some vertical position. Below that position the convection roll sets in, corresponding to stage IV. Thus, the CI appears as a secondary instability following the primary one, the inhomogeneous drainage.
In this paper we present experimental findings on the static deformation of the foam structure and the steady velocity field in the case of CI. Our findings are compared with a linear stability analysis for the growth of inhomogeneities of liquid fraction. This analysis includes the anisotropic drainage, as predicted by Neethling (2006). The identified unstable mode marks the transition from stage II to III. Numerical simulations of the drainage combined with the elastic deformation of the foam and anisotropic drainage (Neethling 2006) reproduce our experimental findings. Finally, we demonstrate the dependency of the critical liquid fraction φ CI for the onset of CI on the initial liquid inhomogeneity as well as on the channel length.

Materials and methods
Experiments are carried out in a vertical foam channel with an effective length of 990 mm and a cross-section T × B = 30 mm × 100 mm. Cox et al. (2006) have demonstrated that a small tilt of the column can lead to CI. Thus, the tilt angle was controlled to be less than 0.05 • . Bubbles are generated at the bottom of the channel by a tube with 20 holes, each 0.5 mm in diameter, submerged in the surfactant solution. The tube is loaded with pressurized air, generating bubbles with radii of R b = 2 ± 0.2 mm at a flow rate of 2 l min −1 .
The surfactant solution is deionized water with 35 mM sodium dodecylsulphate (SDS), generating stable foam. When keeping the foam for 6 hours at 0.5 % liquid fraction, the measured change in bubble size was less than 20 %. Nevertheless, fresh foam is generated for each run, which take less than 1 hour each. At the top of the foam column, a steady liquid flow is imposed through four porous hollow cylinders with outer diameters of 16 mm, yielding a forced drainage configuration similar to that of Leonard & Lemlich (1965). The volumetric liquid flow rate Q is linked to the drainage flow rate q = Q/(B T) by the channel cross-section (B T). The drainage flow rate is also known as the superficial drainage velocity or liquid flux, and equals the volume of liquid that passes through a certain cross-section of the foam channel within a certain time.
Each cylinder is connected to one channel of a four-channel peristaltic pump. This enables each cylinder to be controlled independently. That is, a defined inhomogeneous distribution can be imposed by pumping liquid only through one, two or three of the four cylinders. The notation (oxxx) means that only the leftmost cylinder is charged with liquid. The liquid is taken from the liquid reservoir below the foam. In that way, the total volume of liquid and gas remains constant and independent of the volumetric liquid flow rate.
The liquid fraction is observed with four pairs of electrodes (5 mm × 45 mm electrode area) at a horizontal cross-section 20-65 mm below the porous cylinders (see figure 2a,b), yielding the horizontal profile of the liquid distribution φ( y). From the four measured values φ 1 , . . . , φ 4 , the average liquid fraction φ m = 0.25(φ 1 + φ 2 + φ 3 + φ 4 ) and the first are derived. In our measurements, we never achieved a perfectly homogeneous distribution of liquid fraction. Even when all four porous cylinders are charged, we find M φ of the order of 1 mm. This corresponds to an average horizontal shift of the liquid content by only 1 % of the channel width. We spent a lot of effort to further reduce the inhomogeneity. For example, we measured the flow rate for each channel from our four-channel peristaltic pump. Also, we switched the channels. But, even if we maintained everything constant and only created fresh foam, we found the inhomogeneity changing and sometimes switching to the other side. So we believe that small initial values of M φ can be considered as a random distortion which are always present in a random foam structure. The foam structure inside the transparent channel is observed with a camera taking images with 1920 px × 1200 px at 30 frames per second, yielding a spatial resolution of 3.3 px mm. Backlight illumination with a light-emitting diode panel permits exposure times down to 20 ms. Due to the illumination, the captured foam structure is a cumulative image in the z direction. Consecutive images are analysed with the particle image velocimetry (PIV) algorithm implemented in DaVis 8.0 software. Particle image velocimetry (Adrian, Adrian & Westerweel 2011) compares consecutive images to measure the shift of the foam structure between them.
The PIV analysis is employed in this work in two different cases. In case 1 (stage IV) the network of Plateau borders and vertices moves, while in case 2 (stages I, II and III), the network is static.
In case 1 the shift of the network structure between consecutive images is multiplied with the frame rate, yielding the velocity distribution v of the foam in the x-y plane. This velocity refers to the average velocity of the foam structure, i.e. to the average velocity by which the network of Plateau borders and vertices moves. In this case convection rolls are recorded (stage IV). Typically, the network velocity reaches a steady, non-static state in our experiments. In the case of low velocities, increments of up to 20 frames were used for the PIV analysis. In addition, averaging over 10 image pairs without outlier reduction is employed. This algorithm yields the velocity distribution v(x, y) in the x-y plane with a 5 mm spatial resolution. The velocity of the liquid inside the Plateau borders is not accessible. Due to the backlight illumination, the velocity is automatically averaged in the z direction. For spatially resolved three-dimensional measurements of the foam velocity, other techniques such as ultrasound Doppler velocimetry (Nauber et al. 2018) or radiographic particle tracking (Lappan et al. 2020) are required.
In case 2, the network is static everywhere in the column. In this case, the PIV analysis mentioned above yields vanishing velocity v(x, y) = 0. However, even when the network is static, it has deformed compared with the initial state (stage I) of zero drainage flow. We increase the liquid flow rate stepwise and find at each step a different static network. When we then switched off the drainage flow, the network goes back to its initial shape. It is like a beam balance. For each set of weights (drainage flow) a certain static angle (shear angle) occurs. In order to analyse this static elastic deformation U(x, y) of the foam (figure 3), one image is taken from each drainage flow rate level, and these are combined and fed into the PIV algorithm. By this means, the static elastic shift of the foam network between consecutive levels of drainage flow rate can be computed (see figure 2d). This only works if the foam structure does not undergo a flowing movement, because in the case of flowing foam, the shift between consecutive drainage levels is too large to be tracked.

Stages II and III: static elastic deformation
Each experimental run is prepared by creating fresh foam and letting it rest for 3 minutes. Then, the pump is started, feeding a constant volumetric flow rate Q to the top of the column, resulting in a constant, downward drainage flow rate q = Q/(B T), which is identical to the vertical superficial liquid velocity. In this section, only homogeneous inflow (oooo) is considered. The drainage flow rate q is increased in small steps of 20-50 μm s −1 , giving the foam column 3 minutes after each step to achieve a steady state. This experimental procedure, which avoids undesired drainage fronts, is sketched in figure 1(b).
In stages II and III, the stepwise increase in the drainage flow rate causes a corresponding stepwise static elastic deformation U of the foam. Figure 3(b) shows the measured static vertical displacement U x of the foam at various drainage flow rates. The corresponding initial liquid fraction distributions are depicted in figure 3(a). Figure 3( c) shows U x at selected horizontal lines for different drainage flow rates.When the drainage flow rate is switched off again, the foam relaxes back into its initial shape at stage I. In stage II below q = 65 μm s −1 , the displacement U x is fairly constant in the y direction but increases in the x direction. The relative displacement U x /x is of the order of 0.5 %. This homogeneous displacement is due to the change in the liquid volume inside the foam while the total volume of gas and liquid is kept constant. For example, the column with  a height of 1000 mm displaces the foam-liquid level by 5 mm when the liquid fraction is increased by 0.5 %. At around q = 65 μm s −1 and φ l ≈ 0.6 %, the transition to stage III occurs. A horizontal gradient is seen in the vertical displacement, which is equivalent to a shear deformation ε xy of the foam. This shear ε xy is predicted by Neethling (2006) to generate a horizontal deflection q y,aniso of the vertical drainage flow rate q x , which, in turn, increases the liquid imbalance and shear. This mechanism is discussed in § 5. The liquid imbalance grows with increasing drainage flow rate until, above 160 μm s −1 , CI sets in (see figure 4, predicated on the same run as figure 3).

Homogeneous inflow
The onset of CI is recognized as when foam velocities larger than 0.5 μm s −1 are detected in a steady state at any point of the foam column. Similarly to Hutzler et al. (1998), the critical liquid fraction and drainage flow rate at the onset of the CI are determined. Figure 4 documents the velocity distribution of the foam in a steady state at certain drainage flow rates. In this particular run, the critical drainage flow rate for the onset of CI is between q = 160 and 179 μm s −1 , corresponding to φ CI = 1.05 %. At q = 160 μm s −1 , the foam structure remains static. At q = 179 μm s −1 , the first convection roll is detected at the bottom region, constituting the onset of CI. This roll is steady. It does not grow in size or magnitude as long as the drainage flow rate remains constant. As q increases further, the foam velocity in the convection roll increases, as does the extension of the roll in a vertical direction. At drainage flow rates above approximately q = 500 μm s −1 , the roll 0 covers the entire height of the channel. This relation between the drainage flow rate q and the roll extension supports our mechanism in figure 1: the inhomogeneity of liquid fraction in stage III requires a certain vertical length to grow and reach a critical imbalance. With higher liquid fraction, the inhomogeneity grows faster with respect to the vertical distance.
In the cross-section, where a critical level of inhomogeneity (i.e. the yield strain; see § 6) is exceeded, the steady movement of foam sets in. The resulting convection roll covers the channel height from this critical cross-section to the bottom. At higher liquid fraction the inhomogeneity grows faster and, thus, the convection rolls sets in more close to the top.

Inhomogeneous inflow
In order to further investigate the growth of the inhomogeneity in the liquid fraction with respect to the vertical position, an inhomogeneous liquid fraction is deliberately imposed by feeding only some of the four porous cylinders. The corresponding liquid moment M φ according to (2.1) is derived from the conductivity sensor below the porous cylinders. Figure 5(a) shows how the drainage flow rate q CI and liquid fraction φ CI for the onset of CI depend on the liquid moment. In the case of high initial inhomogeneity, q CI and φ CI can be as low as 70 μm s −1 and 0.65 %, respectively (dash-dotted line in figure 5a). These values are close to the transition from stage II to stage III in homogeneous drainage. This again supports our explanation: if the initial inhomogeneity is very high already, there is only a small increase in inhomogeneity required to reach an imbalance that exceeds the yield stress. In the given channel length, this can be achieved by a small growth rate with respect to the vertical position. As demonstrated in § 5, a liquid fraction just above the critical liquid fraction φ crit will cause such a slow growth of the inhomogeneity with respect to the vertical direction, because it is fed by the anisotropic drainage. However, in the case of a lower initial inhomogeneity (symbol 'o' in figure 5a), even drainage flow rates up to 250 μm s −1 and 1.2 % liquid fraction can result in a static foam structure (stage III). In this case, the liquid inhomogeneity does not grow fast enough to reach a critical imbalance in the provided channel length. Now, the length of the channel is reduced. Figure 5(b) documents the critical drainage flow rate for different channel lengths and initial inhomogeneities. The error bars denote the variation between several runs, each with fresh foam. If the initial liquid fraction exhibits strong inhomogeneity (ooxx), the reduced length has only a small influence on the critical drainage flow rate. Similar to high initial inhomogeneities in figure 5(a) the drainage flow rate has to be just high enough to maintain the inhomogeneity. No large channel length is required to reach a critical imbalance. However, for an initially homogeneous liquid fraction (oooo), the critical drainage flow rate increases significantly the shorter is the channel. Similar to small initial inhomogeneities in figure 5(a) the inhomogeneity requires a sufficient channel length to reach a critical imbalance. If the channel is shorter, the inhomogeneity has to grow much faster to reach a critical imbalance. To accomplish that, the liquid drainage flow rate has to be substantially higher.

Stability analysis
The transition from stage II to stage III marks the primary instability in foam drainage.
In stage II small inhomogeneities in the inflow of liquid, i.e. small M φ , are damped and homogenized. In stage III such inhomogeneities grow in the vertical direction. In order to analyse whether inhomogeneities of the liquid fraction grow with respect to the vertical direction, a stability analysis is performed. It is based on the drainage equation (1.2), which is repeated for the sake of a coherent representation: The foam permeability α also depends on the liquid fraction and the interface mobility. The surfactant used here, SDS, yields a mobile interface. Thus, in the present investigations the dependency is α ∝ φ 3/2 l (Lorenceau et al. 2009). Consequently, the term α/φ 3/2 l in (5.2) becomes a constant and is replaced by K α : Surfactants leading to an immobile interface are considered in appendix B.
To account for the above-mentioned effect of anisotropic drainage, Neethling (2006) suggests using an additional term: A vertical drainage flow rate q x results in a horizontal drainage flow rate q y , which is proportional to the local shear deformation ε xy . Of course, the same mechanism would also cause an additional vertical drainage flow when a horizontal drainage flow is present. However, since gravitationally driven vertical drainage typically exceeds capillary-driven horizontal drainage, the reverse effect is neglected here. The extended drainage equation reads In order to derive the growth of an initial horizontal inhomogeneity of the liquid fraction, the normal mode approach is applied to (5.4). For that purpose, the liquid fraction φ l (x, y) is decomposed into a constant value φ 0 and small perturbations φ k (x) with wavenumbers k which evolve in the vertical direction (cf. figure 6): (5.5) It is important to note that no temporal dependency is incorporated. This is motivated by our experimental observations, where the elastic deformation and the liquid fractions reached a steady state. As already mentioned, even in the case of highly symmetric liquid inflow, we never succeeded in generating a perfectly homogeneous liquid content in the conductivity sensor below the porous cylinders. Randomly, we found that the liquid content was higher at the left or right side of the channel. We assume that the asymmetry of the static foam structure which is in contact with the cylinders always generates small inhomogeneities. These slight inhomogeneities then remain steady over time. Consequently, in our experiments we always had steady but slightly asymmetric inflow conditions, which is reflected by the ansatz function. While the average foam weight is compensated for by a vertical pressure gradient, the perturbations lead to an imbalance in the gravitational force: (5.6) For the stability analysis, we assume an infinitely long (high) channel. Consequently, this imbalance must not accumulate over the channel height. In each horizontal slice, Figure 6. Sketch of the considered system for linear stability analysis. the gravitational imbalance is compensated for by a local shear stress. The vertical force balance on a foam volume stretching between the left wall and the position y =ỹ yields At the vertical boundaries of the channel, a liquid film results in negligible static friction between the foam and wall. Consequently, the shear stress has to fulfil the boundary conditions The shear modulus G of foam depends on the liquid fraction. It can be estimated (Mason, Bibette & Weitz 1995) by (5.10). In the present case, critical liquid fractions of φ 0 ≤ 1 % are considered. In this small range the influence can be neglected, yielding (5.10) Combining (5.10) with (5.9) leads to the shear strain: Considering a steady liquid distribution, the term ∂φ l /∂t in (5.4) equals zero. Feeding the ansatz (5.5) into (5.4), neglecting the influence of capillary pressure on the vertical liquid transport and assuming that only small perturbations φ k φ 0 occur leads to For a detailed derivation we refer the reader to appendix A. This ordinary differential equation is solved by (5.13) The exponent A represents the growth rate with respect to the vertical position x and is given by (5.14) Equations (5.13) and (5.14) show that small imbalances in φ l grow exponentially with respect to the vertical distance x provided that A > 0. Inspecting (5.14), this is the case for sufficiently small horizontal wavenumbers k. The smallest possible wavenumber in a channel that fulfils the boundary conditions (5.8) is k min = π/B. At the critical point, A = 0 and the liquid fraction φ 0 = φ crit . This yields the stability criterion for drainage in an infinitely long vertical foam channel of width B: above which initial inhomogeneities grow exponentially. In the present case, we obtain φ crit ≈ 0.74 %. This value corresponds well with the experimental data for the transition from stage II to stage III. The critical liquid fraction for the onset of the shear deformation of the foam in figure 3(c) is φ crit ≈ 0.6 %.

Simulations
In order to investigate the interaction of stress, deformation and liquid distribution at the transition from stage II to stage III, a two-dimensional, phase-averaging numerical simulation of the experimental set-up is employed. The code computes the unsteady drainage equation, (5.4), in a finite-volume discretization with Euler-explicit time integration. In each time step, the simulation solves iteratively for the linear-elastic deformation U and the stress tensorτ of the foam with zero strain in the third dimension:

4)
∇τ + f g = 0, (6.5) balancing the gravitational load f g from the liquid fraction contained. The shear modulus G depends on the local liquid fraction according to (5.10). The Poisson ratio ν for  incompressible foam would equal 0.5 and the compression module K ν would be infinite.
To avoid numerical problems resulting from an infinite compression module K ν , ν is set to the artificial value of 0.49. It has been found that small variations of ν have negligible effects on the results. The resulting local shear strain ε xy is fed back into the drainage equation. The pertinent boundary conditions are given in table 1. At the sidewalls, zero vertical stress, zero horizontal liquid flux and zero horizontal strain are imposed. At the top wall, zero strain and a prescribed vertical drainage flux q x ( y) = q 0 (1 + q cos( yπ/B)) with inhomogeneity 0 < q < 1 are imposed. At the bottom, 20 % liquid fraction, zero horizontal stress and a hydrostatic vertical stress as a function of the displacement U x are imposed. Figure 7(a) shows the numerical results for a case of unstable drainage marked by a green star in figure 7(b). The initial distortion of the liquid fraction φ l grows in the downward direction (figure 7a). The corresponding weight imbalance causes an inhomogeneous downward displacement U x . Due to the hydrostatic boundary condition, the inhomogeneity of displacement is reduced close to the bottom (figure 7a). The horizontal gradient of U x corresponds to a shear strain ε xy (figure 7a), which is highest in the centre of the channel. The shear drives a horizontal drainage flow q y = q aniso − q cap (figure 7a), which is the anisotropic flow q aniso according to (5.3) minus the capillary-driven drainage flow q cap according to the second term on the right-hand side of (1.2). The horizontal drainage flow q y is positive in regions of high shear, feeding the imbalance. Close to the top and bottom, the boundary conditions inhibit the shear deformation, and a negative horizontal drainage flow occurs due to capillary forces. The ratio of the local strain ε xy to the local critical yield strain ε y = 0.3(0.36 − φ l ) (Saint-Jalmes & Durian 1999) is highest in the lower part of the channel (figure 7a), where experimental investigations have found the onset of CI ( figure 4). Now, the initial inhomogeneity and total drainage flow rate are varied. If the local strain ε xy in the converged simulation exceeds the theoretical limit ε y at any place, the case is considered to show yielding and CI (stage IV). Figure 7(b) depicts the limit for the onset of CI for different channel lengths and initial moments of liquid distribution. It shows very good agreement with the experimental findings in figure 5. For high moments of the liquid fraction M φ , CI occurs at φ CI = 0.65 % liquid fraction and q CI = 70 μm s −1 drainage flow rate, regardless of the channel length. A small moment of 1.5 mm in a long channel causes CI at approximately φ CI = 1.7 % liquid fraction and q CI = 250 μm s −1 drainage flow rate. For shorter channels, the critical values are increased to φ CI = 2.7 % liquid fraction and q CI = 500 μm s −1 drainage flow rate.

Discussion
In this work, we have proven the mechanism for the onset of CI. Our experimental and numerical findings support the concept of an instability in forced drainage due to anisotropic drainage which was initially predicted by Neethling (2006). The dependencies of the critical liquid fraction φ crit on bubble radius, surface tension, gravity and channel size are identical to the findings of Neethling (2006) for the assumption that the critical liquid fraction is negligible compared with the jamming point liquid fraction. However, the stability criterion that Neethling derived differs significantly from ours. In our case it would produce a critical liquid fraction of approximately 0.2 %, which is considerably lower than our experimental findings. We extend his theory by adding the concept of a growing inhomogeneity with respect to the vertical direction and by describing the CI as a secondary instability superimposing on that.
Initial inhomogeneities due to the liquid addition cause a horizontal gradient in the local gravitational force. This gradient is compensated for by a shear deformation of the foam structure. The shear deformation deflects the vertical drainage flow in the horizontal direction, increasing the imbalance. At the same time, capillary forces tend to reduce the liquid imbalance. If a critical liquid fraction is exceeded, the horizontal deflection exceeds the capillary effect and the initial inhomogeneity grows with respect to the vertical distance. This proposed mechanism is supported by the agreement of the critical liquid fraction from experiments, simulations and stability analysis, by the measurement and simulation of the distribution of inhomogeneous vertical displacement and by the dependency of the occurrence of CI on channel length and initial inhomogeneity.
In the experimental results (figure 3) and in the simulation (figure 7), we have observed the occurrence and vertical growth of such a shear deformation above a critical liquid fraction. This process does not necessarily involve a steady flow of the foam. Only if the channel is long and wide enough does the growing imbalance at some point exceed the yield stress of the foam and CI occurs below that point. However, the onset of CI is only a secondary instability of the primary instability, consisting of the shear deformation at φ l > φ crit . Our experimental and numerical observations of the CI have demonstrated that the onset of CI takes place close to the bottom of the channel (figure 4), where the imbalance is highest. Moreover, the liquid fraction φ CI for the onset of CI in the case of an inhomogeneous introduction of liquid in figure 5(a) is φ CI ≈ 0.65 %, which is very close to φ crit ≈ 0.74 %. In line with the discussion above, a reduced channel height increases the critical liquid fraction for the onset of CI (figure 5b). Then again, a strong inhomogeneity in the liquid introduction reduces the channel height required to reach the onset of CI (figure 5b).
From our point of view, the presented findings on elastic deformation also rule out other existing explanations of CI. The idea that dilatancy might cause CI is disproved by the fact that the strongest shear is found in the centre of the channel (figure 3b). Therefore, dilatancy would cause the accumulation of liquid in the centre and not close to the channel side. Also, the idea of a sudden scarification of individual Plateau borders does not explain why the onset of CI is close to the bottom of the channel and why the size of the CI grows with increasing drainage flow rate.
The presented experiments on both the CI with strongly inhomogeneous liquid inflow and on the inhomogeneous elastic deformation together with the elastic simulation all consistently produce a critical liquid fraction of approximately φ crit = 0.65 %. This critical liquid fraction marks the transition from stage II to stage III. In stage II, the static foam stretches homogeneously and initial inhomogeneities in liquid fraction are damped with respect to the vertical direction. In stage III, the static foam stretches inhomogeneously and initial inhomogeneities in liquid fraction grow with respect to the vertical direction. Due to the high scattering between experimental runs there is a certain uncertainty in the critical value resulting in an interval of φ crit ≈ 0.6 % to 0.65 %. Also, there is an uncertainty in distinguishing the different stages. In particular, it is difficult to assign the transition from stage III to stage IV to a fixed value of φ CI , because in some cases the foam network starts to flow and comes to a halt again. In figure 5 in some cases of high initial inhomogeneity, the CI sets in at liquid fractions φ CI even below the theoretical critical liquid fraction φ crit . This might be due to such uncertainties or due to nonlinear effects under high imbalance and liquid fraction gradients.
Carrying out a linear stability analysis, we have established a stability criterion above which initial inhomogeneities in the liquid fraction grow with respect to the vertical position. The resulting critical liquid fraction of φ crit = 0.74 % is close to the experimental result of φ crit ≈ 0.6 % to 0.65 %, supporting the proposed mechanism. Our critical liquid fraction is a factor of five below the critical liquid fraction reported in Hutzler et al. (1998) for the onset of CI in cylindrical channels. One reason for the discrepancy is our significantly larger channel width B, which reduces the critical liquid fraction according to (5.15). Also, Hutzler et al. (1998) investigated the critical liquid fraction φ CI for the onset of CI (transition from stage III to stage IV). By comparing figures 3 and 4 for homogeneous water inflow φ CI ≈ 1.05 % is demonstrated to be larger than φ crit ≈ 0.65 % for the transition from stage II to stage III. And finally, our channel is one order of magnitude longer, giving the inhomogeneity more distance to grow. Thus, smaller liquid fractions can lead to the onset of CI. The influence of the channel length on φ CI is documented in figures 5(b) and 7(b).
We do not believe that the deviations between our critical liquid fraction and previous findings result from the small initial inhomogeneity of the liquid introduction. Despite considerable efforts, a small inhomogeneity of M φ ≈ 1 mm has always been present. Presumably, this is a feature of the water addition to a random foam structure and cannot be avoided. Possibly, a similar inhomogeneity was also present in former studies. But since former studies did not measure the horizontal liquid distribution, this aspect may have gone unnoticed. Despite the different geometry and critical liquid fractions, the dependency of our stability criterion on the bubble size φ crit ∝ R −4/3 b shows the same trend as earlier findings This work takes only monodisperse foam into account. In industrial applications, foam is usually rather polydisperse. Polydispersity modifies the ingredients of the instability mechanism, namely the foam permeability, the anisotropic drainage term, the shear modulus and the capillary pressure term in the drainage equation. But none of these ingredients would be cancelled or negated in the case of polydispersity. Thus, we assume that our mechanism does apply to polydisperse foam as well, but the critical liquid fraction might change considerably.
For the stability analysis, α ∝ φ 3/2 l is assumed. This is only valid for some surfactants. However, as demonstrated in appendix B, the stability criterion is identical for α ∝ φ 2 l . Since permeability acts on all terms on the right-hand side of the drainage equation, its effect is neutralized in steady state.
The stability analysis shows that any mode k becomes unstable at a certain critical liquid fraction, which scales with φ crit ∝ k 4/3 (see (5.15)). But, according to the linear stability analysis, the smallest modes show the highest growth rate in space and presumably suppress higher modes in many cases. Vera et al. (2000) have observed convective structures with higher modes in systems with bubble diameters below 100 μm. In the present study, the case with the liquid being added by the two outermost porous cylinders (oxxo) corresponds to a mode k = 2π/B and shows superior stability (see figure 5a). However, the resulting convection roll then was not a double roll but a single roll filling the full cross-section (see figure 1a), corresponding to k = π/B.
While the one-dimensional drainage equation without anisotropic drainage does not show any signs of instabilities (Verbist, Weaire & Kraynik 1996), we found an instability in two-dimensional steady drainage which presumably also occurs in three-dimensional steady drainage. This has strong implications for any experiment with columns of liquid foam under forced drainage. In such experiments, one should either measure the horizontal liquid fraction distribution or observe the elastic deformation of the foam.
The proposed mechanism of unstable drainage relies on the anisotropic deformation of the Plateau border network. In the case of very small bubbles or even solid grains and particles, (5.3) presumably loses its validity. This should be further investigated by drainage experiments under prescribed shear, and by simulations. It is possible that neutron imaging  or a combination of the Surface Evolver (Brakke 1992) with drainage equation could yield valuable insights into that phenomenon.