Evidence for a universal saturation profile for radial viscous fingers

Complex fingering patterns develop when a low viscosity fluid is injected from a point source into the narrow space between two parallel plates initially saturated with a more viscous, immiscible fluid. We combine historical and new experiments with (a) a constant injection rate; (b) a constant source pressure; and (c) a linearly increasing injection rate, together with numerical simulations based on a model of diffusion limited aggregation (DLA), to show that for viscosity ratios in the range 300–10,000, (i) the finger pattern has a fractal dimension of approximately 1.7 and (ii) the azimuthally-averaged fraction of the area occupied by the fingers, S(r,t), is organised into three regions: an inner region of fixed radius, r < rb, which is fully saturated with injection fluid, S = 1; a frozen finger region, rb < r < rf (t), in which the saturation is independent of time, S(r) = (r/rb)−0.3; and an outer growing finger region, rf(t) < r < 1.44 rf(t), in which the saturation decreases linearly to zero from the value (rf/rb)−0.3 at rf(t). For a given injected volume per unit thickness of the cell, V ≫ πrb2, we find rf = 0.4rb (V/rb2)1/1.7. This apparent universality of the saturation profile of non-linear fingers in terms of the inner region radius, rb, and the injected volume V, demonstrates extraordinary order in such a complex and fractal instability. Furthermore, control strategies designed to suppress viscous fingering through variations in the injection rate, based on linear stability theory, are less effective once the instability becomes fully nonlinear.

If air is injected from a point source into a second immiscible fluid confined within a narrow planar channel, then, with a constant injection rate, the interface eventually becomes unstable and a series of viscous fingers develop 1 . Linear stability theory suggests that the fingering process may be suppressed if the injection rate is varied with time, with the optimal profile involving a linearly increasing injection rate 2 . However, it is not known whether there is any benefit to such an approach once the finger pattern has developed into a fully nonlinear regime. The lack of a simple model for the non-linear fingering process makes it challenging to address such questions. Classical experiments 3 have shown that the nonlinear finger pattern has a fractal dimension of 1.70 ± 0.01, and simulations with diffusion limited aggregation (DLA) models 4 have replicated many of the features of viscous fingering in cases with relatively large flow rates. Guided by this historical work, we now explore the impact of changes in the injection rate with time on the structure of the finger pattern and, in particular, the variation of the saturation with radius.

