4D Neutron Imaging of Solute Transport and Fluid Flow in Sandstone Before and After Mineral Precipitation

In many geological systems, the porosity of rock or soil may evolve during mineral precipitation, a process that controls fluid transport properties. Here, we investigate the use of 4D neutron imaging to image flow and transport in Bentheim sandstone core samples before and after in‐situ calcium carbonate precipitation. First, we demonstrate the applicability of neutron imaging to quantify the solute dispersion along the interface between heavy water and a cadmium aqueous solution. Then, we monitor the flow of heavy water within two Bentheim sandstone core samples before and after a step of in‐situ mineral precipitation. The precipitation of calcium carbonate is induced by reactive mixing of two solutions containing CaCl 2 and Na 2 CO 3 , either by injecting these two fluids one after each other (sequential experiment) or by injecting them in parallel (co‐flow experiment). We use the contrast in neutron attenuation from time‐resolved tomograms to derive three‐ dimensional fluid velocity field by using an inversion technique based on the advection‐dispersion equation. Results show mineral precipitation induces a wider distribution of local flow velocities and leads to alterations in the main flow pathways. The flow distribution appears to be


Introduction
The transport and mixing of reactive fluids through porous media can modify the porosity structure in geological systems (e.g., Baqer & Chen, 2022;Hommel et al., 2018;Kudrolli & Clotet, 2016).A prime example is porosity reduction and clogging due to mineral precipitation (e.g., Dávila et al., 2020;Emmanuel & Berkowitz, 2005), which is commonly observed in the sealing of faults and fractures by pore fluids (e.g., Fisher & Knipe, 1998), microbial activity leading to bio-clogging (Thullner et al., 2002), as well as in industrial applications, such as effluent disposal (Eriksson & Destouni, 1997).However, characterizing the interplay between fluid flow transport and mineral precipitation is challenging as it occurs across a broad range of temporal and spatial scales (Bray et al., 2017;Emmanuel & Berkowitz, 2005).
Among minerals, calcite precipitation (CaCO 3 ) has attracted much attention for its importance in permanent geological storage of CO 2 (e.g., Chang et al., 2017), enhanced oil recovery operations (Wu et al., 2017), geothermal energy recovery (Simmons & Christenson, 1994), and soil stabilization (Mujah et al., 2017).Efforts to understand the kinetics of calcite precipitation have led to the derivation of several rate equations for a range of experimental conditions (Mitchell et al., 2019;Plummer et al., 1979;Zhang & Dawe, 1998).Calcite nucleation and growth at the pore scale depend on several parameters such as pH, temperature, and solute concentration, and require descriptions of coupled processes over a wide range of scales (Noiriel & Soulaine, 2021).Direct imaging of calcite precipitation and solute transport inside porous rocks remains a bottleneck to understand these processes.Since rocks are not transparent to visible light, imaging their interior requires specific techniques that utilize radiations capable of penetrating a rock and providing information about its internal structure.
Recent developments in 4D imaging techniques at the pore to core scales, such as neutrons and X-rays microtomography with spatial resolutions of a few to tens of micrometers offer unprecedented opportunities for highresolution continuum monitoring and modeling of reactive transport in porous rocks (e.g., Mancini et al., 2020).In the past 20 years, X-rays have been widely used to obtain high spatial and temporal resolution of threedimensional microstructures of solid matrix in porous media (e.g., Desrues et al., 2010;Renard et al., 2004).Standard X-ray sources have a much lower attenuation coefficient to water compared to the solid phase of rock; thus, usually, a contrasting agent, such as iodine is used to image the fluid (Andrew et al., 2015;Kim et al., 2013;Scanziani et al., 2020;Vinegar, 1986).Moreover, because X-rays are also adsorbed by the solid matrix, the separation of the fluid from the solid, that is, the segmentation between these phases, may be challenging.Neutron imaging, on the other hand, is very sensitive to a few elements, including hydrogen and cadmium, allowing imaging of these elements with sufficient contrast, whereas most other elements are almost transparent to neutrons (Anderson et al., 2009).Because neutron imaging is sensitive to hydrogen atoms (protium and deuterium), this technique can image minimal amounts of aqueous water in the pore space, and water contained in hydrated minerals, such as clays (Perfect et al., 2014), without imaging the surrounding solid matrix.Conversely, in the case of heavy water where deuterium replaces both hydrogen atoms, the contrast to neutrons is much lower.
Neutron imaging has developed rapidly in the past few years, providing new opportunities for non-destructive characterization of a wide range of geomaterials (e.g., Tengattini et al., 2021).Several recent studies have demonstrated the capabilities of time-resolved neutron imaging to study fluid flow in rocks, for example, to reveal the relationship between wormhole advancement and permeability variation in dissolving rocks (Cooper et al., 2023), or to investigate the interactions between imbibition and pressure-driven flow in microporous deformed limestone (Lewis et al., 2023).Neutron imaging has also been used to study transport properties in other porous media, such as the migration of condensed water vapor through cracked concrete (Gupta et al., 2022), hydro-channel porous transport layer designed for unitized reversible fuel cells (Komini Babu et al., 2023), as well as salt precipitation and water transport in CO 2 electrolysis (Disch et al., 2022).These studies demonstrate the versatility of the neutron imaging technique and its ability to provide a deeper understanding of fluid transport in various porous media.
Several laboratory experiments of carbonate precipitation/dissolution in porous media during reactive transport have highlighted the importance of the dynamics of porosity and permeability evolution in both the pore scale and the Darcy scale (Baqer & Chen, 2022;Beckingham, 2017;Garing et al., 2015).Calcite crystal nucleation and growth in model porous media and fractures have been investigated using X-ray microtomography imaging (Noiriel et al., 2012(Noiriel et al., , 2021)).Experiments with flow-controlled conditions in a Hele-Shaw cell have been performed to observe various calcite precipitation patterns controlled by different injection flow rates (Balog et al., 2019;Izumoto et al., 2022;Schuszter et al., 2016).While a variety of spatial distributions of calcite minerals leading to porosity and permeability evolution have been observed in previous studies, assessing the impact of mineral precipitation on fluid flow properties inside an intact rock matrix remains an important challenge and requires in-situ data.
The main objectives of this study are to use neutron imaging techniques to study fluid flow in porous sandstone and to quantify how the precipitation of calcite into the porous medium changes the flow and transport properties.Here, we aim to (a) demonstrate the feasibility of neutron imaging to explore the mixing of two miscible aqueous fluids during co-injection into porous sandstone and (b) investigate the impact of calcite precipitation on the fluid velocity distribution in porous rock.We use 4D neutron microtomography imaging to track fluids injected into core samples of Bentheim sandstone.Then, we derive the three-dimensional flow velocity field by solving an inverse problem based on the advection-dispersion equation and the contrast in neutron attenuation of the sequential images.Results show that neutron imaging can be used to measure the variability of flow properties in porous rock and quantify how the flow velocity field is modified after calcite precipitation.We also used X-ray imaging at different resolutions, showing the calcite precipitate distribution inside the rock, although it was not directly comparable to the velocity images since the unconsolidated precipitate was transported by flow between the neutron and X-ray scans.

