Experimental Observation of Dissolution Finger Growth in Radial Geometry

Reaction-inﬁltration instability refers to the morphological instability of a reactive ﬂuid front ﬂowing in a soluble porous medium. This process is important for many naturally occurring phenomena, such as the weathering and diagenesis of rocks, dissolution in salt deposits and melt extraction from the mantle. This paper is focused on experiments on dissolution ﬁnger growth in radial geometries in an analog fracture. In the experiments, pure water dissolves a plaster sample forming one of the fracture walls in a Hele-Shaw cell with controlled injection rate and aperture. The ﬂow is directed inwards to the center, and we observe the reaction-inﬁltration instability developing along the relatively long perimeter of the plaster. Our experimental results show a number of features consistent with the theoretical and numerical predictions on the ﬁnger growth dynamics such as screening and selection between the ﬁngers. Statistical properties of the dissolved part evolution with time are also investigated.


INTRODUCTION
In geological systems, dissolution plays an important role in the weathering and diagenesis of Earth's rocks [1,2], chemical erosion of salt deposits [3,4], and melt extraction from the mantle [5]. It is also of fundamental importance in many engineering applications, including dam stability [6] and CO 2 sequestration [7]. The important applications in the oil industry include acidization of petroleum reservoirs [8] in order to enhance oil and gas production by increasing the permeability of the rock [9,10].
Because the reaction-infiltration instability plays an important role in a variety of fields, it stimulates dissolution-related research projects, both theoretical and numerical. Linear stability analysis can be applied to characterize the initial instability in porous media dissolution [11][12][13][14], but after fingering develops, we enter a nonlinear regime, where very few theoretical tools can be applied and one needs to resort to numerical simulations [3,8,[15][16][17][18][19][20][21]. Compared to the theoretical and numerical works on the subject, there are relatively few experimental studies, especially on the observation of dissolution in quasi-2D radial geometry, which is the focus of this article.
In the lab experiments, two different setups are usually used: rock core acidization in Hassler cell and quasi-2D systems in Hele-Shaw cell which are aimed to study the dissolution in quasi-2D porous media and fractures. The number of core-flooding experiments reported in the literature is significantly larger than Hele-Shaw cell studies. There are two reasons for that: First, the core-flooding is closer to the real conditions encountered in the wellbore acidization in petroleum industry. Second, in core flooding it is easier to inject the flow into a porous matrix with negligible boundary effect. Many variables are systematically controlled in core-flooding experiments, including: injection rate [22], sample material [23], system scale [24], pH and temperature [25].
In this work, we have decided to use a 2D Hele-Shaw cell, as quasi-2D systems are easier to visualize and due to a large number of numerical work performed on these systems [18,21,26,27]. On the experimental side, Daccord et al. [28,29] investigated the water/plaster system in radial Hele-Shaw cell with central injection, and Golfier et al. [16] used a water/salt system in rectangular Hele-Shaw cells. The intention of both works was to study flow in a porous matrix, however, because of the difficulty in avoiding wall effects at boundaries, the injected water in most cases was found to flow along the difficult-todetect aperture between the medium and the confining cell plate. This problem can be turned into advantage, if we promote it in a controlled manner instead of avoiding the wall flow. Such a controlled-aperture system can then be considered as an analog of a fracture, and the study is then directed at the investigation of how the fracture aperture evolves in time as a result of dissolution. Within this frame, Detwiler et al. [30] undertook a well-controlled dissolution study of the water/KDP system in rectangular Hele-Shaw cell and systematically measured the evolution of aperture at different flow rates. Osselin et al. [31] have performed experiments on the onset of reactiveinfiltration instabilities in a fracture with a microfluidic setup using a rectangular water/plaster system. However, in many cases the relevant geometry is radial rather than rectangular, for instance, in the oil industry where the acid fluids are injected from a well, and in groundwater protection where pollutants expand with or without dissolution radially from the pollution source. Xu et al. [32] have recently studied dispersion in fractures in radial geometry with a dissolution pattern around the inlet.
The aim of this project is to study the dissolution finger growth in a fracture aperture of radial geometry. In section 2, we describe our experimental setup. In section 3, we present our experimental results discussing the screening effect and the statistical properties of the dissolved part evolution. The conclusions are drawn in section 4.

DESCRIPTION OF EXPERIMENTS
The experimental scheme is illustrated in Figure 1. A Hele-Shaw cell is formed by two circular glass plates which are separated by 1 mm aluminum spacers and held together by clamps. The bottom glass plate (diameter d 1 = 36.0 cm) is larger than the upper one (diameter d 2 = 25.0 cm) and has an external rim to hold the water surface at a fixed altitude level. There is an outlet at the center of the lower glass plate. A lightbox with a homogeneous intensity of light illuminates the system from below. A digital camera (Nikon D7100) records the sample images from the top every 5 min, and the whole dissolution process is thus recorded from the beginning of fluid withdrawing up to the dissolution channel breakthrough at the central outlet. The circular plaster sample between the two plates was prepared as follows: A gypsum saturated water solution was first injected into the Hele-Shaw cell. We made a plaster paste by mixing water and plaster powder with the ratio 2:3 by weight. This paste was then injected from the center of the Hele-Shaw cell. The paste displaced the plaster saturated water and formed a circular plate of radius R 0 = 8.0cm. The hydration of this circular plaster paste requires approximately one hour to complete. During the plaster hydration process, the plaster paste was kept in the cell surrounded by saturated water. Over time, a form of segregation called bleeding takes place, where some of the water in the plaster tends to rise to the top surface of the plaster plate [33]. This process creates a small gap h 0 = 50µm above the upper surface of the plaster. After the completion of the hydration process, we removed the top glass plate, put several plastic films (each film thickness h 1 = 100µm) on the aluminum spacers and put back the top glass plate. The artificial aperture in our experiment h is defined as the distance between the upper glass plate and the surface of the plaster sample. The aperture created in this way is thus h = h 0 + n · h 1 where n is the number of plastic films. In the experiments reported here, we used n = 2 films which gave h = 250µm. When the sample preparation was completed, we started the dissolution experiments.
Because the radial dissolution by injection from the center gives a very short dissolution front around a point-like inlet, it becomes difficult to analyze the evolution of the fingers and the periodic wavelength from the experimental images in such a setup. Therefore, we chose instead to withdraw the water by a syringe pump from the outlet located at the center of the bottom glass plate. In this way, the freshwater flows from the rim toward the center and dissolves the plaster sample from the outer boundary. The instability can then be observed along the external perimeter of the sample. The withdrawing flow rate Q is set as Q = 0.18ml/min and the initial aperture of the artificial fracture is h = 250µm. The permeability of the porous plaster matrix is κ p = 6.0 · 10 −14 m 2 [34] but the permeability of the fracture FIGURE 2 | Experimental photos of the developing fingering pattern at different moments of time. The time interval between the photos is 2 days. In the circular plaster sample, the dark yellow disk is the undissolved part of plaster and the light yellow part represents the dissolved or partially dissolved part. calculated as κ = h 2 /12 is around 5.2·10 −9 m 2 , thus 5 magnitudes larger than the permeability of the porous matrix of our plaster sample. Therefore almost all the freshwater flows through the fracture instead of the porous matrix. The freshwater pumped through the system is distilled water at room temperature T = 22 • C with pH = 7.17. The molecular diffusion coefficient of plaster (gypsum) in water is D m = 1.0 · 10 −9 m 2 /s [35]. One experiment lasts around 10 days. The camera records the entire dynamic evolution process from initial instabilities to fingering formation, then to dissolution finger growth, and finally to the breakthrough of the longest fingers at the outlet.

Characteristic Timescales
The initial aperture h is an important characteristic length scale in our system. A characteristic timescale for diffusion across the aperture is t D = h 2 /Dm = 62.5s and the characteristic timescale for the reaction on the same length scale is t R = h/k = 54.3s, where k is the chemical reaction kinetic constant (k = 4.6 · 10 −6 m/s) [36]. A relevant time scale for convection is the time it takes to flush the system, i.e., t C = πR 2 0 h/Q = 1, 670s. As we see, t C ≫ t D ≈ t R . It means that reaction and diffusion across the aperture happen almost immediately compared with the time it takes for a fluid particle to flow through the system and we therefore expect that the calcium concentration of the water will reach the saturation concentration at the outlet. We performed density measurements to determine the concentration of the effluent solution. The concentration at the outlet is C outlet = 2.5g/L which is consistent with the value of the saturation concentration reported in the literature C sat = 2.53g/L [37].

Dissolution Finger Growth With Screening Effect
As the freshwater flows from the edge of the plaster sample to the center, the plaster begins to dissolve. In this process, the water becomes saturated, thus the dissolution concentrates at the perimeter. After several days, a visible dissolution front appears and it slowly develops into many dissolution fingers as a result of the reactive-infiltration instability. The further growth of these fingers becomes nonlinear, with scarce theoretical results concerning their shapes or growth rates [38,39]. On the other hand, numerical models give a number of predictions for the fingering which can be qualitatively compared with experiments [17,18,21]. In particular, one can study the screening between the fingers, with the longer ones suppressing the growth of their shorter neighbors. As a results, approximately half of the active fingers continue to grow while the other half cease to grow. The process then repeats itself, leading to the scale-free distribution of finger lengths [40,41].
Such a hierarchical growth of the fingers is clearly observed in the experimental images. Four experimental photos are chosen to show the process of the finger growth in Figure 2 (see also Figure S1 and Videos S1, S2).
A dissolution front is extracted from the experimental images by using thresholding. A front position D is calculated as the distance between the point at the dissolution front and the outlet center. The front position varies with the polar angle θ and time t and we define the dissolution length as l(θ , t) = D(θ , t = 0) − D(θ , t), where D(θ , t = 0) is the initial front line. The initial front position D(θ , t = 0) has a small variation with the radius of the plaster sample R because in the experiments, the initial plaster sample is not perfectly circular. The dissolution length evolution with time l(θ , t) is shown in Figure 3. We should notice that the dissolution length function l(θ , t) is not a unique function because one polar angle θ could correspond to more than one dissolution front point, at different radial positions. This is because fingers develop along the radial direction and can also grow wider, with an orthoradial growth component.
In Figure 3, the competitive growth of the fingers is observed, due to the screening effect. The time interval we choose is fixed (1.0 day between neighboring curves), but the growth rate for different fingers varies significantly. One part of the fingering pattern (with θ in between 4 and 5 radians) grows more slowly than the rest over the entire time which means that the fluid flux through this part must have been significantly smaller. Such a flow inhomogeneity could be accidental, because the geometry of the system is not perfectly uniform. It leads to an initial circular symmetry breaking, where one side becomes a freshwater preferential flow path, screening off dissolution at other sides. Due to a positive feedback loop, eventually, some fingers will dominate so that most freshwater concentrates in  these dissolution fingers. These long fingers then continue to grow while the short fingers grow very slowly (see Figure 3B). For the dissolution front propagation with time in a repeated experiment (see Figure S2). The longest fingers grow almost exponentially with time judging from Figure 4A, where we fit the finger length vs. time with an exponential function l = a·e βt . The exponential growth rate β is shown in the legend of Figure 4A. Interestingly, the values of exponential growth rates are the same for all the three different fingers β = 0.32day −1 , which shows that the longest fingers grow largely independent of each other. From Figure 3A, the distance between these longest fingers is comparable to the length of their lengths. This is in agreement with the observation [41,42] that the long finger screens the area of a lateral extent approximately equal to its length. For the analysis of the growth rates of other fingers (see Figure S3A).
In order to perform the statistical analysis of the dissolution finger growth, we find the local maxima of the curves in Figure 3B and define the dissolution finger length as the dissolution length corresponding to these maxima. The dissolution finger length histogram is shown in Figure 4B. The bins are chosen from 0 to 4.0 cm with a bin size of 0.5 cm, which divides the fingers into several types (orders). At an early stage, a lot of small fingers appear, with the length below 0.5 cm; these are the first order fingers. After 3 days, 2nd order fingers become clearly visible and the distribution gets wider as time progresses. With time, the number of small fingers (finger length below 1 cm in Figure 4B) decreases significantly since they are absorbed by the moving dissolution front while the longest fingers (finger length above 2.5 cm in Figure 4B) continue to grow. For the dissolution finger length histogram of other fingers (see Figure S3B).

Statistical Properties of the Dissolution Patterns
We observe from Figure 2 that the interface between the dissolved and the undissolved part is diffuse as the thickness of plaster sample at the dissolution front gradually transitions from the dissolved part to the undissolved part. Therefore, thresholding of the image, although useful for analyzing the dissolution patterns, leads inevitably to the information loss, as partially dissolved region are either interpreted as fully dissolved or as undissolved. In order to measure the aperture variation in the experimental images, a calibration between the thickness profile of the plaster sample and the intensity profile of the photos is performed by an X-ray Thickness Gauging REX-CELL 4X from Flow Capture [43][44][45].
The X-ray measurement for the thickness of dissolved plaster sample is displayed in Figure 5. The data records the photon counts at different positions of the sample by the X-ray measurement. Then the thickness profile is obtained from the intensity data by d = d s · ln(I/I 0 ) ln(I s /I 0 ) , where d s is the initial thickness of the undissolved plaster sample, I is the intensity at a given pixel, I 0 is the average intensity value of the background and I s is the average intensity value of the undissolved part. This formula is derived based on the Beer-Lambert law, I = I 0 e −ζ d [46] where ζ = 0.43mm −1 normalized by a reference intensity I s = I 0 e −ζ d s corresponding to the undissolved sample. The X-ray measurement and the calibration with the images illuminated by white light are shown in Figure 6. We compared the image illuminated by X-ray and the image illuminated by lightbox, we found a linear regression for the scatter plot of two images with y = 0.94 · x − 0.02. The Pearson correlation coefficient [47] between the two images is 0.95, which confirms that we can calculate the thickness of the plaster based on the lightbox measurements.
This allows us to quantify the local volume of dissolved gypsum in the sample. First, we define the aperture growth at a given point as as δh t (r, θ ) = h t (r, θ ) − h where h t (r, θ ) is the aperture at radial position r and polar angle θ and h is the initial aperture. Next, we calculate the dissolved volume per unit angle as V θ (θ , t) = R 0 0 δh t (r) · rdr. Note that the total dissolved volume is an integral of V θ (θ , t) over the polar angle, V diss (t) = 2π 0 V θ (θ , t)dθ . Subsequently, we find it more convenient to use the arc length along the perimeter, p = R 0 θ , instead of θ itself,  to parameterize the experimental data. Then we define a local dissolved part as S D (p, t) = The calculation method is illustrated in Figure 7.
The local dissolved part S D (p, t) evolves in time as the dissolution front propagates, as illustrated in Figure 8A. In order to analyze this function quantitatively, a fast Fourier transform with a Blackman window [48] is applied at different times to obtain power spectrum S(k, t) where k is the spatial frequency or wavenumber (k = 1/λ where λ is the wavelength).
The data in Figure 6A shows another manifestation of the screening of the shorter fingers by the longer ones. The power spectrum in Figure 8B has two main peaks before wavenumber k = 0.1cm −1 . The first peak (Peak 1 indicated by the red circle) at 0.02 cm −1 with wavelength ∼ 50cm comes from the perimeter of the plaster sample and the second peak (Peak 2 indicated by the blue circle) at 0.06 cm −1 with a wavelength ∼ 17cm is connected with the deviation of the initial gypsum disk from the circular shape (see Figures S6, S7 in Supplementary Data). We will ignore these two peaks because they come from the geometric properties of the initial plaster sample and not from the dissolution process. Beyond these two peaks, the maximum of the power spectrum is observed at 0.24cm −1 with wavelength ∼ 4.2cm, indicated by the red star in Figure 8B. This wavelength is related to the characteristic distance between the longest fingers, which is the main contribution to the power spectrum after finger formation. From the semi-log representation in Figure 8C, we see that the . Three main peaks are indicated in the plot. The first two (marked by red and blue circles) are related to the initial geometry of the plaster sample. On the other hand, the peak at 0.24 cm −1 with wavelength ∼ 4.2cm (marked by the red star) is related to the characteristic wavelength of the fingering pattern. (C) Semi-log plot of the power spectrum, which allows to observe power amplitude trends in a larger scope. The amplitude decays almost exponentially with wavenumber at low wavenumber (high wavelength), followed by a flat plateau at higher wavenumbers. The solid black line in (C) is the linear fit, the slope of which gives the characteristic decay length, l D = 0.55cm. amplitude decays almost exponentially with wavenumber at low wavenumbers (high wavelengths), followed by a flat plateau at higher wavenumber. The decay is exponential, S(k, t) ∝ S(0, t) · e −l D ·k with the characteristic decay length l D = 0.55cm as shown in Figure 8C. The Fourier transform of a Lorentzian gives an exponential function [49]. The width of the Lorentzian gives the characteristic decay length of this exponential function. Since the largest fingers are of a similar shape, we expect the decay length l D = 0.55cm to correspond to their characteristic width, which is indeed the case. The amplitude of the power spectrum decays gradually after a crossover at wavenumber k = 20cm −1 corresponding to the wavelength λ = 0.5mm. The part of the spectrum with the wavenumber larger than k = 20cm −1 are considered as noise from the roughness of the dissolution front.
The local dissolved part increases in time, as illustrated in Figure 9. The uncertainty of the measurement of the thickness variation is about 5% according to Figure 6. In Figure 9, we present the evolution of S D (p, t) with error bars at several different points along the perimeter of the sample, the positions of which are marked by circles with the corresponding color in Figure 9A. Figure 9B further confirms the crucial role played in the dynamics by the competition between the fingers-the longer fingers speed up at the expense of the shorter ones.
To validate the calibration between the aperture profile and intensity profile, we have also calculated the growth of the total dissolved volume in time. The results, presented in Figure 9C, show that V D increases linearly in time with the slope of 528mm 3 per day. This is close to the theoretical estimate of the growth of the dissolved volume based on the mass balance of the reactant where Q is the injection rate, c sat is the solubility of gypsum in pure water at 20 • C and ρ b is the bulk density of experimental plaster sample. For Q = 0.180ml/min, c sat = 2.53g/l [37] and ρ b = 1.1g/ml we get the rate of growth of the dissolved volume of 596mm 3 /day with an uncertainty of 9.1%, a reference line of this slope with error bars is shown in Figure 9C.
The possible factors influencing the difference between the theoretical slope and the one measured experimentally come from inaccuracy in theH measurement of the bulk density ρ b , measurement of the thickness variation δh and the linear fitting for V D (t) based on the image analysis. Importantly, the second run of this experiment results in similar statistical properties of the dissolved pattern (see the Supplemental Data in Figures S4, S5). We now introduce a calculation method by Vinningland et al. [50] who characterized the mean wavenumber of a growing interface. We again ignore the low frequencies since they are related to the initial geometric properties of the sample and only consider the wavelength λ < 10cm i.e., k > 0.1cm −1 . An average wavenumber < k > from the power spectrum S(k) for k > 0.1cm −1 is defined as: where the sum is over all k > 0.1cm −1 . The standard deviation σ f is defined as: We calculate the temporal evolution of < k > and σ k , shown in Figure 10. As observed, the average wavenumber decreases with time, i.e. the average wavelength increases. This is another manifestation of the competition between the fingers: longer fingers develop and screen off nearby shorter fingers, which leads to a larger wavelength. The average wavenumber ranges from 0.4 to 1.8 cm −1 , thus the average wavelength ranges from 0.6 to 2.5 cm which fits the observation of characteristic wavelength from the experimental images.

CONCLUSIONS
Understanding of the dissolution of fractures is important both for basic science (e.g., for studies of speleogenesis) as well as for technological applications, particularly in petroleum industry. However, the dissolution experiments in quasi-2d, radial geometry are relatively seldom performed. We dissolve a plaster disk in a Hele-Shaw cell by withdrawing water from the center, thus creating the inward flow pattern, and for the first time to our knowledge, we report the reactiveinfiltration instability and finger growth along the perimeter of the circular plaster sample. The perimeter is 50.3 cm long which is sufficient to perform a statistical study of the reactive-infiltration instability and the dissolution fingers growth with time. We observe strong competition between the growing fingers with the longest fingers growing exponentially with time in our experiments. We measure the thickness variation of the plaster sample by X-ray gauge and quantify the evolution in time of the dissolved volume. We then analyze it using statistical methods. The characteristic wavelength of the perturbations of the front is measured by a fast Fourier transform of the local dissolved volume. The power spectrum shows exponential decay with a characteristic decay length l D = 0.55cm. On the other hand, the average wavelength increases linearly with time.
Our experimental setup allows us to adjust both the aperture thickness and flow rate. In the future, we will systematically control these two parameters to study a dissolution phase diagram in a radial geometry.

DATA AVAILABILITY
The raw data supporting the conclusions of this manuscript will be made available by the authors, without undue reservation, to any qualified researcher.

AUTHOR CONTRIBUTIONS
LX and KM designed the experiment. LX performed the experiments, analyzed the data, and authored the paper. PS, RT, EF, and KM assisted with the interpretation and data analysis, and editing of the manuscript.

ACKNOWLEDGMENTS
We thank Atle Jensen, Bin Hu for technical support of X-ray measurement and Marcel Moura, Fredrik K. Eriksen, Mihailo Jankov, Monem Ayaz, Guillaume Dumazer, and Florian Osselin for useful discussions.

SUPPLEMENTARY MATERIAL
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphy. 2019.00096/full#supplementary-material Video S1 | The video of the experiment described in the main text shows the sequence of the experimental photos of the developing fingering pattern at different moments of time.
Video S2 | The video of the experiment described in the Supplementary Material shows the sequence of the experimental photos of the developing fingering pattern at different moments of time.