Viscous fingering in a radial elastic-walled Hele-Shaw cell

We study the viscous-fingering instability in a radial Hele-Shaw cell in which the top boundary has been replaced by a thin elastic sheet. The introduction of wall elasticity delays the onset of the fingering instability to much larger values of the injection flow rate. Furthermore, when the instability develops, the fingers that form on the expanding air–liquid interface are short and stubby, in contrast with the highly branched patterns observed in rigid-walled cells (Pihler-Puzović et al., Phys. Rev. Lett., vol. 108, 2012, 074502). We report the outcome of a comprehensive experimental study of this problem and compare the experimental observations to the predictions from a theoretical model that is based on the solution of the Reynolds lubrication equations, coupled to the Föppl–von-Kármán equations which describe the deformation of the elastic sheet. We perform a linear stability analysis to study the evolution of small-amplitude non-axisymmetric perturbations to the time-evolving base flow. We then derive a simplified model by exploiting the observations (i) that the non-axisymmetric perturbations to the sheet are very small and (ii) that perturbations to the flow occur predominantly in a small wedge-shaped region ahead of the air–liquid interface. This allows us to identify the various physical mechanisms by which viscous fingering is weakened (or even suppressed) by the presence of wall elasticity. We show that the theoretical predictions for the growth rate of small-amplitude perturbations are in good agreement with experimental observations for injection flow rates that are slightly larger than the critical flow rate required for the onset of the instability. We also characterize the large-amplitude fingering patterns that develop at larger injection flow rates. We show that the wavenumber of these patterns is still well predicted by the linear stability analysis, and that the length of the fingers is set by the local geometry of the compliant cell.