Experimental Set-Up and 4D Imaging
We performed three flow-through experiments on cylindrical core samples of Bentheim sandstone, 5 cm diameter and 5 cm height, cored perpendicularly to sedimentary bedding from the same block.This reservoir rock contains mainly quartz grains with sizes in the range 200-400 μm and the pore sizes are in the range 50-500 μm (Klein et al., 2001).The rock has a porosity of 21% and is considered to be homogeneous in the continuum scale.Experiments were performed at the ambient temperature of 24°C.They are labeled BS0 for the fluid tracer test, and BS1 and BS2 for the calcium carbonate precipitation experiments.
We built core holders made of polytetrafluorethylene (PTFE) and a jacket of fluorinated ethylene propylene surrounds the rock sample within the core holder.The chosen materials do not contain hydrogen atoms and are therefore almost transparent to neutrons.The core holder contains two inlets connected to a syringe pump injecting fluid at a given flow rate from either the top extremity of the core sample (experiments BS0 and BS1) or its bottom extremity (experiment BS2).An outlet at the opposite side of the core holder collects the fluid after flow-through in the sample and is at atmospheric pressure (Figure 1a).Porous diffusers made of sponge or porous glass were used at the inlets and outlet points to diffuse the fluid, allowing almost homogenous injection and outflow conditions over the circular surface area of the sample.At the inlet side, a barrier separates the two injected fluids in co-flow experiments, preventing their mixing before reaching the sample.The core holder was installed on the rotating stage of the ICON neutron imaging beamline of the SINQ neutron spallation source (Paul Scherrer Institute, Switzerland), with a source neutron flux of 10 14 n/cm 2 /s (Kaestner et al., 2011).To minimize the penumbra effect arising from the slight divergence of the neutron beam in blurring the image at the geometric irregularities, we used a Gadox screen on the camera box.By using this device, we can estimate the penumbra blurring at varying distances from the detector, and together with the lower Gadox screen thickness, the spatial resolution can be improved (Kaestner, Kis, et al., 2017).
A laboratory X-ray tomograph was installed on the neutron beamline at the Paul Scherrer Institute (Kaestner, Hovind, et al., 2017) and enabled X-ray imaging at different stages of experiment BS2 (Figure 1b).Throughout the rotation of the stage during the tomography acquisitions, the injection points (inlets) remained stationary within the core holder, which itself underwent full rotation.Time-resolved neutron three-dimensional tomography images of the rocks were acquired with a voxel size of 62 μm, and each radiograph of the samples contained 800 × 800 pixels.The voxel size is on the same order of magnitude as the mean sandstone pore size (50-200 μm) and therefore the neutron measurements are averaged over around one pore size.The three-dimensional tomograms were obtained from 626 projections covering the 360-degree rotation of the core holder, and the average acquisition time of a scan was 32 min.The X-ray Along with all the acquired scans, some dark field, and several flat field images were acquired as well.The neutron and X-ray volumes were reconstructed with scattering correction using the open-source software MuhRec developed by Kaestner (2011).For each reconstruction, the intensity values were normalized using the dark and flat field images, as described by Kaestner et al. (2008) to represent the normalized attenuation of the beam instead of an absolute intensity reading.Sample BS2 after precipitation was also imaged using X-ray microtomography at the European Synchrotron Radiation Facility, 1 year after the experiment was performed.This duration was necessary to let the radioactivity of the sample decay sufficiently enough after imaging it with neutrons.During this year, the sample, saturated with heavy water, was stored horizontally in the experimental hutch of the SINQ neutron beamline at the Paul Scherrer Institute.At the European Synchrotron Radiation Facility, we acquired a high-resolution X-ray tomogram of the entire sample on beamline BM18 with a voxel size of 5.85 μm, allowing imaging details within pores.Although high-resolution X-ray imaging was useful to verify the presence of calcite precipitate inside the rock, the distribution of unconsolidated precipitate inside the rock was likely affected by the horizontal storage of the sample.

Fluid Mixing Test (Experiment BS0)
The experiment BS0 is a test to demonstrate the use of three-dimensional time-resolved neutron imaging in investigating fluid mixing phenomena within porous rock.We used two miscible fluids, a solution of deuterium oxide heavy water (D 2 O) and a solution of heavy water with dissolved cadmium (CdCl 2 ), an element that provides good contrast for neutron imaging and that was used in a previous study (Cordonnier et al., 2019).After fully flushing the sample with heavy water, we injected heavy water from one inlet at the top of the core sample and 1 M of CdCl 2 (aq) (dissolved in heavy water) from the other inlet at the top.Figure 2 shows the time series of the neutron tomograms during the D 2 O and CdCl 2 injections.For all experiments, chemicals were provided by Sigma Aldrich.
With this experiment, we investigate the use of neutron imaging to analyze fluid mixing by measuring transverse dispersion in sandstone rock.Assuming that the linear attenuation coefficients of injected fluids are constant through time, we can estimate the concentration of cadmium solution at each voxel of the three-dimensional volume and its relative evolution with time, based on the Beer-Lambert law (Swinehart, 1962).Taking the last scan of this test where the mixture reached a steady state, a CdCl 2 solution was injected continuously and formed a plume that extended downward via advection and grew in the lateral dimension through dispersion.Assuming the lateral velocity is small compared to the longitudinal velocity, the transverse dispersion equation for the steadystate plume is: where u z [L/T] is the average downward velocity, I[ ] is the normalized gray level intensity at each pixel (i.e., normalized neutron attenuation coefficient), and D T [L 2 /T] is the transverse dispersion, perpendicular to the interface between the two fluids.The solution of Equation 1 with unbounded boundary conditions is (Gautschi & Cahill, 1964): where I 1 and I 2 are the relative concentrations of the two solutions that correspond in our case to the normalized neutron attenuation coefficient of heavy water and cadmium solutions, x 0 [L] is the initial position of the interface between two solutions, and σ = ̅̅̅̅̅ αz √ , where solute dispersivity is defined as