Results
We have analysed the output from (a) a series of new DLA calculations, (b) many of the classical radial viscous fingering experiments 3,5-8 and (c) a series of new experiments, some of which have injection rate which varies with time: each experiment analysed is listed in the Table in the supplementary information. In Fig. 1a, we illustrate a typical nonlinear viscous fingering pattern: this was obtained from an experiment we performed with constant-flow injection. In Fig. 1b we show the pattern from an experiment with constant-pressure injection 3 , which involves much more tip-splitting than in Fig. 1a. Measurement of the area of the fingers as a function of time from data supplied by the authors 3 shows that in this case the injection rate evolves approximately as = .
Q Q t 0 0 6 . In Fig. 1c we illustrate the pattern produced by a DLA calculation using an off-lattice algorithm 9 , which replicates many features of the nonlinear fingering process.
In order to characterize the fingering pattern, for each experiment we have computed the area covered by fingers, 2π r S(r,t) dr, in the annular ring bounded by the circles with radii r and r + dr, with dr ≪ r. This reveals the variation of saturation S(r,t) with radius at a series of times t. Figure 2 shows the saturation profile for three experiments, corresponding to: (a) constant flow injection, (b) constant pressure injection and (c) approximately linearly increasing injection as well as (d) DLA calculations. In each case, the data are accompanied by a model www.nature.com/scientificreports www.nature.com/scientificreports/ saturation profile (black broken lines), which we now develop and explain. We note that since the oil preferentially wets the glass in favor of the air, there is a thin film of liquid which remains on the walls of the cell as the fingers advance, but this film is much thinner than the gap width 10,11 . The thickness of the film scales with b(μU/σ) 2/311,12 , where b is the cell thickness, μ and U the fluid viscosity and speed, and σ the surface tension, and this has value of order 0.01b, suggesting that the residual liquid film is very thin relative to the gap width. This is consistent with estimates of the injected volume of air based on the product of the gap width and the finger area, which agree within measurement error with the actual volume injected.
As may be seen in Fig. 2 (panels a-f), we find that three distinct regions develop: (i) a central fully flooded region, r < r b , in which S = 1, where r b depends on the injection rate but is approximately independent of time, (ii) an inner frozen finger zone, r b < r < r f (t), where the saturation is independent of time and gradually decays with radius, and (iii) an active outer finger zone, r f (t) < r < r t (t), in which the fingers advance in time but also thicken and undergo tip-splitting events. Historical data have shown that the viscous fingering pattern is fractal with dimension D = 1.7. If the saturation in the frozen finger region decays with radius as then for r f  r b , integration shows that this saturation profile is consistent with the region covered by the frozen finger pattern being fractal (see eqn 2 below).
Observations of the saturation in the finger patterns shown in Fig. 2 suggest that in the evolving outer finger zone, r f (t) < r < r t (t), the saturation decreases to zero approximately linearly with radius. This suggests a simple model for the saturation in this region given by here, the outer radius of the finger zone, r t (t), corresponds to the minimum radius at which S(r, t) = 0, and we define Λ according to the relation r t = Λr f . Integration of the model saturation profile across the frozen finger zone and the advancing front region leads to the approximate relation which is valid provided that V(t), the invaded area per unit thickness of the cell at time t, satisfies V(t) ≫ πr b 2 . In order that this expression represents a fractal shape for all values of r f (t), with a constant value of D, we expect that Λ is also a constant.
By comparing the saturation profile for the DLA calculations with eqn (1) at a series of times, we find that the best fit values for the model with fractal dimension D = 1.7 3,4,9 are Λ = 1.43 ± 0.04 and r b = 0.14 ± 0.01, where lengths have been scaled with the radius of the particle. To find the best fit, we minimise the quantity: www.nature.com/scientificreports www.nature.com/scientificreports/ To illustrate the sensitivity of the model to the values of Λ and r b we show the contour on which E = 0.1, as a function of Λ and r b with D = 1.7 (Fig. 3a). We have carried out a similar analysis for each of the viscous fingering experiments to determine the best fit values of Λ and r b, and in Fig. 3b illustrate the contour E = 0.1 for a typical experiment. Repeating this analysis across all our experiments, we find that Λ = 1.44 ± 0.01 for D = 1.7.
We have also estimated the value of r b for each experiment as determined from (i) the point at which the saturation first decreases below unity (crosses, Fig. 3c), and (ii) the best fit of the data to the model saturation profile (circles; Fig. 3c). The dependence of r b on the source parameters is more subtle, but for each experiment a www.nature.com/scientificreports www.nature.com/scientificreports/ convenient length-scale with which to scale r b is the radius r 0 at which the interface first becomes unstable according to linear stability analysis 1 , as given by  Figure 3c shows the dependence of this dimensionless radius, r b /r 0 , on the non-dimensional injection parameter, C*, as given by where the areal injection rate Q(t) = Q 0 t α , μ is the viscosity, b the cell thickness, and σ the surface tension. It is seen that r b /r 0 increases with C* corresponding to there being a larger fully saturated region in a more intense flow. It is interesting to note that r b  r 0 , which indicates that the fully saturated zone extends some distance beyond the radius at which the linear instability initially sets in. This may be expected since at the onset of instability, the interface is still migrating outwards, and as the radius increases, the wavenumber of the most unstable mode also changes 1,2 .
As a further test of the model, Fig. 4 illustrates the variation with time of the outer radii of the frozen finger zone, r f , and the growing finger zone, r t , scaled with the best-fit values of r b for each experiment. Data are given for a large number of experiments from both the literature and our new experiments. We have also evaluated r f /r b (eqn 2) using the parameters obtained from the DLA calculations D = 1.7, Λ = 1.43 and r b = 0.14 (black broken lines), where V is defined in terms of the number of particles added as V = nπ. Figure 4 illustrates that, for all experiments analysed, eqn (2) provides a good approximation for the variation of the radii r f and r t as a function of V/r b 2 , irrespective of whether the experimental injection rate varies in time or is constant, and also for the DLA simulations. This is consistent with Fig. 2 where the model (broken black lines) appears to give good agreement with the evolving saturation profiles for the various injection conditions.

Discussion
The universality of the model for the variation of the saturation profile with radius, and the recognition that there is a transition from the frozen non-linear finger region to the growing non-linear finger region in the leading part of the flow provide new insight about and a succinct description of the evolution of such nonlinear viscous fingers. In contrast to predictions of linear stability theory 1,2,13 , in which varying the injection rate with time but injecting a fixed volume leads to differences in the amplitude of the most unstable mode 2 , the apparent www.nature.com/scientificreports www.nature.com/scientificreports/ universality of the non-linear finger pattern, as characterised by the saturation profile shown herein, does not seem to depend on the injection rate history for a given total volume injected.

Methods
In our experiments, two World Precision Instruments AL1000-220 pumps were used to inject air into a Hele-Shaw cell of radius 30 cm, with plate thickness 1 cm, and gap width of either 1 mm or 1.5 mm set by positioning micrometers. The cell was initially filled with rapeseed oil (μ = 50 mPa s; σ = 31 mN/m 2 ). Another experiment, marked C4 in Fig. 4, was provided by I. Bischofberger et al. and used the same setup 10 but is previously unpublished. It involved the displacement of silicone oil (μ = 298 mPa s; σ = 27 mN/m 2 ) by water. In the DLA simulations, we used code supplied by the co-authors of publication 9 and used the mean of 7 saturation profiles for each aggregate size.