Introduction
When a viscous fluid is displaced by a gas bubble in the narrow gap between the rigid walls of a Hele-Shaw cell, small perturbations to an initially circular interface readily develop into fingers whose tips subsequently become unstable themselves and split. Repeated tip splitting combined with the near arrest of the bases of the fingers, i.e. the portions of interface left behind the developing fingers, ultimately creates highly branched patterns (Paterson 1981;Chen 1989;Lajeunesse & Couder 2000). This fingering instability is an archetype for pattern-forming phenomena in a large variety of systems (Saffman & Taylor 1958;Mullins & Sekerka 1964;Ben Jacob et al. 1992;Clanet & Searby 1998;Hull 1999;Couder 2000).
In this paper, we study the development of viscous fingering in a radial elasticwalled Hele-Shaw cell whose top boundary has been replaced by an elastic sheet as shown in figure 1(a,b). The axisymmetric expansion of the interface in this set-up has been characterized both experimentally and theoretically by Lister, Peng & Neufeld (2013), Pihler-Puzović et al. (2013 and Peng et al. (2015). One of the features of this flow is the presence of an exponentially damped oscillatory perturbation to the thin liquid layer ahead of the interface, which is common to a broad range of systems, see, e.g. Gaver III et al. (1996), Heil (2000), Salez et al. (2012) and Dixit & Homsy (2013). In elastic-walled cells the injected gas works against elastic, viscous and capillary forces to form a blister (Juel, Pihler-Puzović & Heil 2018), thereby creating a converging gap into which the interface propagates (see figure 1b). Pihler-Puzović et al. (2012) found that, in this system, the fingering instability is weakened or even suppressed in the sense that the instability only arises at much larger injection rates than in an equivalent rigid cell. When the interface becomes unstable in an elastic-walled cell, the entire interface (i.e. the tips and the bases of the fingers) continues to propagate outwards. As a result, a large number of relatively short and stubby fingers develops on the interface as shown in figure 1(c,d). Similar patterns have been observed in other examples of two-phase flows, where a gas-liquid interface exists in a converging gap between solid boundaries, e.g. in roll coating (Ruschak 1985;Rabaud, Couder & Michalland 1991;Couder 2000;Weinstein & Ruschak 2004), tape peeling (McEwan & Taylor 1966) and models of airway reopening (Ducloué et al. 2017a,b).
The suppression of fingering in radial elastic-walled cells was explored numerically by Pihler-Puzović et al. (2013), who showed that the key mechanisms are the slowing of the interface and the formation of a convergent gap ahead of the interface due to the inflation of the sheet. These features have been exploited to control viscous fingering in rigid Hele-Shaw cells, by varying e.g. the injection rate (Li et al. 2009;Dias & Miranda 2010;Dias et al. 2012), the separation of the rigid bounding plates (Zheng, Kim & Stone 2015) and the cell geometry (Al-Housseiny, Tsai & Stone 2012).
In radial elastic-walled cells, neither the unperturbed circular interface nor the fingering instability itself grow at a constant rate when the injection volume flux is fixed. Therefore, the fingering is transient in nature, similarly to other instabilities with time-evolving base states (Trevelyan, Almarcha & De Wit 2011;Haudin et al. 2014). Viscous fingers grow initially but they do not saturate to a constant length, and instead ultimately decay again.
In contrast, viscous-fingering instabilities in rectilinear geometries are not transient. One of the earliest studies that reports the development of saturated finite-amplitude fingers in such geometries and involves fluid-structure interaction is McEwan & Taylor's (1966) study of tape peeling. The paper shows photographs of saturated, short Viscous fingering in a radial elastic-walled Hele-Shaw cell  and stubby fingers on the air-liquid interface but the authors did not systematically explore these patterns. Peng & Lister (2018) studied theoretically the fingering instability in an idealized rectangular elastic-walled Hele-Shaw cell, which could support a steadily propagating infinitely wide flat front. They used linear stability analysis to show that the main stabilizing effect of the elastic sheet is due to the fact that the deforming sheet creates a converging geometry with a taper angle α i , shown schematically in figure 1(b), into which the interface advances. Approximation of the geometry as a triangular wedge allowed the growth rate of the perturbations to be expressed analytically in terms of α i and the base-state capillary number Ca = µU/γ , where γ is the surface tension at the interface, U is its steady speed and µ is the dynamic viscosity of the viscous fluid. Ducloué et al. (2017a,b) studied the two-phase displacement flow in a finite-width elastic-walled rectangular Hele-Shaw channel; see also McCue (2018). In a rigid channel this flow always results in the formation of a single Saffman-Taylor finger which propagates steadily and is linearly stable (Saffman & Taylor 1958;McLean & Saffman 1981;Combescot et al. 1986). However, if one of the cell boundaries is replaced by an elastic sheet then the propagating interface is prone to instabilities which emanate from the tip of the finger. At high capillary numbers Ca, the elastic-walled channel can support a modified Saffman-Taylor finger with a locally flat front, which is in turn susceptible to Peng & Lister's (2018) fingering instability. The emerging finite-amplitude fingers are similar to the short and stubby fingers shown in figure 1(c,d), but their average width and length do not vary in time. They resemble the fingers observed in the printer's instability that develop on the meniscus between two cylinders co-rotating with the same speed (Couder 2000). Ducloué et al. (2017a) characterized the fingers as a function of the base-state capillary number Ca and the peeling slope α i . It was not possible to compare the growth rate of the instability to the predictions from Peng & Lister's (2018) linear stability analysis because the advancing front of the modified Saffman-Taylor finger was flat only at sufficiently high Ca. This meant that the onset of fingering could not be captured because the instability was strongly nonlinear for all values of Ca investigated.
In this paper, we extend the axisymmetric models considered by Peng et al. (2015) and Pihler-Puzović et al. (2015) to include non-axisymmetric perturbations and perform a linear stability analysis similar to that of Peng & Lister (2018) to obtain predictions for the growth rates of these perturbations in an elastic-walled radial Hele-Shaw cell. We compare the predictions to experimental results and obtain good quantitative agreement with the linear stability analysis for injection rates slightly above the critical flow rates required for the onset of the instability. We also characterize the large-amplitude fingering patterns at higher injection rates. We find that while nonlinear interactions do not affect the pattern selection (i.e. the number of fingers that form on the interface), they do control the maximum length to which the fingers grow before they start to decay again.
The paper is organized as follows. We start by describing the experimental method in § 2 and continue with an outline of the theoretical model and linear stability analysis in § 3. The results are presented in § 4, where we compare the predictions from the linear stability analysis with experimental measurements ( § 4.1), and characterize the large-amplitude fingering patterns in the nonlinear regime ( § 4.2). Conclusions are given in § 5.

Experimental set-up
A schematic of the experimental set-up is shown in figure 1(a) (see also Pihler-Puzović et al. (2012) and Pihler-Puzović et al. (2015)). The experiments were performed in an elastic-walled Hele-Shaw cell with a rigid glass bottom and an elastic top boundary separated by a latex spacer with a circular cutout of radius R cell = 180 mm. The top boundary was a thin square elastic sheet, made either of latex (Young's modulus E = 2.1 MPa and Poisson's ratio ν = 0.5) or polypropylene (E = 3.7 GPa and ν = 0.44) with side length 2R lid = 400 mm and uniform thickness d. The cell was initially filled with viscous liquid (silicone oil of viscosity µ = 0.962 Pa s, density ρ = 961 kg m −3 and surface tension γ = 0.021 N m −1 for the latex cell, or a glycerol-water mixture with µ = 0.305 Pa s, ρ = 1235 kg m −3 and γ = 0.065 N m −1 for the polypropylene cell) to a uniform height h 0 equal to the thickness of the spacer and radius R fluid ≈ 150 mm. Prior to the experiment, a small amount of gas was injected through a central hole in the base to form a bubble with radius R init ≈ 5 mm. Nitrogen was then injected at a fixed flow rate Q in the range 100-1300 ml min −1 , which caused the bubble to expand radially as well as inflate the sheet vertically. The evolution of the expanding, initially approximately circular gas-liquid interface and the development of the fingering instabilities were captured by a camera with a top-down view of the experiment. Retraction of the interface due to dewetting in the case of the glycerol/water solution on polypropylene occurred on time scales of the order of minutes, which were much longer than the experimental time scales ( 15 s). Hence, its influence on the fingering instability was negligible.
The values of d, h 0 and Q used in the experiments are shown in table 1. In the experiments with latex sheets, we either varied d while h 0 was held constant at 0.56 mm (cells named D33-D97, respectively) or varied h 0 while d was held constant at 0.56 mm (cells named H46-H79, respectively). The cell with h 0 = d = 0.56 mm is listed as both D56 and H56 in table 1. In the experiments with polypropylene sheets, we investigated only one cell with h 0 = 0.56 mm and d = 0.03 mm (named PP).
In order to quantify the variability inherent in the system, we performed at least six experiments for each value of Q in each latex cell and three experiments for each value of Q in the cell with the polypropylene sheet. In total, more than 650 experiments were performed, during which we observed the development of the instability. Results for lower values of Q, for which the interface did not develop fingers, are not shown in table 1 but were reported in Pihler-Puzović et al. (2015). More details of the experimental protocol can be found in appendix A.
The typical time evolution of the patterns observed in the latex and polypropylene cells is shown in figure 1(b,c), using cells D56/H56 and PP with a volume flux Q = 500 ml min −1 . At the start of the experiment, the perturbation amplitude is not uniform across modes and it is not a priori constant between experiments. Therefore, each pattern includes a wide range of finger widths -even at the very early stages of the system's evolution. This is particularly visible in the innermost interface shape in figure 1(c), indicating that the system is very sensitive to initial perturbations. Tip splitting, most visible in the top right corner of the pattern in figure 1(c), occurs as the mean radius (circumference) of the pattern increases.

Image processing and extracted quantities
All images were processed in MATLAB 2013a by first converting each true colour picture to a binary image using an input value for the intensity threshold, resulting in a white image with black pixels representing the interface. Occasionally this procedure resulted in discontinuities in the interface of up to two pixels, particularly for the thicker (less translucent) latex sheets. Such discontinuities were filled in using a morphological closing algorithm (bwmorph in MATLAB). An active contour algorithm was then applied to extract the interface as a sequence of points x i , with 1 i M, going around the gas bubble. (The total number of points, M, increased as the interface expanded.) In regions where the interface was thicker than a pixel, the active contour algorithm detected a meandering path, which we smoothed using a Savitzky-Golay filter with a 31-point window and a third-order polynomial. A typical result of the image processing is shown in figure 2(a), where the extracted interface is superimposed onto a top view of the instantaneous interfacial pattern.
Once the interface of an instantaneous pattern had been identified, the number N(t) of fingers was determined. This required a criterion to distinguish the local minima of the interfacial radius associated with the base of a fully formed finger from those arising from a slight indentation at the tip of a wider finger that was only just beginning to split, and hence could be ignored (see e.g. figures 1-3). We used the criterion that a minimum x j is counted if it lay closer to the centre than the line segment connecting the points x j+λM and x j−λM , which are positioned λM points along the circumference to either side of x j , as shown in figure 2(b). We used λ = 0.015 in the experiments with latex sheets and λ = 0.0315 for the experiments with polypropylene sheets. Hence, the point x i in figure 2(b) was counted as the base of a finger, whereas the point x j was not. We confirmed that this method worked well over the whole range of experiments performed with the same elastic sheet. Alternative criteria based on Fourier decomposition of the patterns or local curvatures of the interface were found to be less robust.
The local minima of the interfacial radius that were deemed to indicate the presence of a finger were then connected with a spline to yield an inner envelope of the pattern (the dot-dashed curve in figure 2a). A corresponding outer envelope was constructed by fitting a spline through every local maximum; here we included multiple maxima on the same finger.
From the outer and inner envelopes of the pattern, the corresponding radial coordinates, R i outer and R i inner (for 1 i M) were found at the M sampling points, and the average radius of the interface was calculated as . (2.1) Throughout this paper, we will plot the evolution of the instability as a function of R(t) as a proxy for time t. This is because the interface undergoes rapid initial growth followed by a much slower evolution, so that a nonlinear scale would be required if t was used as the independent variable. Moreover, data could only be collected while the interface was inside the field of view of the camera, so each data series ends at approximately the same value ofR but at different values of t. The amplitude of the perturbations to the interface is quantified using the average and standard deviation, (2.2a,b) We will refer toR as the 'perturbation amplitude' or the 'finger length' (although the average distance from the tip of a finger to its base is 2R). The standard deviation δR is often a significant fraction ofR (see e.g. figures 1-3 and 11).
Finally, the time-series data forR,R and δR were all smoothed using the same Savitzky-Golay filter used to smooth the raw interface. A typical result is shown in figure 2(c). The error bars show one standard deviation δR above and below the  average valueR. Similarly, experimental results compiled in later figures are shown as mean points with error bars indicating one standard deviation to either side. As shown in figure 2(c), the perturbation amplitudeR initially increases rapidly before reaching a maximum and subsequently decreases upon further expansion of the interface. This behaviour is observed in all experiments. We defineR peak to be the maximal value of the perturbation amplitudeR reached in each experiment (i.e. its value before the instability starts decaying) andR peak and N peak to be the values of the average radiusR and the number N of fingers at this point.
The variability in finger length observed within a single pattern and shown with error bars in figure 2(c) is comparable to the variability between different experimental runs performed under the same control parameters. This can be seen in figure 3(a), which compares results forR in three identical experiments, both in the D56/H56 latex cell and in the polypropylene cell. BothR peak andR =R peak vary significantly between runs as shown in figure 3(c). This is in contrast with the time evolution of the average radiusR(t) between different experimental runs which shows far less variability (see Pihler-Puzović et al. (2015)).
The fixed field of view of the camera used to record the interface in experiments with latex sheets posed a problem for determiningR peak because any maximum inR appearing afterR ≈ 70 mm could not be detected. In some cases, such as run 1 in figure 3(a),R reached a clear maximum and started decaying well beforeR = 70 mm. In other cases, such as run 3, the values ofR appear to plateau aroundR peak , but a well defined global maximum was not seen implying that decay occurred for values ofR > 70 mm. In these cases, the plateau value provided a good estimate ofR peak , but the value ofR at which the maximum occurred could not be detected. In the experiments with polypropylene sheets,R always reached a maximum value within the field of view.
We describe the motion of the viscous fluid occupying the narrow gap underneath the deformed elastic sheet using the lubrication approximation, which yields the equations for the vertically averaged velocity u(r, θ , t) and the fluid pressure as where R(θ , t) is the bubble radius.
In the region r < R(θ, t), the gas bubble has a uniform pressure p g (t). We note that since the liquid wets the walls of the cell, thin films of liquid are left behind the advancing bubble tip. We assume that the films have no dynamical effect in this region, so that the pressure acting on the elastic sheet is given by (3.4) At r = R(θ , t), we impose the kinematic and dynamic conditions where U n is the speed of the gas-liquid interface in the direction of the unit normal, n, to that interface. The functions f 1 and f 2 model the effect of the films that are deposited on the upper and lower walls. Comparisons against full free-surface Navier-Stokes simulations by Peng et al. (2015), Pihler-Puzović et al. (2015) showed that, for the capillary numbers encountered in our experiments, the presence of these films has a significant effect on the axisymmetric base flow. The choices which were obtained by a fit to Reinelt & Saffman's (1985) computational data for the width of the finger penetrating a viscous fluid between parallel plates and the pressure drop across the tip of that finger as a function of the capillary number, produced results that were in excellent agreement with solutions of the full Navier-Stokes equations. We neglect the small corrections to the expressions (3.6) due to the walls of the cell being non-parallel, and refer to papers by  and Jensen et al. (2002) which study the asymptotic structure of the flow in similar geometries with non-parallel walls. The effect of the films deposited on the upper and lower walls can be neglected by taking f 1 = 0 and f 2 = 1. Alternative models which also account for the effect of the films in rigid cells can be found in Martyushev & Birzina (2008), Anjos & Miranda (2013) and Jackson et al. (2015).
In the experiments, the sheet rests freely on the viscous layer and there is no significant deformation far ahead of the gas-liquid interface. We mimic this behaviour in the outer boundary conditions of the elastic sheet, but use simplified conditions for the fluid, by imposing with R out = 200 mm. We checked that changing R out from 150 mm (edge of the fluid layer in experiments) to 300 mm (corresponding to the diagonal of the elastic sheet) changed the data reported in § 4 by less than 2 %. Initially, the cell is unperturbed and, as in the experiments, contains a small bubble of radius R init = 5 mm, so we imposed the initial conditions: Finally, mass conservation requires that 3.2. Linear stability analysis We investigate the early states of the fingering instability using linear stability analysis. The base state is axisymmetric, so we employ a Fourier decomposition of the linear perturbations in the azimuthal direction and consider each mode separately. Each quantity f is thus decomposed into a sum of an axisymmetric but time-dependent base solutionf and a perturbation with small amplitude 1 and azimuthal wavenumber n. The resulting equations are given in appendix B. At O(1) we recover the governing equations for the axisymmetric base state (B 1), which were analysed by Peng et al. (2015), and at O( ), we obtain the governing equations for the perturbations (B 2), which we shall analyse here. Henceforth, we use primes and dots to denote differentiation with respect to r and t, respectively.
The evolution of the perturbation was followed from the initial conditionŝ which corresponds to an initial perturbation to the interface position only. We solve the governing equations, (B 1) for the base flow and (B 2) for each linear perturbation with wavenumber 1 n 50, simultaneously using a fully implicit backward-Euler finite-difference scheme with adaptive spatial grid and adaptive time stepping (for details see appendix B). We define σ for each azimuthal wavenumber n as in (2.3). We also compare the instantaneous growth rates to the growth rates obtained using the 'frozen-time' approximation (the eigenvalue analysis of (B 2), conducted at each time step using the numerical solutions to (B 1)).
In figure 4(a) we show the time evolution of a representative base-state height profileh and interfacial position. The injected gas deflects the sheet upwards and causes the gas-liquid interface to spread outwards. Ahead of the interface, the displaced liquid accumulates in a wedge whose size increases with time. In figure 4(a) below the axish = 0 the symbols show the interface observed in experiments, with each point corresponding to one of the 12 repeated experimental runs and the horizontal bars indicating the rangeR ±R observed in them. The model can be seen to capture the mean position of the interface in the experiments satisfactorily.
To analyse the characteristic features of the perturbations, we focus on the mode that is predicted to grow to the largest amplitudeR in figure 4, namely n = 25. The height perturbationsĥ are shown in figure 4(b) using solid lines, for the same three times as in figure 4(a). They have multiple undulations in the vicinity of the wedge region, but are otherwise very small. A more detailed inspection of the data shows that the perturbations to the fluid pressure and to the position of the air-liquid interface are much larger than the perturbations to the deflection of the sheet. Specifically, (ĥ/h)/(R/R) never rises above a value of O(10 −6 ) in the time interval shown in this figure. In contrast, the corresponding pressure perturbationsp, illustrated in figure 4(c) using solid lines, are much more significant, with (p/p)/(R/R) being in the range O(10 2 -10 3 ) over the same time interval. We will exploit this observation in § 3.3 below.
In figure 4(b), we also show the results obtained from the 'frozen-time' approximation using dotted curves. We find that within the experimental field of  (d) The evolution of amplitudesR for a range of different wavenumbers, comparing the results from the full analysis (solid curves) to those obtained when using the rigid-wedge approximation (dashed curves). The maximum for each wavenumber is marked by a vertical bar.
view the maximum instantaneous growth rate σ obtained using this approximation is no more than 1 % larger than that obtained from the direct calculations for all wavenumbers between 5 and 50.
In figure 4(d), we show the logarithm of the perturbation amplitudeR as a function of the mean radiusR for wavenumbers between 10 and 40 in increments of 5. A large number of modes can be seen to have comparable (positive) growth rates and the linear analysis predicts all of these modes to grow to a maximum amplitude before decaying. Which one of these modes will actually become dominant in a specific experiment depends on the relative initial amplitudes of the perturbations and on nonlinear effects that start to play a role when the amplitude of the perturbation becomes sufficiently large, an issue we will discuss in more detail below. It is, however, interesting to observe that in figure 4(d) a single mode (n = 25) maintains the largest instantaneous growth rate throughout the system's evolution and therefore also grows to the largest amplitude. This behaviour was found to be typical for all the cases considered in this study (see table 1). 3.3. The rigid-wedge approximation and the physical mechanism of fingering suppression The results in figure 4 show that the most rapidly growing modes have relatively large wavenumber, n 10, and hence comparatively short azimuthal (lateral) length scales. As discussed in § 3.2, we also find that the pressure perturbationsp only generate a small perturbationĥ to the sheet (i.e. gap height), suggesting thatĥ can be neglected. Moreover, since the perturbations were found to be mostly confined to the wedge region (see figure 4b,c), we approximate the perturbation flow ahead of the interface as a flow inside a wedge-shaped domain with straight converging walls. These are precisely the 'rigid-lid' and 'wedge' approximations, first introduced by Peng & Lister (2018) for the idealized rectangular elastic-walled Hele-Shaw cell. Using these approximations, we will now derive an analytical expression for the instantaneous growth rate σ of the different modes in terms of the (numerically calculated) base-state quantities at the interface. This will then allow us to identify the physical mechanisms by which the presence of the elastic sheet weakens (or even suppresses) the fingering instability.
In order to simplify the notation, we denote the base-state height, slope and (compressive) longitudinal velocity gradient at the interface (subscript 'i') by (3.12a−c) By assumingĥ = 0 (the 'rigid-lid' approximation), substitutingṘ = σR and using (B 1d), the interfacial conditions (B 2e), (B 2f ) become ( whereû andp are evaluated at r =R + . The dynamic condition (3.13b) describes how the pressure perturbation is generated. The first term on the right-hand side describes how the base viscous pressure drop associated with the mean advance of the interface results in a relative increase in pressure at the tip of the fingers. This is the driving force of the instability, as the increased pressure drives a perturbation flow that causes growth of the fingers. The term includes a stabilizing correction due to the wetting films left behind the interface, as these films reduce the amount of fluid that needs to be displaced to achieve a given interfacial speedṘ and hence reduce the base viscous pressure drop. There are also two stabilizing surface-tension terms in (3.13b). The first term is due to the perturbations to the in-plane curvature; this is the main restoring effect in the classical case of viscous fingering in a rigid-walled Hele-Shaw cell. The second term is due to surface tension acting on the vertical curvature of the interface as observed in a radial cross-section of the cell. Due to the converging geometry of the cell, the tip of the finger sits in a narrower gap and hence has a larger vertical curvature. The result is a lowering of the pressure at the tip, which has a restoring effect. This is the 'taper' mechanism described by Al-Housseiny et al. (2012).
In the kinematic condition (3.13a), the first term on the right-hand side describes how perturbations to the interface grow or decay due to advection of the interface by the perturbation velocity, which is driven by the pressure perturbation from (3.13b). The second term describes how spatial variations in the base flow affect the advection of the finger tip. The simulations show that typically c i = −ū > 0, so this is a stabilizing effect and we call it 'kinematic compression'. (In Juel et al. (2018), the term 'geometry' is used instead.) The correction factors in the parentheses on the left-hand side account for the effects of the thin wetting films left behind the tip.
In order to solve for σ in (3.13), we need to calculate the ratioû/p = −h 2 ip /(12µp) by solving the perturbation lubrication equation (B 2c) (while neglectingĥ). We introduce a local variable x =R + h i /α i − r, and approximate the base-state height profile as a wedgeh = α i x with slope α i and height h i at the interface r =R. Assuming that the wedge is short compared with the radius, x R , allows us to simplify the lubrication equation further by neglecting the terms arising from the radial geometry to obtain where the boundary condition is analogous to (B 2i). The solution isp ∝ I 1 (nx/R)/x, where I 1 is the order-1 modified Bessel function of the first kind. Following Peng & Lister (2018), we rescale the ratioû/p (evaluated at r =R + or equivalently The scaling is chosen so that Y = 1 for an infinite cell with rigid and parallel walls, while the restricted flow in a converging wedge gives a smaller value ofû/p and hence Y < 1. Using the admittance, we can eliminate the unknownsp/R andû/R from (3.13) to obtain (3.16) The evolution of the perturbation amplitude is then given bŷ The denominator in (3.16) is always positive, so the stability of the perturbation is determined by the sign of the numerator. The perturbation decays, i.e. σ < 0, provided that These four terms represent the stabilizing effects of films deposited behind the interface, taper, horizontal surface tension and kinematic compression, respectively, and they give a quantitative measure of the stabilization due to each physical mechanism.
Using the base-state equations (B 1b) and (B 1d), we can further decompose the longitudinal compression into These three terms represent three mechanisms that stretch the base flow near the interface in the azimuthal and vertical directions, and cause the base flow to compress longitudinally through mass conservation. In order to understand these three terms, we consider the local frame of reference moving with the interface at velocityṘ so that the interface is stationary. The first term in (3.19) is the same as in the unidirectional steady peeling case (Peng & Lister 2018), and describes how the flow towards the interface, with velocitȳ u −Ṙ = −f 1Ṙ , experiences a vertically widening cell at the relative rate α i /h i per unit length in the longitudinal direction. The second term is due to the radial geometry of the cell, and is the relative expansion rateṘ/R of the circumference multiplied by the correction factor (1 − f 1 ) to account for the deposited films. Finally, the third term is the relative rate of change of the interfacial gap height (in the frame of reference moving with the interface). Typically, the height is (slowly) increasing with time, so this represents a stretching in the vertical direction. All three terms in (3.19) are positive and therefore contribute to the suppression of the instability.
For n 1, the expression (3.16) for σ reduces to the expression found by Peng & Lister (2018) if we identify the lateral wavenumber as k = n/R. The only difference is in the compression term involving the velocity gradient c i = −ū . The kinematic compression was also found to be stabilizing in the cell studied by Peng & Lister (2018), but in their two-dimensional time-independent system the peeling process did not involve the additional mechanisms provided by the second and third terms on right-hand side of (3.19). Hence, we conclude that the present radial time-evolving system is more stable than the one studied by Peng & Lister (2018).
Our numerical results (e.g. figure 4) show that, in the regime under consideration, the perturbation pressure profiles obtained using the full calculation and the rigidwedge approximation agree well overall, particularly near the interface (see figure 4c). Consequently there is excellent agreement between the models regarding the evolution of the perturbation amplitudesR and growth rate σ (figure 4d). Thus, for simplicity, we focus henceforth on results from the rigid-wedge approximation only.

Results
We discuss the development of fingers as the interface expands by direct comparison between the experiments and the rigid-wedge model in § 4.1. Nonlinear pattern formation in the experiments is then explored in § 4.2.

Small-amplitude perturbations
The curves in figure 5(a) show the theoretically predicted instantaneous growth rate σ as a function of the wavenumber n from (3.16), at three different instants when the average radius isR = 25, 35, 45 mm, for the D56/H56 latex cell with injection rate Q = 400 ml min −1 . Positive growth rates are observed for a range of intermediate wavenumbers, while modes with smaller or larger wavenumber are seen to decay. The experimental data from each of 12 repeated runs with the same parameters are plotted with symbols. AtR = 25 mm, there is a rather large spread in both wavenumber and growth rate, likely due to the sensitivity to initial conditions discussed in § 2.1. However, nearR = 35 mm andR = 45 mm the model satisfactorily captures the development of the fingering pattern, although the wavenumbers observed in the experiment are generally smaller than that of the mode for which the linear stability analysis predicts the maximum instantaneous growth rate. A typical example of the evolution of the growth rate σ with the average radius R is shown in figure 5(b), for the same cell as in figure 5(a). The experimental data are obtained by averaging 12 repeated runs from figure 5(a) and the shading indicates the standard deviation. We find that the instantaneous growth rate of the most rapidly growing mode (over all wavenumbers) calculated numerically provides an upper limit on the growth rates obtained in the experiments. Initially it is very large (see also figure 4d), but at later times, the growth rate decreases and eventually becomes negative (figure 5b). Henceforth, we refer to the maximum (i.e. dominant) growth rate across all wavenumbers in the numerical simulation as simply the growth rate σ .
In figure 6, we show σ 35 , the value of the growth rate measured atR = 35 mm, as a function of Q for all the experimental parameters given in table 1. We compare the experimental data (symbols) with predictions from the linear stability analysis (solid lines) for Q ranging from 50 to 1300 ml min −1 in increments of 50 ml min −1 or less. The values Q = Q 35 , for which σ 35 = 0 in the linear stability analysis, are shown in table 1.
In each case the value of σ 35 increases with Q, because a larger injection rate leads to faster spreading and hence a stronger viscous force driving the instability. The growth rates observed in the experiment are in excellent agreement with the linear stability analysis when Q is only slightly larger than Q 35 (see the regions where symbols overlap with the lines in figure 6). The deviation of the experimental growth  rates from the theoretical predictions for larger values of Q suggests that nonlinear effects become important within the experimental field of view. A closer analysis of the data shown in figure 6 reveals that the value of the growth rate σ 35 predicted by linear stability analysis increases with the sheet thickness d (cells D33-D97) and decreases with the initial cell depth h 0 (cells H46-H79). Increasing d stiffens the sheet so that it deforms relatively less. This allows the interface to propagate faster (i.e. with a larger capillary number Ca = µṘ/γ ) and also reduces the stabilizing taper angle α i . Both of these effects result in a faster growth of the instability, as discussed in § 3.3. On the other hand, we find that an increase in the value of h 0 for a given volume flux results in the interface propagating more slowly, which reduces σ 35 , but also in a decrease in the taper angle α i , which increases σ 35 . Therefore, an increase in h 0 could in principle be stabilizing or destabilizing depending on the relative magnitude of these opposing effects; see also Pihler-Puzović et al. (2012).
We quantify the balance of stabilizing and destabilizing effects in our system using the criterion (3.18). Figure 7 shows the evolution of the four stabilizing terms from (3.18) in the D56/H56 latex cell with Q = Q 35 (see table 1) for the mode n = 13, which is the last one to decay. We observe that the sum of the terms initially lies below 1, indicating that the mode is unstable. Following a short initial decay, the sum then grows steadily and crosses 1 whenR = 35 mm, where the mode stabilizes. At the time of stabilization, we find that all four physical mechanisms have a significant contribution. However, at late times, when Ca 1, the taper mechanism has the largest contribution towards stability as was argued by Al-Housseiny, Christov & Stone (2013) and Peng & Lister (2018), although taper alone would only stabilize the system at R ≈ 70 mm. Similar behaviour is observed for all other cells in table 1. Figure 8 shows the path taken by the system in Ca-α i space, for each cell in table 1 with Q = Q 35 . In each case, the system initially has a large capillary number and small interfacial slope. However, immediately after the volume flux is imposed the slope quickly increases, and then remains approximately constant as the capillary number decreases due to the interface slowing down. AtR = 35 mm, marked with vertical bars in figure 8, the system becomes stable, while the capillary number continues to Viscous fingering in a radial elastic-walled Hele-Shaw cell . (Colour online) The number of fingers N peak observed at the maximum ofR in experiments (points with error bars indicating standard deviation) and the wavenumber N peak of the mode that grows to the largest amplitude before decaying in the linear stability analysis (solid curves), as a function of the injection volume flux Q, for each cell in table 1. The theoretical data were calculated for Q ranging from 50 ml min −1 to 1300 ml min −1 in increments of 50 ml min −1 and the corresponding integer values of N peak were connected with straight lines, which appear as steps in the curves. decrease. The shaded region in figure 8 is the stable region calculated by Peng & Lister (2018) for steady unidirectional-peeling solutions. It underestimates the stability of the present time-dependent axisymmetric system, but only by a small amount, as can be expected from our discussion in § 3.3.

Nonlinear evolution and pattern formation
Next we study the pattern formation beyond the early development of fingering. As discussed above, a key feature of the fingering instability in the radial Hele-Shaw cell is that the instability is always transient: small-amplitude perturbations to the axisymmetrically expanding interface start to decay once its radius has exceeded a certain value. Thus even the linearized analysis of § 3.2 predicts that fingers grow to a maximum length before decaying; see, e.g. figure 4(d). However, the maximum finger length and the time of its occurrence predicted by the linear analysis differs significantly from experimental observations, indicating that the saturation and the ultimate decay of the fingers occurs in a nonlinear regime.
Nevertheless, figure 9 shows that the wavenumber N peak observed atR peak in the experiments is close to the wavenumber of the mode that is predicted to grow to the largest amplitude before decaying by the linear stability analysis. Recall that, as discussed in § 3.2, this mode also tends to have the largest instantaneous growth rate throughout the system's evolution, particularly during the early stages of the instability when the linearized analysis can be expected to be accurate. This suggests that the number of fingers is set during the early development of the instability when the linear stability analysis applies. Comparing the results for different cells with latex sheets in figure 9, we find that the number of fingers N peak decreases with increasing initial cell depth h 0 (H46-H79) and with increasing sheet thickness d (D33-D97) (or equivalently with extensional stiffness Ed), although it is much less sensitive to the latter. Indeed, increasing d by a factor of three from D33 to D97 reduces N peak only by a small amount. A significant change in N peak was only seen when comparing the latex sheets with the polypropylene sheet, whose extensional stiffness is two orders of magnitude larger than that of the latex sheets. However, we note that a different liquid was used in the experiments with polypropylene sheets, which may be an additional contributory factor to the observed variation.
The dependence of the maximum finger lengthR peak on the injection flux Q is shown in figure 10 for latex sheets of different thickness and h 0 = 0.56 mm; the same behaviour was observed for different values of h 0 .R peak initially grows with increasing Q before it reaches an approximately saturated plateau. This apparent saturation correlates with the taper angle created by the inflation of the elastic membrane reaching an approximately constant value (see below); see also Ducloué et al. (2017a). We refer to this range as the 'saturation regime' with respect to the parameter Q rather than to the time t. In figure 10, we fit a straight line through the region of growth inR peak and a horizontal line in the region of approximate saturation inR peak minimizing the overall mean-squared error to find the approximately constant valueR sat ofR peak observed for each cell in the saturation regime. The typical fingering patterns observed in the saturation regime are shown in figure 11. The finger length increases with both the sheet thickness (figure 11, D33-D97) and the initial cell depth (figure 11, H46-H79). We quantify this increase by plotting the variation ofR sat with d and h 0 in figure 12. Since the sheet thickness d affects the system mainly by altering the extensional stiffness Ed of the sheet, we plot the data for both latex (D33-D97) and polypropylene (PP) sheets with h 0 = 0.56 mm on a semi-logarithmic scale in the inset to figure 12(a). The trend of increasing finger length for the latex sheet experiments, for which Ed varies by a factor 3, is consistent with the polypropylene data point, for which Ed is two orders of magnitude larger.
Looking back at figure 4(a,b), we notice that at late times the finger length observed in the experiments is comparable to the lengthL w of the liquid wedge ahead of the interface in the numerical simulation of the axisymmetric base state. (We defineL w as the distance from the interface r =R to the first point ahead of the interface whereh(r, t) crosses h 0 .) Since the maximumR peak is attained in each experiment at a slightly different time, and the saturation regime spans different ranges of Q, we simply calculated the minimum and maximum values ofL w obtained in the numerical simulations across all times corresponding to the range 60 mm R 80 mm, for all flow rates in the range 400 ml min −1 Q 1300 ml min −1 . We found thatL w was comparable to 2R sat , i.e. the distance from base to tip of the fingers. Hence, we plot the results forL w /2 in figure 12, using boxes that span the range between the minimum and maximum values, for comparison withR sat . The correlation between L w /2 andR sat is an indication that the nonlinear saturation process is acting on a length scale comparable to the length of the wedge.
Our results are consistent with those reported by Ducloué et al. (2017a). Their experiments studied the fingering of a flat front that propagated steadily below an elastic sheet, thus yielding truly saturated fingers (of constant lengthR peak ) both in time as well as in the control parameter Ca. Firstly, they found a linear dependence of the finger lengthR peak on the initial layer thickness h 0 , which agrees with our results in figure 12(b). Secondly, their data indicate thatR peak is approximately proportional to 1/tan(α i ), which in their set-up is equivalent toR peak increasing with the wedge length. This is also consistent with our results in figure 12. Similar behaviour is observed in the printer's instability. We conclude that cells with stiffer walls have a smaller peeling angle and hence a longer wedge, which in turn leads to longer saturated fingers. In the experiments of Ducloué et al. (2017a) this implied that the ridges of viscous fluid that separate adjacent fingers terminated (at the upstream end) when they reached a critical height within the tapered channel. In the limit of infinite stiffness (rigid cell with parallel walls), the finger length does not saturate (Saffman & Taylor 1958;Paterson 1981).

Conclusion
We have studied viscous fingering in an elastic-walled Hele-Shaw cell, from the early-time linear growth to the late-time nonlinear saturation and decay. A direct comparison between experiments and a linear stability analysis confirms that the fingering initially develops as a linear perturbation to an axisymmetric base state.
The linear analysis elucidates several physical mechanisms by which the elasticity of the sheet has a stabilizing effect on the viscous-fingering instability. Firstly, since the instability is driven by the injection of gas at a constant flow rate Q, the inflation of the sheet reduces the spreading rateṘ of the interface and hence the capillary number Ca, which is the main parameter controlling the strength of the instability. Secondly, the inflation of the elastic sheet implies that the interface propagates into a converging cell gap, which introduces additional stabilizing effects compared with a classical parallel-walled cell, even at the same value of Ca. Our 'rigid-wedge' approximation reveals that the elasticity of the top wall affects the instability mainly in this indirect way via the change in cell geometry, while the non-axisymmetric perturbations to the elastic wall can be neglected.
The following physical mechanisms affect the instability, and all have stabilizing effects. The films of liquid deposited on the cell walls behind the advancing interface reduce the viscous driving force of the instability. Surface tension interacts with the horizontal (in-plane) curvature as well as the radial variation of vertical (transverse) curvature in the converging cell, and provides a restoring force. The kinematic compression in the base flow (i.e. the presence of a negative radial velocity gradient) advects the perturbed interface back to its mean position; it is due to the diverging radial flow, the converging cell gap and the lifting of the sheet. In the parameter regime covered by this study, all of these mechanisms have a significant effect on the stabilization of the system. In general, a more compliant wall results in a larger change in the cell geometry, and hence yields a smaller capillary number Ca and a larger slope α i , both of which have a stabilizing effect on the system.
The linear stability analysis does not apply to the nonlinear late-time evolution observed in the experiments. However, the wavenumbers observed in the experiments agree well with those predicted by the linear theory. This suggests that they are set by the linear evolution at the early stages of development and are not significantly altered by the nonlinear effects. Similar conclusions were reached by Ducloué et al. (2017a), who found that the number of fully grown fingers in an elastic-walled channel followed the scaling predicted by the linear stability analysis of a planar front in the Saffman-Taylor problem (Saffman & Taylor 1958). Unlike the saturated patterns found by Ducloué et al. (2017a), the fingers in the radial cell continuously evolve as the cell geometry (such as interfacial circumference, slope of the sheet and depth of the cell) changes, ultimately resulting in the decay of the fingers. However, the maximum length of fingers during their evolution does reach a saturation in the control parameter Q. As in the study by Ducloué et al. (2017a), this saturated length varies approximately linearly with the initial depth of the layer and scales with the length of the liquid wedge ahead of the interface.
The present study serves to bridge the gap between the fingering instability in a rigid-walled Hele-Shaw cell and the printer's instability in coating flows. While the similarity between the two has been noted before (Couder 2000), we are not aware of any quantitative studies of pattern formation relevant to both classical instabilities, apart from the recent work by Ducloué et al. (2017a). As the stiffness of the elastic sheet increases, the stabilizing taper of the cell decreases and we expect the short and stubby fingers to eventually transition into the highly branched fingers observed in classical viscous fingering.

Appendix A. Experimental protocol
Here we present more details of the experimental protocol described in § 2.1. The set-up is shown schematically in figure 1(a,b). The bottom boundary of the cell was a 15 mm thick float-glass plate accurately levelled to within 0.1 • . The upper boundary of the cell was a square elastic sheet made from latex (Supatex) or polypropylene (Film Products Ltd), and the initial depth of the cell was set by a spacer with a circular cutout, which was made from the same type of latex sheet. Their properties are listed in table 1. The fluids were injected from below through a small nozzle with an internal diameter of 2.28 mm embedded in a brass cylinder with a diameter of 19 mm, which was mounted flush with the bottom boundary and appears as the central black circle in figures 1(b,c), 2(a), 3(c) and 11.
Prior to each experiment, the spacer was laid down on the glass plate and covered with the elastic sheet. The cell was then filled with a viscous fluid. For the experiments with latex sheets, silicone oil from Basildon Chemicals Ltd was used, because it wets both glass and latex surfaces. However, we found that if silicone oil was used in experiments with polypropylene sheets, significant electrostatic forces developed resulting in an erratic distortion of the gas-liquid interface. To avoid this, a mixture consisting of 90 % glycerol and 10 % purified water (by volume) was used for these experiments. This mixture only partially wetted the polypropylene sheet, but the expanding bubble deposited a continuous film of viscous fluid on the cell boundaries, and dewetting was only observed on time scales longer than that of a typical experiment.
During the filling procedure, a heavy glass plate was placed on top of the elastic sheet to keep it from deforming, thus ensuring that the fluid layer had a uniform initial depth, h 0 , set by the thickness of the spacer separating the top and bottom boundaries of the cell. The fluid occupied a circular region of radius R fluid < R cell , in order to leave an air gap into which the fluid could propagate, with the displaced air escaping through the unsealed sides of the cell between the latex separator and the elastic sheet resting on top. (In practice the experiments were terminated long before any significant motion had occurred near the edge of the fluid domain, and we found that changing R fluid from 150 to 170 mm had no effect on the experimental results.) At the end of each experiment, the cell was taken apart and cleaned carefully before the whole procedure was repeated again. All experiments were performed in a temperature-controlled laboratory at 21 • C.
It was necessary to inject fluids through the bottom boundary to avoid interfering with the elastic wall, so some viscous fluid inevitably leaked during the filling procedure and clogged the injection nozzle. Therefore, in order to reduce the initial perturbations to the system, a small circular gas bubble of radius R init ≈ 5 mm was injected into the viscous fluid layer through the inlet and allowed to relax before each experiment was initiated.
The displacement flow was driven by injecting nitrogen from a compressed gas cylinder via a mass airflow meter (Rd-Y Smart Meter PCU1000, Icenta Controls Ltd), a fine needle valve and a three-way pneumatic solenoid valve. Before each experiment, the three-way valve was set to discharge the nitrogen into the atmosphere while the fine-needle valve was adjusted manually to set the flow rate Q to a desired value indicated in table 1. Then, the three-way valve was switched to divert the flow into the cell to start the experiment. Careful initial adjustment allowed us to achieve a value of Q within 3 % of the target value, and the flow rate remained constant to within 0.5 % during the experiment.
The evolution of the interface shape was recorded by a CCD camera (JVC KY-F1030, 1360 pixels × 1024 pixels, 7.5 frames per second, for the experiments with latex sheets; Dalsa Genie HM1400, 1400 pixels × 1024 pixels, 12.5 frames per second, for the experiments with polypropylene sheets) looking down on the experiment. The polypropylene sheets were transparent and allowed a clear view of the position of the interface, while backlighting was used to obtain a clear image through the translucent latex sheets. While the time intervals between frames were known, the exact start time of each experiment was not recorded. Instead, it was extrapolated from the images before and after the start of the experiments. Therefore, the maximum error on the start time was 0.13 s in the experiments with latex, and 0.08 s in the experiments with polypropylene. Each time series of top-down images showing the gas-liquid interface was then processed to extract the evolution of various quantities with t andR as discussed in § 2.2. variable spacing using a conservative finite-difference scheme, and the grid is adapted regularly to ensure that the spacing at each point is at most 1 % of the length scale on which the solutions are varying. Time stepping is performed using the implicit backward-Euler method: Given the solution at the current time step, the solution at the next time step is found using Newton iteration. (Each iterative step involves solving a banded matrix equation, which is computationally inexpensive.) For each time step, one full step is taken in parallel with two half-sized steps, and the step size is adjusted to ensure that the discrepancy between the two results is at most 1 %. The solutions to the perturbation equations (B 2) are calculated in parallel with the base-state solutions for 1 n 50, and, due to linearity, the Newton iteration always converges after a single step. This numerical method was verified by direct comparison to the (different) finite-difference method presented in Pihler-Puzović et al. (2013). For the frozen-time approximation, the solution to (B 1) is held fixed while, for each n, (B 2) is continually evolved forward in time (while the solution is regularly rescaled to keep it O(1)) until the structure of the solution no longer changes, indicating that the dominant eigenmode has been found (and all eigenmodes with smaller growth rates have become negligible).