Flow Prior to and After Calcium Carbonate Precipitation (Experiments BS1 and BS2)
The objectives of experiments BS1 and BS2 were to measure and assess the spatial distribution of velocity before and after two different calcite precipitation modes: the two solutions used to precipitate calcite were injected in the sample either one after each other from the same inlet in a sequential mode (experiment BS1), or simultaneously in two inlets in a co-flow mode (experiment BS2).We used several aqueous fluids.Distilled water (H 2 O) and heavy water (D 2 O) were mixed at different proportions in the preparation of the solutions to control the neutron imaging contrast.The steps of fluid injections in the rocks and the injection flow rates are indicated in Figure 3.To effectively detect and monitor D 2 O injection, the core sample initially had to be saturated with a mixture of D 2 O and H 2 O, ensuring an optimal signal-to-noise ratio.Following a series of calibration tests, similar to the procedure followed by Cordonnier et al. (2019), we chose a combination of 50% D 2 O and 50% H 2 O as a reasonable mixture, which is neither too transparent and nor too opaque to neutrons.The flow rates were calculated to ensure that at least one pore volume was injected in the sample over a duration sufficient to acquire the desired number of complete scans.The experiment BS1 (injection from the top of the core sample) contained two steps of mineral precipitation, and heavy water was injected before and after the second step (Figure 3a).For experiment BS2 (injection from the bottom), one step of mineral precipitation was performed with heavy water injections before and after this step (Figure 3b).For the heavy water injections, the characteristic advection time is approximated as  condition of 24°C.The combination of 75% heavy water and 25% distilled water was used for the CaCl 2 (aq) solution, and a combination of 25% heavy water and 75% distilled water was used for the Na 2 CO 3 solution.These proportions of distilled and heavy water provide good contrast in the neutron images to detect the miscible fluids separately.Because the precipitation is quasi-instantaneous when the two fluids mix into the rock, several types of calcium carbonate could form, such as amorphous calcium carbonate, vaterite, aragonite, and calcite (Brečević & Nielsen, 1989;Morse et al., 2007).Calcite is the most stable phase under our experimental conditions and therefore should be the dominating mineral at the end of each experiment.
The initially dry BS1 sample was filled with D 2 O at a constant flow rate, with water replacing air (Figure 3a, from 0 to 368 min).This imbibition step was followed by a step of injection of CaCl 2 (aq) and Na 2 CO 3 (aq) in parallel from the two inlets (co-flow from 368 to 462 min).However, the flow stopped due to the clogging of one of the inlet tubes with calcium carbonate precipitates.The experiment was continued as a sequential injection of CaCl 2 (from 462 to 940 min) and then Na 2 CO 3 (from 940 to 1463 min) solutions successively from the same inlet.After these steps, we injected D 2 O to measure the fluid velocity distribution.Then, two successive steps of sequential injections of CaCl 2 (from 1950 to 2,086 min) and Na 2 CO 3 (from 2,118 min to 2,309 min) induced a second episode of precipitation before another step of D 2 O injection was performed.It is important to note some experimental artifacts where we observed side flow (between the jacket and the sample rock as shown in Figure 1a) during injections, which was attributed to insufficient confining pressure in the system.Additionally, we faced failures of the camera of the neutron tomograph, which resulted in some incomplete scans indicated as yellow stars on the time axis in Figure 3.The duration of this experiment was 47 hr.
We performed a second experiment with calcium carbonate precipitation on sample BS2 (Figure 3b).Here, we injected the fluids from the bottom of the sample.Before the experiment, we first fully flushed the sample (∼two pore volumes) with a mixture of 50% of D 2 O and 50% of H 2 O.Then, we injected D 2 O (from 30 to 695 min) to measure the fluid velocity distribution.Following, we performed a co-flow injection of the two reactive solutions from the two inlets for an episode of calcium carbonate precipitation (from 775 to 1216 min), and finally injection of D 2 O (from 1246 min to 1955 min) to measure the velocity distribution after precipitation.In the course of the co-flow injection of CaCl 2 and Na 2 CO 3 , we experienced clogging of the tubes due to the formation of calcium carbonate precipitates in the porous media that redirected the flow of one solution toward the other.Therefore, this co-flow experiment was unsteady and effectively consisted of a series of injections of one and the other solution as the tube alternatively clogged and unclogged.The pre-mixing before the injection into the rock and also the side flows due to insufficient confining pressure contributed to a less homogeneous distribution of precipitation within the rock.This was interesting to test the flow imaging method in a heterogeneous sample.Two X-ray scans were acquired when the rock was saturated with D 2 O, one before and one after the precipitation phase.The second X-ray was subsequent to the final D 2 O injection, and this additional flow likely altered the unconsolidated precipitate distribution during the D 2 O injection process.The duration of this experiment was 34 hr.
In addition to the systematic errors associated with the flow condition in the experimental setup as mentioned above, other random error sources may exist, such as beam fluctuations and camera failure, which amplify the overall error magnitude.Precise quantification of the uncertainty pertaining to the measured velocity field, as presented in Section 2.3, is complicated.However, one of the main contributors to the error can be measured from random noise and uncertainty in the gray scale value of the pixels.To minimize this error, several postprocessing steps are taken, including scattering correction and spatial median filtering.The residual error is then quantified through an analysis of two successive scans conducted during the injection of heavy water prior to precipitation.During these scans, the flow in the core sample is in a stationary state.Comparison of two successive reconstructed images labeled 1 and 2, is performed using the root mean squared (RMS) difference as: where I 1(x,y,z) and I 2(x,y,z) are the normalized intensity of the voxels in the two reconstructed images of size M × N × L, when the sample rock is saturated with heavy water and flow is in a stationary state.The calculated RMS error in pixel intensity for sample BS1 is 4.6%, and for sample BS2 it is 4.2%.

Calculation of the Flow Velocity Field
The distinctive grayscale levels of the various miscible fluids utilized in our study enable us to identify and track the interface between them over time.The interface is seen as a strong gradient in the intensity of the neutron tomograms.Due to imposed fluid flow conditions, the interface is transported mainly in the direction of the cylindrical core axis.In order to compute the velocity of the interface, we solve an inverse problem based on the advection-dispersion equation by using the contrast in neutron attenuation between different fluids from timeresolved tomograms.The advection-dispersion equation is: where where u * is the velocity field transported by advection at each voxel for each pair of sequential images.We have combined this advection equation with a global smoothness term to constrain the estimated velocity field (Horn & Schunck, 1981) by minimizing the function: where Ω is the domain of the three-dimensional volume, and λ is a Lagrange multiplier to smooth out the sharp velocity changes.There is no rigorous methodology for determining the Lagrange multiplier in the variational formulation; however, the solution of this equation remains relatively insensitive to the value of λ within a broad range (Liu, 2017).We used a value λ = 50 empirically selected as suggested by Liu and Shen (2008).As shown in Figure S1, modifying λ in the range 10-50 does not have a significant impact on the velocity distribution, and when λ is larger, the variation will be smoothed out.We calculated the generalized form of the corresponding Euler-Lagrange equation proposed by Horn and Schunck (1981) to minimize J (u * ,v * ,w * ) and obtain the velocity field.The iterative numerical solution for the convergence is: where I x , I y , I z , and I t are intensity derivatives, and u * n , v * n , w * n are averages of the velocity at each voxel of 4 × 4 × 4 neighborhood voxels at iteration n.The estimation of the partial derivatives of the normalized intensity at each voxel is shown in the Supplementary Information S1 (Text S1).For the voxels on the boundaries of the image, the u* estimates are the averages of the neighboring velocity estimates.We use a maximum of 50 iterations to obtain the solution, and the iteration stops when the flow field has converged, that is, when at iteration n, where TOL is the tolerance level chosen here equal to 10 9 m/s.
The values of the cells in each u* matrix is obtained based on the intensity gradient which is coming from the interface of the injected fluid.After obtaining matrices of u* for each pair of sequential images, the final velocity field u = (u,v,w) transported by both advection and dispersion can be estimated by incorporating these matrices of u* into Equation 10.We use Equation 5 multiplied by a weighting function W = |∇I| as u * includes a contribution from the gradient of intensity.Therefore, the actual velocity field can be derived by solving the set of equations: where the dispersion coefficient D is approximated by calculating the divergence of Equation 5

Fluid Mixing During Co-Injection of Water and a Cadmium Solution (Experiment BS0)
This fluid mixing test is conducted to measure transverse dispersivity from neutron images.Figure 4a displays the last scan of the neutron tomogram series following the mixing of D 2 O and CdCl 2 shown in Figure 2. Slices at some selected heights of this three-dimensional volume are shown in Figure 4b.The top slices of the rock exhibit a distinct sharp interface between D 2 O and CdCl 2 solutions.As the fluids move toward the lower part of the sample, this sharp interface becomes more diffuse and widen owing to lateral dispersion.Figure 4c presents a series of normalized neutron attenuation profiles along the x-direction (perpendicular to the interface and the axis of the sample) at different heights in the sample.Each curve in Figure 4c is plotted along 45 mm in the x-direction where each point is calculated by averaging 20 × 20 voxels in the yz-plane and covers the central area (along the y-axis) of the slices at several heights z.By comparing the curves obtained from the upper (z = 7 mm) and lower (z = 42 mm) sections of the rock, we observe that the mixing zone has become wider toward the lower part of the sample (Figure 4c) since the change in attenuation with lateral distance has become more gradual.This widening is attributed to the transverse dispersion of solutes within the rock increasing with distance from the inlet.
The function described by Equation 2 can be used to fit the normalized attenuation curves obtained at various heights (Figure 4c).This fitting enables the determination of all the parameters in Equation 2 including the interface width σ at different heights z.The values obtained are shown in Figure 5.By fitting a square root dependence σ = ̅̅̅̅̅ αz √ to these data, we estimate the dispersivity α = 7.0 × 10 4 m.Considering the average downward velocity u z = 1.92 × 10 6 m/s (flow rate divided by the pore area of the surface), we estimate the transverse dispersion coefficient D T = αu z = 1.3 × 10 9 m 2 /s.This dispersion coefficient is of the same order of magnitude as the diffusion coefficient of the considered solutes, which is consistent with the considered Pe < 1 condition.The front progression of the injected D 2 O in the neutron tomograms allows the velocity distribution of the flow in samples BS1 and BS2 to be calculated using the inversion method based on the advection-dispersion equation (Equation 4and Section 2.3). Figure 8a shows the probability density function of the derived velocity fields of D 2 O flow in sample BS1.Note that these velocities correspond to coarse-grained velocity averaged over the neutron image pixel size and are therefore not measured at the pore scale.This result is attributed to the limited voxel size of the neutron images being of the same order of magnitude as the mean pore size of the rock.The fluid velocity magnitude at three horizontal cross-sections of the second D 2 O injection (before the second step of calcium carbonate precipitation) and the third D 2 O injection (after the second step of calcium carbonate  The probability density function of the velocity magnitude for D 2 O injection before precipitation (Figure 8a) shows a fluid velocity range of 0.14-0.26mm/min and a mean value of 0.19 mm/min.This value is in agreement with a velocity of 0.18 mm/min calculated from the injected flow rate of 0.073 ml/min and the available pore area of the surface (412 mm 2 surface area of the horizontal cross-section of the rock multiplied by the porosity equal to 21%).After precipitation, the mean value of fluid velocity is 0.20 mm/min.Comparing the second and third D 2 O injections shows that, for the same flow rate, the average fluid velocity after this step of precipitation is higher than before it, consistent with porosity reduction due to mineral precipitation.We also note that mineral   4) before and after the second step of mineral precipitation.

10.1029/2023WR036293
precipitation is associated with a widening of the velocity distribution, increasing the standard deviation from 0.018 mm/min to 0.023 mm/min.
For sample BS2, Figure 9 shows the fluid velocity magnitude of D 2 O injections before and after mineral precipitation and the probability density functions of fluid velocity.Before precipitation, the average fluid velocity is 0.14 mm/min, and after precipitation, it is 0.16 mm/min.The velocity magnitude images in Figure 9b show that the central part of sample BS2 hosts lower fluid velocities after precipitation in comparison to the sides of the core sample.The lower velocity in the center, compared to the sides, could be attributed to insufficient confining pressure in the system during the experiment which resulted in some flow along the side of the sample.We observe that the enhancement in flow heterogeneity is more visible in BS2 than in BS1.The probability density functions of fluid velocity (Figure 9a) clearly show the influence of mineral precipitation between the two injections of D 2 O. Similar to the results obtained for sample BS1, fluid flow in sample BS2 shows a higher average velocity and a wider range of fluid velocities after a step of calcium carbonate precipitation.In this case, the standard deviation increased from 0.019 mm/min to 0.026 mm/min.In the calculation of the velocity field, we assume that the estimated dispersion coefficient derived from Equation 11 is constant during each injection.While this assumption appears somewhat restrictive, we verified that the estimated values of D are nevertheless too small to significantly influence the velocity distributions.The derived values of dispersion coefficient D are equal to 5.6 × 10 11 and 9.3 × 10 11 m 2 /s for the second and third D 2 O injections into BS1.For the first and second D 2 O injections into BS2, the values are 8.0 × 10 11 and 3.4 × 10 10 m 2 /s, respectively.
One of the primary sources of error arises from random noise and the uncertainty associated with the grayscale values of voxels.To quantify how this error propagates into the velocity measurements, we added random noise with Gaussian distribution to each voxel of the scans during heavy water injection before precipitation to achieve the same errors we measured in Section 2.2.2 (4.6% for the scans of sample BS1 and 4.2% for sample BS2).Subsequently, flow velocity fields were computed based on the images with added noise.The correlation between velocities from the noisy and original images is shown in Figure S2.The RMS error between the original velocity and the velocity derived from images with added noise is 0.018 mm/min for sample BS1 and 0.013 mm/min for sample BS2.Dividing these errors by the average fluid velocity during heavy water injection before precipitation (BS1: 0.19 mm/min and BS2: 0.14 mm/min) gives the estimated propagated error percentages in the velocity field which are estimated at 9.5% for BS1 and 9.2% for BS2.Notably, this approach may potentially overestimate errors, considering the additional noise introduced and the non-linear nature of propagated errors in the velocity field.Figures 10a and 10b display density maps illustrating how local fluid velocities after a precipitation step, compare to local velocities before the precipitation step for samples BS1 and BS2, respectively.In the case of BS1, where the precipitation step occurred through sequential injections, there is no correlation between these velocities, suggesting more random changes in the velocity field in response to the precipitation event.For sample BS2, where precipitation was induced by unsteady co-flow injections, we observe a slight positive correlation, where higher initial velocities tend to further increase following precipitation.This observation is substantiated by measuring the average velocity after precipitation as a function of the initial velocities, as shown by the red line along with the variance in the same figure.In the case of BS1, we observe a consistent average velocity after precipitation.However, for sample BS2, there exists an upward trend, signifying a positive correlation.This highlights the distinct behavior exhibited by these two samples in response to the precipitation event.

X-Ray Microtomography Imaging
For experiment BS2, we acquired two X-ray tomography scans at the Paul Scherrer Institute, one before and one after mineral precipitation.Figure 11a shows the difference between these two scans, displaying the areas where  calcium carbonate precipitate appeared.Note, however, that some unconsolidated precipitate likely moved during the heavy water injection.Instead of a straight plane of precipitation in the middle of the rock, which would be expected for a steady co-flow mixing experiment, we observe a complex region of precipitation likely due to the unsteady injection of CaCl 2 and Na 2 CO 3 .The amount of precipitation is higher at the bottom near the inlet points where the solutions of CaCl 2 and Na 2 CO 3 were injected and also at the top where the calcium carbonate precipitates accumulated at the outlet.The spatial resolution of these scans was too limited to observe the details of the calcite crystals, and only some clusters calcite accumulated near the boundaries are noticeably visible.Noting that the signal-to-noise ratio of the image falls within an acceptable range, some precipitation occurred in the middle of the sample.Figure 11b shows the average of the difference of gray values after and before mineral precipitation (X-ray attenuation coefficients) along the height of the sample BS2.The curve demonstrates a sharp variation in density near both the inlet and outlet points.We also see increased X-ray adsorption due to the precipitation of calcium carbonate in the center part of the sample, albeit slightly less pronounced than near the inlet/outlet regions.
Segmentation of the pore space in the high-resolution three-dimensional tomogram acquired at the European Synchrotron Radiation Facility shows the porosity after precipitation (Figure 12a).Figure 12b displays four selected slices perpendicular to the axis of the core sample, at different heights.The pores appear in dark gray because of low X-ray adsorption, while the rock grains and precipitates have lighter gray colors.Figure 12c shows the profile of porosity along the height of the core derived from pore segmentation, which we calculate as the summation of segmented pores' surface area divided by the slice surface area along the axis of the core sample.Porosity at both the top and bottom regions of the rock sample is lower than the average value.Figure 12d plots the normalized X-ray attenuation profile by averaging gray-scales along the core sample's axis.These observations from Figures 12c and 12d indicate that a higher level of precipitation occurred at the top and bottom sections of the rock sample, a result that agrees with the observations shown in Figure 11 for the X-ray data acquired at the Paul Scherrer Institute.However, it should be mentioned that these results were acquired using different camera specifications with different spatial resolution, and there was a considerable time gap of several months between the acquisitions that may have impacted the density contrast and attenuation range observed in the images.Therefore, we focus on the qualitative aspect of the similarity between Figures 11b and 12d, which are indicative of calcite precipitation distribution over the sample rock.

Neutron Imaging of Fluid Flow in Porous Media
Recent advances in neutron imaging technique offer the capability to acquire three-dimensional time series of fluid flow in opaque rocks with acquisition times in the range 5-30 min and spatial resolution of several tens of micrometers (Kaestner et al., 2007;Tengattini et al., 2021;Tötzke et al., 2011;Tudisco et al., 2019).The distinctive characteristics of a neutron beam, such as its isotope sensitivity (e.g., D 2 O/H 2 O) or its strong attenuation in the presence of cadmium, make neutron imaging a useful technique to study fluid flow and contaminant transport in porous rock.Heavy water (D 2 O) exhibits significantly reduced neutron attenuation compared to regular water (H 2 O), allowing for the investigation of fluid flow dynamics.Recent advances in neutron tomography techniques have facilitated quicker imaging and enabled scientists to capture fluid dynamics with enhanced spatial and temporal resolution.As an example, the CONRAD-2 neutron source (Helmholtz-Zentrum, Berlin) enables the monitoring of fluid front movement in saturated samples by utilizing the distinctive contrast in neutron attenuation (Etxegarai et al., 2021).The Institut Laue-Langevin in Grenoble, France is a neutron source where high-speed neutron tomographies with a rapid acquisition time of around 1 minute have been obtained to study water invasion into air-filled samples of rock (Tudisco et al., 2019).The neutron source SINQ at the Paul Scherrer Institute in Switzerland is another facility where we show in the present study the capability of in situ analysis of the local changes in fluid flow properties of rocks, expanding the possibilities for studying fluid dynamics in rock samples (Cordonnier et al., 2019).In the present study, we quantify the spatial fluctuations of fluid velocity before and after calcium carbonate precipitation in porous sandstone.We solved the inverse advection-dispersion equation using the difference in neutron attenuation between two miscible fluids, a technique that to the best of our knowledge has not been previously explored.The resolution of neutron images implies that the fluid velocities correspond to average pore velocities, which do not resolve the velocity field inside pores.

Dispersion of Solutes Quantified by Neutron Imaging
Recent advancements in understanding transverse dispersion dynamics in natural rock media have been made using X-ray tomography to study the steady-state dispersion of solutes (Boon et al., 2016(Boon et al., , 2017)).The potential of utilizing neutron imaging to measure dispersion is evaluated in this study.In experiment BS0, the width of the mixing zone between the heavy water and cadmium solutions increases from the top to the bottom of rock sample BS0, as expected, showing the applicability of the neutron imaging technique in exploring the transverse dispersion in rock.The derived value of the transverse dispersion coefficient was obtained by fitting Equation 2 to the normalized attenuation curves (Figure 4c) collected at different heights in sample BS0 is D T = 1.3 × 10 9 m 2 /s.Assuming the medium is homogeneous and the flow rate is low (Pe = 0.5), the transverse dispersion would be close to the molecular diffusion.The value of the dispersion coefficient we calculate is within a factor two comparable to the diffusion coefficient of cadmium into water (7 × 10 10 m 2 /s) published in other studies (Furukawa et al., 2007;Lide & Kehiaian, 1994).

Fluctuations of Fluid Velocity in Sandstone Before and After Mineral Precipitation
The injection steps of D 2 O in experiments BS1 and BS2 produce a visible fluid front in the tomograms.The shape and velocity of this front change before and after mineral precipitation steps, indicating local modifications of porosity that we interpret to be related to the clogging of some pores by calcium carbonate precipitation.In both Water Resources Research 10.1029/2023WR036293 samples BS1 and BS2, the velocity distribution curves observed after a precipitation step in Figures 8a and 9a exhibit a higher and broader range of velocities.This effect is more pronounced in experiment BS2, where a higher amount of calcium carbonate precipitation occurred following the co-injection of Na 2 CO 3 and CaCl 2 .
Our results suggest that non-uniform precipitates can readily form at the pore scale.The decrease in permeability is associated with a reduction in porosity (Reis & Am, 1994), due to the decrease in the size of the pore throats.This decrease in permeability induces an increase in the average fluid velocity under conditions of imposed flow rate.The formation of new minerals may also result in the development of roughness in the pores (Noiriel et al., 2016), which, in turn, may lead to a reduction in fluid velocity near the grain surface.This reduction is controlled by both the location and size of the newly formed crystals.The impact of local precipitation on fluid velocity may have substantial implications for the hydrodynamics and transport of reactants and products in the vicinity of the fluid-rock interface (Steefel & Maher, 2009).
For both experiments BS1 and BS2, we observed fluid velocity increases in certain areas and decreases in other areas of the rock sample.This effect is related to the porosity reduction due to localized calcium carbonate precipitation.To quantify the velocity field, we first calculated the path-averaged velocity of D 2 O transported by advection, and then, we added a dispersion term with the approximation described in Equation 11.The derived values of the dispersion coefficients in samples BS1 and BS2 before precipitation.5.6 × 10 11 m 2 /s and 8.0 × 10 11 m 2 /s (see Section 3.2), are smaller than the value of 1.3 × 10 9 m 2 /s calculated for sample BS0.One possible explanation is the use of a smoothing Lagrange multiplier in the variational formulation (Equation 7), which acts like a diffusion term (Liu & Shen, 2008).However, we have shown that varying the value of the Lagrange multiplier between 10 and 50 does not modify significantly the distribution of the velocities (Figure S1).Another explanation could be attributed to a bias in the three-dimensional volume reconstruction process, which may have arisen from the relatively long 32-min temporal resolution between each scan, during which fluid flow was occurring within the rock.These factors combined with the assumption of considering the constant dispersion coefficient throughout the injections have led to the underestimation of this coefficient in samples BS1 and BS2.

X-Ray Imaging of Local Calcium Carbonate Precipitation and Porosity and Comparison With Neutron Imaging
X-ray imaging has been used to investigate the precipitation of calcium carbonate and porosity distribution in rocks.One study has shown that X-ray micro-tomography can provide a detailed view of the morphology and distribution of calcite precipitates in fractured rocks (Aben et al., 2017).X-ray imaging has also been used to observe how calcite precipitation affects the permeability, porosity, and fluid flow in porous rocks (e.g., Hebert et al., 2015;Noiriel, 2015;Noiriel et al., 2016).In this study, X-ray tomography acquisition of the sample BS2 allowed us to verify that calcite precipitation occurred in the rock and have a general estimation of its distribution.However, these X-rays cannot be used to explain the flow distribution since the unconsolidated precipitate likely moved in the pore space between the neutron and X-ray scans.Figure 13 presents vertical cross-sections of the sandstone samples acquired using different imaging methods.Figure 13a shows a section obtained using neutron imaging, which, despite its lower resolution, provides valuable insights to image the solutions present within the samples.Figures 13b and 13c display sections acquired using low and high-resolution X-ray imaging.Notably, the high-resolution X-ray imaging provides enhanced visualization of the pore space and calcite distribution within the sample.At the pore scale, calcite-filled pores can be identified.As discussed above, it is not yet possible to fully connect the three-dimensional X-ray image with the three-dimensional velocity field obtained from the neutron tomogram series in order to determine the relationship between mineral precipitation and the fluid velocity at that region.Developing imaging techniques to allow for such developments should be the topic of future studies.

Conclusions
This study uses 4D neutron imaging to investigate flow and transport before and after mineral precipitation in core samples of Bentheim sandstone.Neutron imaging is shown to be a useful tool to track the transport of solute fronts in porous rocks since the rock and sample holder are mostly transparent to neutrons, whereas the aqueous fluid can be imaged under in situ conditions.
Hence, we investigated the feasibility and precision of neutron imaging to measure transverse solute dispersion induced by flow in porous rock.Furthermore, neutron imaging was used to visualize the flow of two aqueous miscible fluids within porous rock and its alteration by calcium carbonate precipitation.We developed an inverse method that uses the advection-dispersion equation to quantify the spatial fluctuations of the velocity of injected fluid by tracking the moving interface between two miscible fluids.While this technique did not provide porescale resolution, it was useful to monitor the distribution of the average pore velocity within the rock.After precipitation, the average velocity was larger and the flow field was more heterogeneous.This effect was more pronounced in the co-flow experiment compared to the sequential experiment.In the former, larger initial velocity tended to increase, explaining the broadening of the velocity distribution.X-ray imaging at different resolutions allowed us to confirm the presence of calcite precipitate inside the rock and estimate its spatial distribution in the rock.However, the transport of unconsolidated precipitate between the neutron and X-ray scan did not allow a direct correlation between the flow and precipitation fields.This could be achieved in the future by optimizing the coupling of neutron and X-ray imaging sequences.
These findings therefore open new opportunities to use neutron imaging, coupled with other tomography techniques, to investigate flow, transport, and reaction processes in porous media.Our results suggest an interesting coupling between flow and precipitation, particularly in the co-flow experiment.Following this proof of concept, the technique could be largely improved by better controlling the injection of reactive fluids, investigating reactive transport regimes where the precipitate is not transported, and optimizing the coupling of neutron and Xray imaging.
the Paul Scherrer Institute have a voxel size of 32 μm with an acquisition time of 75 min.

Figure 1 .
Figure 1.Schematic of the experimental set-up at the ICON beamline at the Paul Scherrer Institute, Switzerland: (a) Sample rock (5 cm diameter, 5 cm height) inside the core holder with two inlets for fluid co-flow injection and one fluid outlet.Here, the sketch shows injection from the bottom of the core sample, against gravity.Rotating this set-up upside down allows injecting fluids from the top, in the direction of gravity.(b) Sketch of the rotating core holder on the rotating stage and image acquisition using both neutron and X-ray sources.
100s, where L = 200 μm is the pore size, v = 2 × 10 6 m/s is the fluid velocity.The characteristic diffusion time is approximated as τ d = L 2 /D = 50s where D = 7.5 × 10 10 m 2 /s.Thus, the Péclet number is Pe = τ d /τ a = 0.5.While the characteristic reaction time varies locally within the porous media and changes over time, it is approximated τ r = 5s byIzumoto et al. (2022), obtained from a millifluidic experiment of calcite formation leading to a Damköhler number Da = τ a /τ r = 20.The precipitation of calcium carbonate is controlled by the constant injection of two aqueous solutions containing CaCl 2 and Na 2 CO 3 at the inlets (CaCl 2 (aq) + Na 2 CO 3 (aq) → CaCO 3 (s) + 2NaCl (aq) ) under ambient temperature

Figure 2 .
Figure 2. Time series of neutron tomograms acquired during the co-injection and mixing of D 2 O (transparent) and CdCl 2 (orange color) solutions injected simultaneously from the two inlets at the top of sample BS0 (gray color).Initially, the rock was saturated with D 2 O.The time t = 0 corresponds to the start of the injection.

Figure 3 .
Figure 3.The fluid injection steps of the two mineral precipitation experiments performed on Bentheim sandstone.For experiment BS1 (a), the fluids were injected from the top to the bottom, along the direction of gravity, and for experiment BS2 (b) the injection was from the bottom to the top.Calcium carbonate precipitation was performed by either injecting simultaneously the Na 2 CO 3 and CaCl 2 solutions in the two inlets of the core holder (co-flow, pink) or by injecting one solution after the other from the same inlet (sequential flow, orange).Short periods where water was injected to flush the tubes and avoid calcium carbonate precipitation in them are also indicated in bright blue.The core samples were imaged with three-dimensional neutron tomography during the main injection steps, and two-dimensional radiographs were acquired during the flushing steps.The red stars on the time axis show the complete neutron tomography scans and the yellow stars indicate incomplete scans (due to camera failure), from which we used two-dimensional radiographs.Due to their longer acquisition time, two X-ray scans were acquired before and after calcium carbonate precipitation for sample BS2.
φ [ ] is the porosity, I [ ] is the normalized gray level intensity at each pixel, t [T] is the time between two consecutive tomograms, u = (u,v,w) [L/T] is the fluid velocity, and D [L 2 /T] is a dispersion coefficient.Defining the velocity vector u * = (u * ,v * ,w * ) as u * = u D ∇I I (5) allows deriving the three-dimensional velocity field of the injected heavy water by tracking how the interface between heavy water and distilled water moves through the sample.The path-averaged velocity of fluid transported by advection only is calculated by: φ ∂I ∂t + u * • ∇I = 0 and considering mass conservation (∇ • u = 0) as

Figure 4 .
Figure 4. (a) Neutron tomogram with a vertical cross-section of the core sample in experiment BS0 after co-flow injection of D 2 O and CdCl 2 .The orange color shows the cadmium solution injected from the right inlet located at the top of the sample while injecting heavy water from the left inlet.(b) Two-dimensional slices perpendicular to the core axis at different heights.(c) Normalized attenuation along transects perpendicular to the interface of D 2 O and CdCl 2 mixing at different heights of sample rock BS0.These lines are calculated as an average of 20 pixels wide transects located at the position of the red arrow in panel (b) and crossing the interface between the two fluids perpendicularly.The broadening of the intensity profiles when crossing the interface between the two fluids from the top (z = 7 mm) to the bottom (z = 42 mm) of the sample is due to flow dispersion.
Distribution Before and After Calcium Carbonate Precipitation (Experiments BS1 and BS2) The difference in neutron absorption and scattering cross-sections between D 2 O and the initial fluid inside the samples (mixture of D 2 O, H 2 O, and different solutes) is used to image the flow of D 2 O. Figures6 and 7show the time series of the neutron tomograms during the injection of D 2 O into samples BS1 and BS2 before and after calcite precipitation.The boundary separating the two aqueous solutions is visible but is not a sharp front due to the dispersion occurring within the interface of the fluids.It should also be noted that the exposure time needed to image the flow interface during injection may blur the moving front in the three-dimensional images; however, injecting with a low flow rate helped reduce this effect.Visual comparison between the first row (before precipitation) and the second row (after the precipitation) within Figures6 and 7indicates dissimilarities in the shape and the velocity of the front progression.This result suggests that calcium carbonate precipitation modified flow pathways in both samples.

Figure 5 .
Figure 5. Derived values of interface width σ at different heights z in experiment BS0 (filled squares).The solid line shows the best fit of the function σ = ̅̅̅̅̅ αz √ to the data, which allows us to estimate the dispersivity α.Injection of the two fluids occurs at the height z = 0.

Figure 6 .
Figure 6.Neutron tomograms during injection of D 2 O (dark color) from the top with flow rate 0.073 ml/min into sample BS1.Prior to each injection, the rock was saturated with the calcium carbonate growth solution (Na 2 CO 3 and CaCl 2 ), which contained a 25%-75% mixture of D 2 O and H 2 O (green color).Times t correspond to those indicated in Figure 3a, and Δt is the time since the beginning of each injection step.

Figure 7 .
Figure 7. Neutron tomograms during injection of D 2 O (dark color) from the bottom with flow rate 0.073 ml/min into sample BS2.Before injection, the rock was saturated with a 50%-50% mixture of D 2 O and H 2 O (blue color).Times correspond to those indicated in Figure 3b, and Δt is the time since the beginning of each injection step.

Figure 8 .
Figure 8.(a) Probability density functions of fluid velocity magnitude during two injections of D 2 O (2nd and 3rd D 2 O injections) in sample BS1 with an error confined within 9.5% margin.The second injection (blue curve) was performed before the second step of calcium carbonate precipitation (Figure 3a) and the third one (green curve) after it.(b) Cross-sections perpendicular to the core axis at different heights of local fluid velocities (Equation4) before and after the second step of mineral precipitation.

Figure 9 .
Figure 9. (a) Probability density functions of fluid velocity magnitude during the injections of D 2 O in sample BS2, before (blue curve) and after (green curve) calcium carbonate precipitation.The error is confined within 9.2% margin.(b) Cross-sections perpendicular to the core axis at different heights of local fluid velocities before and after the mineral precipitation.

Figure 10 .
Figure 10.Density maps depict the local velocity after a precipitation step in relation to the velocity before the precipitation step for samples BS1 (a) and BS2 (b).The color bar indicates the density of the scatter data.The red lines indicate the average post-precipitation velocity at a given initial velocity, along with the variance.

Figure 11 .
Figure 11.(a) X-ray tomography three-dimensional rendering of precipitation of calcium carbonate in sample BS2 acquired at the Paul Scherrer Institute, with a voxel size of 32 μm.The image was produced by subtracting the X-ray tomogram after precipitation from the X-ray tomogram before precipitation.(b) Effect of calcium carbonate precipitation: The differential Xray attenuation profile along the axis of the core sample is calculated by subtracting the three-dimensional volume after precipitation from the three-dimensional volume before precipitation and averaging along the axis of the core sample.

Figure 12 .
Figure 12.(a) Rendering of the X-ray tomogram of sample BS2 acquired at the European Synchrotron Radiation Facility, with a voxel size of 5.85 μm.Pores are shown in blue.(b) Two-dimensional slices of porosity distribution perpendicular to the core axis at different heights including X-ray tomograph zooms of small squares, where empty pores appear in black.(c) Profile of porosity distribution along the axis of the core sample calculated by segmentation of the pores (sum of segmented pores divided by slice area), as shown in panel (a).(d) Normalized X-ray attenuation profile calculated by averaging gray values along the axis of the core sample.