Creep Burst Coincident With Faulting in Marble Observed in 4‐D Synchrotron X‐Ray Imaging Triaxial Compression Experiments

Faults in carbonate rocks show both seismic and aseismic deformation processes, leading to a wide range of slip velocities. We deformed two centimeter‐scale cores of Carrara marble at 25°C and imaged the nucleation and growth of faults using dynamic synchrotron X‐ray microtomography. The first sample experienced a constant confinement of 30 MPa and no pore fluid. The second sample experienced confinement in the range 35–23 MPa and water as a pore fluid at 10 MPa pore pressure. We increased the axial stress by steps until creep deformation occurred and imaged deformation in 4‐D. The samples deformed with a quasi‐constant or increasing strain rate when the differential stress was constant, a process called creep. However, for both samples, we also observed transient events that include the acceleration of creep, that is, creep bursts, phenomena similar to slow slip events that occur in continental active faults. During these transient creep events, strain rates increase and correlate in time with strain localization and the slow development of system‐spanning fault networks. In both samples, the acceleration of opening and shearing of microfractures accommodated creep bursts. High‐resolution time‐lapse X‐ray microtomography imaging and digital image correlation during triaxial deformation quantify creep in laboratory faults at subgrain spatial resolution. This work demonstrates that transient creep events, that is, creep bursts or slow slip events, correlate with the nucleation and slow growth of faults and not only with slip on preexisting faults.


Introduction
Fault slip can reach velocities on the order of a meter per second during earthquakes or much slower velocities when displacements occur over hours to weeks. Events of such slower velocities, known as creep transients, or slow slip events when applied to active faults, indicate a permanent deformation whose rate ranges between the tectonic loading rate and earthquake rates (e.g., Bürgmann, 2018). At least two physical mechanisms have been linked to creep in the Earth's upper crust. Brittle creep can occur by the chemically activated slow growth and coalescence of microfractures (Brantut et al., 2013;Scholz, 1968). Alternatively, pressure solution creep arises from the coupling between mechanical and chemical forces at the grain scale such that mass transfer through dissolution and precipitation processes controls volumetric deformation (e.g., Gratier et al., 2013;Rutter, 1976). The rheological laws of these two mechanisms describe strain rate as a function of a series of mechanical, chemical, and petrophysical parameters, such as stress, temperature, rock composition, grain size, and fluid composition.
When the system is subjected to less than 90% of failure stress, the creep rate is usually quasi-constant through time (e.g., Lockner, 1993). When approaching failure, the coalescence of microfractures may lead to an exponential or power law increase of creep rate, until catastrophic failure (Amitrano & Helmstetter, 2006;Main, 2000;Reches & Lockner, 1994). Exponential or power law increase of porosity is also observed in experiments where the loading stress is not constant but increases by steps when approaching failure (Cartwright-Taylor et al., 2020;Renard et al., 2018). In creeping faults, the creep rate may be constant through time, corresponding to steady-state slip. However, transient stages of the acceleration of slow displacements in continental faults have been observed in borehole strain meters (Linde et al., 1996) and in geodetic data (Crescentini et al., 1999;Jolivet et al., 2013;Rousset et al., 2016). Such creep bursts, or slow slip events, are observed as periods of the increase of slip rate before the fault either becomes locked again or continues to creep at a lower rate. Creep rates measured in major continental faults using time-lapse satellite interferometry indicate that slow slip events can have a wide range of slip surface areas and durations (Jolivet et al., 2015). Some slow slip events can occur on a fault tens of years after a major earthquake (Aslan et al., 2019). Slow slip events can trigger seismicity (Lohman & McGuire, 2007) and may control the nucleation process before some major earthquakes (Bouchon et al., 2011). Creep transients induced by anthropogenic fluid injections can also trigger seismicity at the meter to kilometer scales (Guglielmi et al., 2015;Wei et al., 2015).
Because creep deformation can precede seismic failure in rocks, an acceleration of creep may indicate an approaching catastrophic event (e.g., Kranz, 1980). This recognition led to the concept of predicting the time to failure of major earthquakes, landslides, and volcanic eruptions (Main, 1999;Voight, 1989). The creep evolution of some landslides before failure follows this concept (Carlà et al., 2019). However, transient slip acceleration does not always indicate the onset of a catastrophic failure, as for the Mud Creek landslide in California, for example (Handwerger et al., 2019).
Observations of rock creep in the crust have motivated the development of laboratory experiments to measure the process and propose rheological laws. Since early creep experiments on sedimentary rocks (Griggs, 1939) and granodiorite and gabbro rocks (Lomnitz, 1956), series of laboratory experiment studies have characterized brittle creep in various rocks, such as granite (Kie et al., 1989;Kranz & Scholz, 1977;Lei et al., 2000;Lockner, 1993;Lockner & Byerlee, 1977;Ross et al., 1983), basalts (Heap et al., 2011), amphibolite (Satoh et al., 1996), marble (Fredrich et al., 1989;Liu & Shao, 2017;Quintanilla-Terminel & Evans, 2016;Tal et al., 2016;Yang et al., 2015), and sandstone (Baud & Meredith, 1997;Ngwenya et al., 2001;Shengqi & Jiang, 2010;Tsai et al., 2008). In all of the experiments that reached failure, an acceleration of creep strain rate occurred before failure, and for those for which accoustic emission recording was performed, an increase of acoustic emissions was observed. These results suggest the predictability of the time to failure. However, experimental techniques that use acoustic emission recording are blind to aseismic deformation mechanisms at the grain scale because they can only detect the seismic component of deformation. This limitation challenges attempts to estimate the time to failure from microstructural parameters, such as the evolving fracture network geometry, and the validation of theoretical studies on creep in rocks.
In creep experiments under constant stress conditions, the macroscopic axial strain rate may increase or decrease, often with three main stages until macroscopic failure, the so-called primary, secondary, and tertiary phases (e.g., Brantut et al., 2013;Lockner, 1993). First, primary creep occurs with an initial nonlinear increase in strain with time. Next, the onset of secondary creep occurs as a decrease in strain rate, producing a quasi-linear relationship between strain and time. Finally, tertiary creep occurs as an acceleration of strain with an exponential or power law dependence that ends with catastrophic failure. More recent literature on brittle creep has moved away from describing the creep curve in terms of three phases (e.g., Figure 6 and associated discussion in Brantut et al., 2014a). The curve of strain versus time may be divided into two stages: a deceleration and an acceleration of strain with time. Between these two regimes, the creep rate is nearly constant and the strain increases quasi-linearly with time ( Figure 1). In the present study, we show in laboratory experiments that the stage of accelerating creep may contain transient accelerations of strain, that is, creep bursts, superimposed on the exponential or power law trend of strain and time. Strain transients were previously observed in creep experiments in Darley Dale sandstone (Figure 1 in Heap et al., 2009) and basalt from Mount Etna (Figure 1 in Heap et al., 2011). Here we demonstrate that strain transients in Carrara marble correspond to the slow nucleation and growth of faults.
A key question in studying creep processes is how mechanisms at the grain scale, such as microfracturing, interact at mesoscopic and then macroscopic scales to produce the observed permanent creep or transient creep bursts observed in the upper crust (e.g., Brantut et al., 2013). The rock deformation experiments in the creep regime with dynamic X-ray microcomputed tomography imaging (4-D XCT), described in the present study, enable linking the macroscale behavior and the microscale processes under in situ conditions of the upper crust. This technique images the interior of a rock sample, while it deforms, with spatial sampling below the grain scale. Here, we show that the nucleation and growth of microfractures accommodate creep under constant stress conditions. A creep burst detected during each experiment coincides in time with the nucleation and propagation of system-spanning fault networks. We track the microstructural geometric properties of the microfractures and the local incremental strain deformation with digital volume correlation. The results characterize the evolution of porosity detected in the X-ray images, dilation, compaction, and shear strain before, during, and after the creep burst event.

Dynamic Synchrotron In Situ Experiments
The rock samples used here are cylindrical cores, 5 mm in diameter and 10 mm in height, drilled from a block of Carrara marble. They come from the same block used in the experiments of rock failure performed in Kandula et al. (2019). Carrara marble is a coarse-grained, nearly pure calcite rock with grain size in the range of 100-200 μm and initial (undeformed) porosity less than 1%. Each core sample was inserted into a Viton jacket and between two stainless steel pistons. The interfaces between the samples and the piston were not lubricated. This sample assembly was mounted into the Hades triaxial rock deformation apparatus. The details of the apparatus, including sketch and operating conditions are described in Renard et al. (2016). This apparatus is installed on the beamline ID19 at the European Synchrotron Radiation Facility (ESRF) and is used to perform 4-D XCT imaging during rock physics experiments at in situ conditions (Renard et al., 2017). Table 1 describes the experimental conditions for the two samples. Figure 2 shows the axial stress, confining pressure, differential stress, and axial strain as a function of time. Experiments were performed at the room temperature of the hutch of beamline ID19 at ESRF, in the range 23-25°C.
The 4-D XCT data were acquired by rotating the Hades rig, with the sample inside, over 180°and taking either 1,800 (Sample M83) or 1,600 (Sample M84) radiographs using the full white beam of the synchrotron. In these experiments, the Hades rig acts as a filter for X-rays so that an equivalent energy of 85 keV crosses the sample. For Experiment M83, each scan duration lasted between 1.5 and 2 min with an average duration of 1.7 min. When the axial stress was increased between two sequential scans, the scans were acquired 2 min apart. When the axial stress was not modified between two successive scans, each scan was acquired immediately after the other without a rest period between scans. For Experiment M84, the complete series of scans Note. P c : confining pressure; P p : pore fluid pressure; P O-ring : differential pressure due to friction in the O-rings of the rig; T: temperature; σ CB : differential stress during the creep burst; ε y : macroscopic axial strain at yield; nb. XCT: number of 3-D microtomography tomograms acquired during the experiment; duration: total duration of the experiment from the onset of loading to unloading.

Figure 1.
Sketch of the evolution of strain until sample failure during a laboratory creep test, with two stages: a first stage where strain rate decelerates with time (creep rate decreases) and a second stage where strain rate accelerates (creep rate increases) until failure. At the transition between the two stages, the creep rate is quasi-constant and reaches a minimum value. The inset above the main curve shows a transient acceleration of strain, also called a creep burst. The inset below the main curve shows the strain rate as a function of time, the separation between the two regimes, and a creep burst. The present study demonstrates that creep bursts may be coincident with the nucleation and growth of faults in laboratory experiments on Carrara marble. had 2 min scan durations and 2 min rest periods between scans. Tomographic reconstruction is performed using the program PyHST2 (Mirone et al., 2014). The voxel size is 6.5 μm 3 . The tomograms map the three-dimensional X-ray attenuation in the samples, with the calcite grains having a stronger attenuation (light gray values) than air-filled voids (dark gray values). Each tomogram of the rock sample contained around 7.1·10 8 voxels.
Due to phase contrast, one can detect features that are smaller than the voxel size. For example, crack interfaces produce a phase contrast. However, this identification in the phase contrast data does not indicate that one can resolve the crack volume. On beamline ID17 at ESRF, and using the same detector and optics as on ID19, a study has quantified a spatial resolution of the order of 20 μm (see Table 5 in Mittone et al., 2017). However, these authors did not use the phase contrast imaging technique we used on ID19, and they used an energy (50 keV) lower than the energy used here (85 keV). In summary, the property of phase contrast allows separating objects with dimensions of the order of the voxel size. The true spatial resolution is of the order of 20 μm, and the sampling distance (i.e., voxel size) of our images is 6.5 μm.

Stress Control and Axial Strain Measurement
For both samples, the axial stress, σ 1 , confining pressure, P c , and pore pressure, P p , are controlled by independent pumps connected to the Hades rig. The differential stress in the sample is  For Experiment M83, around 20 tomograms were acquired during fault network nucleation and formation (i.e., creep burst and slower continuous creep before and afterward). For Experiment M84, only three tomograms were acquired during macroscopic fault formation (i.e., creep burst) and four more scans were acquired before this sample failed. The largest principal stress, σ 1 , is parallel to the axis of the cylindrical sample.  Figure S1. (Table 1). For the first sample, M83 (dry), the confining pressure was constant and equal to 30 MPa, with no pore fluid pressure. The axial stress was increased by steps of 5 MPa below a differential stress of 100 MPa, then by steps of 2 MPa between differential stresses of 100 and 120 MPa, then steps of 1 MPa until a differential stress of 154 MPa was reached. The average equivalent stress rate corresponds to the slope of the red curves (differential stress) in Figure 2. The average stress rate was 1 MPa/min from 0 to 100 MPa differential stress and then 0.8 MPa/min between 90 and 154 MPa differential stress.
The sample was left under a constant differential stress of 154 MPa and then started to creep (left gray rectangle in Figure 2a). A first series of 3-D radiographs was then acquired (Scans 90-130, Figure 3a). Then, the differential stress was increased to 159 MPa, cycled four times between 158 and 159 MPa to enhance creep, and left constant at 159 MPa for~30 min (Scans 139-162, right gray rectangle in Figure 2a), during which a fault network propagated across the sample ( Figure 3b). The slow propagation of the fault network corresponds in time with an acceleration of creep strain rate between two periods of slower creep under constant differential stress, and we define this acceleration as a creep burst. It is observed macroscopically (Figures 3b and 3c) and occurred 4 min after the differential stress was increased from 158 to 159 MPa. Ten tomograms were acquired during this event, between Scan Numbers 141 and 150, then the samples continued to creep until Scan 162 where the experiment was stopped before a final failure could have been reached (Figures 2b and 3b).
For the second sample, M84 (wet), the initial confining pressure was set to 35 MPa during most of the experiment and then reduced to 23 MPa at the end of the experiment. The pore fluid pressure was constant and equal to 10 MPa. The pore fluid was deionized water equilibrated with calcium carbonate powder for 24 hr and then filtered. The sample was saturated in water under vacuum for 24 hr before the experiment. The axial stress was increased by steps of 5 MPa until a differential stress of 93 MPa, then by steps of 2 MPa between differential stresses of 93 MPa and 133 MPa, then steps of 1 MPa until a differential stress of 165 MPa was reached. At this stress level, the axial stress was 200 MPa, the maximum available on the Hades rig. We acquired a series of tomograms and measured negligible creep of the sample. Because of the limited time available at the ESRF for this experiment, the confining pressure was reduced by steps of 1 to 23 MPa to increase creep rate. After a period of 4 min at 23 MPa confining stress, corresponding to 172 MPa differential stress, several faults nucleated in the sample and propagated, during a creep burst that occurred after 9.2 hr and that lasted for around 8 min ( Figure 4). After this creep burst, the sample continued to creep. Similar to the creep burst observed in Experiment M83, we define an acceleration of creep between two periods of slower creep under constant differential stress as a creep burst. Then, the confining stress was reduced by 2 MPa, leading to a differential stress of 170 MPa. At this stress level, the sample continued to creep until a macroscopic failure occurred after 9.6 hr. During this macroscopic failure, the sample collapsed by slip along several faults, leading to a final axial strain of 0.6 (gray rectangle in Figures 2c and 4). Because the final deformation occurred after more than 9 hr with constant pore fluid pressure and the sample was water saturated before the experiment, we assume that the fluid pressure was equilibrated in the sample before the creep burst occurred. We did not measure the permeability of our samples. However, the permeability of Carrara marble is low, in the range 10 −18 to 10 −17 m 2 (e.g., Delle Piane et al., 2015). Therefore, the sample likely experienced effectively undrained conditions throughout most of the experiment. Consequently, the water did likely not escape from or enter into the sample during the creep burst and the final failure event.
Only three tomograms were acquired during this slow faulting episode, from Scan Numbers 94 to 96 ( Figure 2d). Because of sample deformation during these scans, the 3-D volumes are blurred and so we could not robustly quantify the void space ( Figure S2 in the supporting information). Despite the low number of  Figure 2b. The transient acceleration of creep (i.e., creep burst between 9.2 and 9.3 hr), observed macroscopically, corresponds to the nucleation and slow growth of a fault network that hosts system-spanning conjugate faults, observed with X-ray tomography. Each circle indicates the acquisition of a 3-D tomogram. Insets show 3-D views of the sample. Local macroscopic axial strain rates are indicated above the axial strain curve. After 9.6 hr the sample failed. A vertical slice of the sample during the creep burst is shown in Figure S2. scans and the blurring bias, this sample is significant because a creep burst is also observed in a second sample, under different experimental conditions (presence of water and confining pressure of 35-23 MPa for Sample M84 and dry condition and confining pressure of 30 MPa for Sample M83). Together, these two samples show the reproducibility of the existence of creep bursts in the creep regime and that these creep bursts are related to the nucleation and slow growth of faults. Movies S1 and S2, provided as supporting information, display time-lapse 3-D rendering of the samples during the experiments.
Our experimental setup contains an important improvement compared with a previous study where Carrara marble was deformed to failure under increasing stress conditions . In the present study, a LVDT sensor was added on the Hades rig. For both samples, we measured the macroscopic axial strain by two independent techniques. First, the displacement sensor measures the axial displacement of the upper piston, which records the shortening of the sample (LVDT data in Figures 2b, 2d, and 3c). These data are corrected from the elastic deformation of the rig. Second, from the 3-D tomograms we measured the height of the sample as a function of time (XCT data in Figures 2b, 2d, 3a, and 3b). Both measurements techniques show similar results. The experiments were stopped after axial strains of 0.11 (Sample M83) or 0.6 (Sample M84, after failure) were reached. After the experiment was stopped, the samples were removed from the Hades rig. For both experiments, sample extraction led to their destruction and only a powder could be recovered.

Porosity Imaging and Quantification
Because of the contrast in X-ray attenuation between the calcite grains of the Carrara marble and the air-filled voids, we can separate the microfractures and pores from the solid grains ( Figures 5 and 6). This procedure (i.e., segmentation) extracts the porosity of each scan. In our samples, the initial porosity detected by the segmentation procedure is less than 0.05% and is lower than the porosity of the sample of the order of 1%. The actual porosity was estimated by (1) measuring the weight of the dry samples, (2) water impregnation of the samples under vacuum for 24 hr and measuring the weight again, and (3) calculating the porosity assuming that Carrara marble contains 100% calcite with a density of 2.71 g/cm 3 . The difference between porosity measured by water impregnation and porosity measured in the 3-D images arises because some pores have dimensions smaller than the spatial sampling of the X-ray tomography images. In the following analysis, we report the porosity calculated from the X-ray computed tomography images, which we call XCT  Figure 3a shows the loading conditions and axial strain of each scan. The gray scale of scans indicates X-ray attenuation, with darker shades corresponding to air-filled voids. Between these two scans, the XCT porosity increased from 0.02% to 0.1% (Figure 3a) and the sample dilated via the nucleation and growth of microfractures, shown by darker gray levels. Scale bars are both 2 mm. porosity (ϕ XCT ). This porosity only includes voids with dimensions larger than the spatial sampling of the images that is constant for all scans. Thus, most of the XCT porosity that subsequently develops arises from microfracture propagation. We used segmentation to extract the microfracture network of Sample M83. First, each tomogram was filtered using a nonlocal mean filter (Buades et al., 2005) to reduce noise in the data and thus (1) enhance boundaries between grains and voids and (2) homogeneize the gray level to reduce noise in subdomains with locally similar gray levels (calcite or pores). Then a cylindrical mask was used to select a subvolume centered in the rock sample and comprising 63% of the core sample volume, thereby reducing boundary effects by removing data near the jacket and the pistons ( Figure S1). Then, a threshold in gray level of 13,800 was applied. All voxels below this value are considered voids, and all voxels above it are considered grains. The histogram of the gray levels and the thresholding result are shown in Figure S1. We chose this segmentation workflow, rather than more elaborate procedures, such as Hessian filter (Voorn et al., 2013) or machine learning (e.g., Andrew, 2018;Zhao et al., 2020), because it is well established and contains only one parameter (the grayscale threshold). Several studies have estimated the effect of using different segmentation techniques on the variability of the segmented porosity or fracture surface area (e.g., Voorn et al., 2013;Zhao et al., 2020) and found that these techniques could estimate that the difference in the segmented fracture surface area between different techniques could be up to 20% (e.g., Figure 4 in  Zhao et al., 2020). In most studies, the segmented fracture surface area or porosity cannot be compared with a ground-truth value, so it is not possible to provide an actual error of the porosity for a given segmentation procedure, but only a variability of possible values. In one study, the variability among different segmentation techniques was assessed quantitatively because a ground-truth porosity in a synthetic data set was known exactly (Andrew, 2018). This work shows that differences between segmentation procedures are negligible when the level of noise in the data set is low (Figure 3h in Andrew, 2018). For this reason, we filtered the data with a nonlocal mean filter to reduce noise and to increase the edge contrast prior to segmentation.
In the final step, we used the thresholded image to identify microfractures as individual clusters of voids. A microfracture corresponds to a set of void voxels that are connected to each other by a face and edge or a corner and surrounded by voxels identified as solid. When a cluster of voxels is identified, it is given a unique label (Figure 7). Individual clusters of voxels correspond either to a single microfracture or to an undetermined number of connected microfractures that resemble a chocolate cornflake and that we call a microfracture cluster. The procedure of thresholding of the 3-D volume to extract microfractures and estimate the XCT porosity selects microfractures that have an aperture of several voxels and thus of the order of the spatial sampling of the images. The sample contains microfractures with aperture smaller than this resolution size ( Figure S1). Because we report results of a time series of the same sample, the error in the XCT porosity and missing smaller fractures is the same for all 3-D images.
We then calculated the evolution of the XCT porosity ( Figure 8a) and the evolutions of the number of microfractures, mean volume of microfractures, and mean distance between microfractures during the creep burst (Figure 9), as well as the histograms of microfracture clusters surface are and volume ( Figure S3). We also calculated the time evolution of the volume of the largest microfracture cluster ( Figure S3). To characterize the spatial variability of the XCT porosity along the axis of the core sample, we divided the tomograms into a series of vertical layers, each 130 μm thick (20 voxels), calculated the XCT porosity in each layer, and plotted the evolution with time ( Figure 8b).
To characterize the shape of the microfracture clusters (Figures 10  and S4), we calculated the covariance matrix of each microfracture cluster, from which the three eigenvalues λ 1 , λ 2 , and λ 3 , represent the three principal axes of the microfracture cluster (inset in Figure 10c). The mean values of these eigenvalues, the evolutions of the angle β 1 , between λ 1 and the vertical direction, and the angle β 2 , between λ 3 and the vertical direction are displayed in Figures 10a and 10b. We also calculated the distribution and the average of the microfracture flatness, defined by the ratio λ 3 /λ 2 , with flat objects having a flatness close to 0 (Figures 10c and 10d).

Digital Volume Correlation
To characterize qualitatively the volumetric and shear strain evolution in Sample M83, we calculated the incremental 3-D internal displacement vectors between successive scans using digital volume correlation analysis, implemented in the software TomoWarp2 (Tudisco et al., 2017) and following the same  Renard et al., 2019). Digital volume correlation calculates the displacement vectors inside the sample at the locations of points called nodes and uses a cubic correlation window around each node to calculate the incremental displacement done between two successive 3-D scans. The node spacing was 20 voxels (130 μm) and the correlation windows size was 10 voxels (65 μm).
Using the incremental displacement fields calculated between two tomograms, we examine both the volumetric and shear components of the strain by calculating the divergence and curl of the displacement fields, respectively. The divergence is proportional to the first invariant of the incremental strain tensor and thus represents a measurement of local volume changes. Positive divergence indicates local dilation and negative divergence indicates local contraction. The curl is a vector that characterizes the rotational component of the displacement field. The norm of this vector is used here as a proxy for incremental shear strain. In a Cartesian coordinate system where the z axis is vertical and parallel to the core sample axis, and the x and y axis are perpendicular, horizontal, and arbitrarily selected, a positive curl indicates right-lateral shear strain, and a negative curl indicates left-lateral shear strain with respect to the coordinate system. Figure 9. Evolution of the number of microfracture clusters and total microfracture volume (a) and the mean microfracture volume and mean distance between microfractures (b) as a function of scan number. During the phase with a quasi-constant creep rate, the total number of fractures (single and clusters) increases, whereas the total fracture volume remains constant. During the creep burst (Scans 141-150, shaded rectangle), the number and the mean volume of microfracture clusters increase whereas the mean distance between microfracture clusters decreases. Following the creep burst, the number of fractures, mean fracture volume, and mean distance between fractures all tend to plateau, indicating the slowing or stalling of fracture network development (i.e., fault locking) that occurs as the largest microfracture cluster starts developing and occupies up to 50% of the porosity (Figures 7 and S4).

Journal of Geophysical Research: Solid Earth
To visualize and quantify the strain increments we proceed in two steps. First, we show the localization patterns of the high (>95th percentile) incremental strain values (Figure 11 and Movie S3). Second, we show the evolution of the cumulative values of these values with time, as deformation progresses (Figure 12). We show the cumulative values, rather than the incremental strain values, in order to compare this evolution to the trends in XCT porosity, which is a cumulative property.

Macroscopic Strain and Deformation Pattern
The macroscopic differential stress and axial strain relationships show mechanical behavior typical of a triaxial deformation test during the initial stages when axial stress was increased (Figure 2). In the first step, the differential stress versus axial strain is nonlinear, concave upward, due to the closure of voids and the settling of the sample between the two pistons. The second step includes a linear increase of axial strain with differential stress. The third stage begins when a yield point is reached after around 2-3% axial strain Figure 10. Statistics of the geometry of microfracture clusters between a period of slower creep when the differential stress was constant and equal to 154 MPa (Scans 90-130) and faster creep (i.e., creep burst and continuous creep after the creep burst) when the differential stress was constant and equal to 159 MPa (Scans 139-162). (a) Evolution of the average lengths of the three principal eigenvalues of the covariance matrix of all microfractures in the sample, λ 1 > λ 2 > λ 3 . During the creep burst, all three values increase indicating that microfractures grow along each dimension. Figure S3 shows the probability density function of the microfracture clusters at each scan. (b) Evolution of the average angles β 1 and β 2 between the horizontal axis and λ 1 and λ 3 , respectively. The evolution of the angle β 2 during the creep burst indicates that the smallest axis of the microfractures and microfracture clusters becomes oriented at an angle of 60°with respect to σ 1 , indicating that the microfractures' planes tend to become oriented at 30°, a direction more favorable for shear slip.  (Table 1), followed by a fourth stage during which microfractures grow under constant differential stress (Figures 2b and 2d). We defined the yield point when a deviation of 3% from linearity occurs in the stress-strain curve. The fourth and final stage includes an acceleration of macroscopic creep driven by microfracture propagation and coalescence under either constant differential stress (Sample M83) or slightly increasing differential stress (M84). For both samples, conjugate sets of faults developed (insets in Figure 2 and Movies S1 and S2). The main difference between the two samples is that for Sample M83, deformed without pore fluid, the conjugate faults developed in about 20 min. In contrast, in Sample M84, including water, the faults became system-spanning in 8 min. Because we could acquire only three scans of Sample M84 during fault initiation and propagation (Movie S2), and the deformation of the sample during the scan acquisition blurred the images ( Figure S2), we do not characterize the XCT porosity and incremental strain evolution of this sample in the following analysis. We focus this analysis on Sample M83 for which we acquired twenty tomograms during creep deformation at constant differential stress (Movie S1). No significant fracture propagation or deformation was detected during each scan acquisition on this sample. Figure 3 shows two periods of deformation of Sample M83 (dry), corresponding to the two gray rectangles indicated in Figure 2a. Under a constant differential stress of 154 MPa (Scans 90 to 130), the sample creeps at a quasi-constant strain rate close to 1.5·10 −6 s −1 (Figure 3a), between Scans 90 and 130. The nucleation and growth of microfractures produces the macroscopic creep ( Figure 5). Twinning at the grain scale could accommodate part of the deformation. However, we could not identify this mechanism with the X-ray microtomography data. The quasi-linear trend in the increase of strain with time suggests that this stage can be defined as the transition between the stage of strain rate deceleration and acceleration (Figure 1).
Under a differential stress of 159 MPa (Scans 139-162), an increase of axial creep rate from 2.4·10 −5 to 5.9·10 −5 s −1 , defined as a creep burst, is observed from Scans 141-150 (Figure 3b). The creep acceleration correlates in time with the nucleation and growth of a fault network that spans the rock core, including smaller conjugate faults ( Figure 6 and Movie S1). From Scans 139-142, the XCT porosity increases, signaling an acceleration of fracture dilation (Figure 8a). Then the fractures that later comprise the system-spanning network begin to develop at Scan 144, and continue to propagate until Scan 150 (Figures 6 and 7). Between Scans 150 and 161, the main fault network spans the system and the macroscopic strain rate decreases from 5.9·10 −5 s −1 to 1.9·10 −5 s −1 .
In Sample M84 (wet), the strain rate at 170 MPa differential stress is equal to 5.0·10 −5 s −1 and deformation occurs by the nucleation of microfractures in the volume (Figure 4). After 9.2 hr, a creep burst starts with a strain rate of 6.4·10 −4 s −1 and lasts for 8 min. During this period, several conjugate faults develop and accommodate most of the deformation. After this creep burst, the strain rate decreases to 1.4·10 −4 s −1 until the sample failed at 9.7 hr (Figure 4).

Evolution of XCT Porosity
The XCT porosity in Sample M83 (dry) evolves during deformation (Figures 6-8, S3, and S4). Microtomography images reveal that fault growth occurs as the formation of a zone of increased XCT porosity in the top half of the sample (Figures 6, 7, and 8a). This dilating volume propagates across the

10.1029/2020JB020354
Journal of Geophysical Research: Solid Earth sample until it reaches the boundary of the rock core. The XCT porosity of the entire sample increases from less than 0.5% before deformation to 2.5% after the fault network spans the core (Figures 7 and  8a). At the beginning of the creep burst, the porosity is made mainly of individual microfractures (Figure 7b), which then coalesce to form microfracture clusters (Figure 7a). The largest cluster occupies up to 50% of the total XCT porosity at the end of the creep burst when a fault network has crossed the sample ( Figure S4).
After the fault network has crossed the sample, the strain rate decreases. The time evolution of the XCT porosity along a vertical profile in the sample shows that XCT porosity increases mainly in the middle of the sample, at some distance from the two pistons ( Figure 8b). The XCT porosity in the sample increases with time during fault propagation (Scans 141-150) and below the fault network, after the fault has stopped propagating. In the upper part of the sample, porosity stops to increase of the sample once the fault has developed (Scans 150-162). The rock sample continues to creep at a lower rate during this last period.
Segmenting the voids (fracture and pores) from the solid rock produce evolving statistics of microfracture cluster geometry (Figures 9  and 10). During the quasi-constant creep rate stage (Scans 90-130 in Figure 3a), the mean values of the three eigenvalues of the covariance matrix, λ 1 , λ 2 , and λ 3 , increase slowly (Figure 10a). Idealizing a microfracture or microfracture cluster as a perfect ellipsoid, the three eigenvalues correspond to the three main axes of the ellipsoid, with the smallest eigenvalue, λ 3 , corresponding to the fracture aperture. During Scans 90-130, a phase with a quasi-constant creep rate before the creep burst, the average angle between the largest principal axis of the microfractures, λ 1 , and the horizontal plane, β 1 , remains constant, around 64° (  Figure 10b). The angle between the smallest principal axis of the microfractures (λ 3 ) and a horizontal plane, β 2 , also remains constant, above 70°. This high value of the angle β 2 indicates that the smallest dimensions of the microfractures is almost perpendicular to the vertical direction and thus that the long axes of the microfractures are oriented almost parallel to σ 1 . The flatness, corresponding to the ratio between the smallest and the intermediate eigenvalues of the microfractures, is constant and equal to 0.2 ( Figure 10c). Thus, the microfractures are generally flat objects with a penny-shaped geometry.
During the stage of the experiment with the fastest rate of macroscopic creep (Scans 142-150 in Figure 3b), that is, the creep burst, the number, and mean volume of microfracture clusters increase ( Figure 9) and the histograms of volume and surface of microfracture clusters also show an increase ( Figure S3). As a consequence, the mean distance between microfracture cluster decreases (Figure 9b). During the creep burst, microfractures coalesce into a large microfracture cluster that underline the main fault zone and contains around 50% of the porosity of the sample (Figures 7 and S4). The three eigenvalues of the microfracture clusters and their flatness increase (Figures 10a, 10c, and 10d). This observation indicates that the microfracture clusters grow in length and width, explaining the increase of XCT porosity and thus dilatancy. Concurrently, the mean orientation of the microfractures evolves: the angle β 1 increases to 68°, while the angle β 2 decreases to 60° (Figure 10b). Because the flatness of the voids is small, less than 0.3, these voids have a flat shape with λ 1 , λ 2 ≫ λ 3 . With such geometry, the longest axis of the voids could be either vertical or horizontal, with an angle β 1 that can therefore be either close to 90°or close to 0°for a vertical microfracture. Therefore, the angle β 2 of the smallest axis, λ 3 , is the most relevant geometrical parameter to estimate the orientation of penny-shaped microfractures. High angles of β 2 indicates that the aperture of the fracture is perpendicular to σ 1 or that the penny-shaped fracture is close to being vertical (inset in Figure 10b). The evolution of β 2 observed in Figure 10b shows a rotation of the microfracture orientation, with the smallest dimension oriented at 60°to σ 1 . This orientation corresponds to one of the other two axes of the microfractures oriented at 30°to σ 1 , the optimal orientation for shear faulting when considering an internal friction coefficient equal to 0.6 (Anderson, 1905). During the creep burst, the angles β 1 and β 2 of the largest microfracture cluster

10.1029/2020JB020354
Journal of Geophysical Research: Solid Earth evolves as well from orientation of microfractures almost parallel to σ 1 (β 1 close to 90°and β 2 close to 0°, Figure 7) to an increase of β 1 to 40°and then 60°while β 2 increases and stabilizes close to 50°when the fault zone is fully developed at Scan 162.

Local Strain Accumulation
Digital volume correlation analysis characterizes the incremental volumetric and shear strains during creep deformation. In Sample M83 (dry), the highest magnitudes of the incremental local strain (>95th percentile) occurred pervasively in the volume and did not become localized until Scan 141. In particular, between Scans 90 and 140, high magnitudes of dilation and shear strain concentrated in the lower half of the sample (Figure 11a). During this period, the number of microfractures increased but the total volume of fractures did not increase significantly ( Figure 9a). As the fault network nucleated, a band of high magnitudes of dilation, and a more localized zone of shear strain, formed in the upper part of the sample (Figure 11a, between Scans 139 and 141, Movie S3). Then, larger volumes of high dilation and shear strain developed along the system-spanning fault network between Scans 145 and 149 (Figure 11b). Then, slip along the fault slowed as well as the rate of new fault propagation (Figure 9a), and the highest magnitudes of incremental deformation occurred mostly in the lower part of the sample (Figure 11c) between Scans 150 and 162. This increase of incremental shear strain (i.e., creep burst) corresponds to the development of the sample-spanning fault network ( Figure 6) and the increase of XCT porosity in the fault zone (Figure 8b). At the same time, the number of microfractures increased as well as their volume (Figure 9). The mean distance between microfracture centroids increased as well, indicative of a global dilation in the fault zone. Once the fault locked, after Scan 150, the XCT porosity continued to increase, mainly in the upper part of the sample (Figure 8b). We infer the fault locked from the absence of high magnitudes of shear strain near the fault zones after the creep burst (Figure 11c).
To characterize the accumulation of incremental strain in the sample, we calculate the evolution of the cumulative values of dilation, compaction, and curl using digital volume correlation analysis (Tudisco et al., 2017). First, we sum all these values in each digital volume correlation calculation. Then, we sum cumulatively between successive scans. This procedure allows comparing the digital volume correlationresults with XCT porosity evolution. These cumulative sums over time of dilation and shear increments show a similar evolution (Figure 12). From Scans 90 to 141, the cumulative dilation increases in the sample, while the shear strain magnitude remains low. Then the cumulative incremental dilation increases, as well as compaction and shear strain between Scans 141 to 145, concurrent with the XCT porosity increase in the sample. From Scans 145 to 149, both incremental dilation and shear strain accelerate, corresponding to the growth and slip on the sample-spanning fault. After Scan 150, slip on the fault decreases and both cumulative volumetric and shear strain continue to increase at a smaller rate in the lower part of the sample (Figure 11c).
These data show that the acceleration of macroscopic creep observed between Scans 141 and 150 occurred by (1) dilation of the sample along a nascent shear zone, (2) propagation of this shear zone as a quasi-planar structure that crossed the sample and slipped, and (3) the slowing of slip on this fault network, slowing of the rate of fracture nucleation and diffuse deformation in the rock volume around it. This microscopic evolution characterizes the macroscopic creep bursts observed in the strain data (Figure 3b).

Creep in Carrara Marble
Our experimental results share similarities with previous experiments of brittle creep in Carrara marble. At high temperature, above 400°C, several mechanisms produce creep of Carrara marble including the formation of twins, grain boundary sliding, and activation of dislocation displacements (e.g., Quintanilla-Terminel & Evans, 2016). However, at lower temperature, the creep mechanism in marble is mainly due to the growth and coalescence of microfractures (Tal et al., 2016). Triaxial creep experiments with cyclic loading on marble at room temperature have shown that the macroscopic axial strain increases with time with a slightly nonlinear trend, under a differential stress loading equal to 60% of the short-term strength of the material (Figure 7 in Yang et al., 2015). Two different viscous dissipation processes were proposed by these authors that could produce this nonlinearity: a viscoelastic term with an exponential dependence in time and a viscoplastic term with a linear dependence in time.

Journal of Geophysical Research: Solid Earth
Below 50% of the differential stress at macroscopic failure, creep in marble was below the detection limit of triaxial compression laboratory experiments (Liu & Shao, 2017). At room temperature, 30 MPa confining pressure and 150 MPa differential stress, the creep rate of Jinping Bed marble in China was close to 10 −7 s −1 (Liu & Shao, 2017). For other marble samples from the same area, the creep strain rate measured at 35 MPa confining pressure, 145 MPa differential stress, and room temperature, was smaller, around 2·10 −8 s −1 (Yang et al., 2015). In our experiments, the stage with quasi-constant creep rate before the creep burst is in the range 1.5·10 −6 to 5.0·10 −5 s −1 under constant differential stress and room temperature (Figures 3a and 4). These values are higher than those reported by Yang et al. (2015) and Liu and Shao (2017), probably because our creep rates are measured at a differential stress closer to failure than in their experiments. However, the main difference is that we observed creep bursts in our data, with strain rates in the range 5.9·10 −5 to 6.4·10 −4 s −1 (Figures 3b and 4), several orders of magnitude larger than creep rates previously reported in marble (e.g., Liu & Shao, 2017;Yang et al., 2015).
In series of laboratory experiments performed on Carrara marble, the confining pressure was varied to study the deformation of either pristine rock cylinders (Fredrich et al., 1989) or prefaulted cylinders (Meyer et al., 2019). A transition between brittle and semibrittle deformation was observed for confining pressures in the range 30-50 MPa. Under these confining pressures, both microfractures and plastic micromechanisms of deformation, such as twinning and movement of dislocations in the calcite crystals, accommodated deformation categorized as semibrittle. Creep bursts in our experiments are observed when the confining pressures are 30 MPa (Sample M83) or~23 MPa (Sample M84), close to the lower limit of confining pressure where the brittle-ductile transition of Carrara marble was measured at room temperature (30 MPa, Meyer et al., 2019). Because our samples were destroyed during unloading from the Hades rig, we cannot accurately identify the micromechanisms that operated, such as twinning, dislocations, or pressure solution creep. Future experiments in which the sample is preserved after a creep burst would be necessary to quantify the partitioning between ductile and brittle micromechanisms because other experiments on carbonate rocks have shown that pressure solution creep, crystal plasticity, and twinning could participate to macroscopic creep at room temperature (Brantut et al., 2014a;Nicolas et al., 2017). However, our results suggest that the creep bursts are mainly due to the nucleation, growth, and coalescence of microfractures, and if other plastic deformation mechanisms exist, they accommodate a minor fraction of the total strain. Therefore, our results agree with results on prefaulted Carrara marble samples deformed near the brittle-ductile transition, which show that at confining pressures below 30 MPa, brittle deformation and fault slip accommodate more than 50% of the total strain, and pervasive ductile deformation of the matrix accommodates less than 50% of the total strain ( Figure 2 in Meyer et al., 2019).

Brittle Creep and Faulting in Rock Experiments
Other laboratory experiments performed on crystalline low porosity rocks have measured strain rates during brittle creep. Lockner (1993) deformed granite samples at room temperature, under constant confining and pore pressure conditions, as the axial stress was cycled to enhance creep. Macroscopic axial strain increases with a linear dependence of time when the applied differential stress is at some distance (80%) of the failure stress, that is, constant creep rate (Figure 1). This constant creep strain rate follows an exponential dependence with the applied (constant) stress (Lockner, 1993). From the recording of acoustic emissions and postmortem observations of microfractures, the nucleation and subcritical growth of microfractures have been identified as dominant microstructural mechanisms that accumulate irreversible strain (Lockner, 1993;Ross et al., 1983;Scholz, 1968). In our experiments, constant creep rate occurs under a differential stress of 154 MPa in Sample M83 (Scans 90-130) and coincides in time with the nucleation and growth of microfractures (Figures 3a and 5). During this deformation stage, the XCT porosity increases slightly (Figure 3a). The average shape of microfracture remains similar, and their orientation is consistently near parallel to σ 1 with an angle β 2 close to 70° (Figure 10). These dilatant vertical cracks formed in the sample produce irreversible damage. With continuing creep, when approaching fault formation, these microfractures coalesced to form a fault network that extended across the sample. These stages of fracture coalescence were also observed in 2-D experiments on Carrara marble (Tal et al., 2016).
In a series of experiments on crystalline rock performed under constant differential stress, Lei et al. (2000) measured a nonlinear increase of the rate of acoustic emissions that occurred within 1 min in a rock sample loaded near failure. The spatial distribution of the acoustic emissions indicates that they were initially 10.1029/2020JB020354 Journal of Geophysical Research: Solid Earth distributed throughout the rock and then concentrated in a quasi-planar structure that evolved into a more localized fault. The zone with tensile fracturing located in front of the propagating fault, identified from acoustic emission recording, supports the existence of a dilatant process zone. The fault then propagates within this process zone, a mechanism proposed in Reches and Lockner (1994). In porous rocks, such as sandstone, acoustic emission recording show that creep remains diffuse in the volume, without localizing into a subplanar structure (Heap et al., 2009). In the present study, strain localization occurs during a creep burst, defined as an acceleration of creep between two periods of slower creep rate and before a final failure may occur (Figure 1). In Sample M84 (final axial strain of 60%), this final failure was reached whereas it was not reached in Sample M83 because the experiment was stopped. During the creep burst in Sample M83, the dilatant (high porosity) zone forms in the middle of the sample (Figures 6, 7, and 11, Scans 143-145), with low magnitudes of cumulative shear strain ( Figure 12) and higher local porosity increase. Then high magnitudes of shear strain localized into a system-spanning fault (Scans 143-147 in Figure 11). After the fault network crossed the sample, slip along it decreased, the rate of fracture propagation slowed, and creep continued in the lower half of the core (Figures 9b and 11c, Scans 150-162). This evolution characterizes the creep bursts reported for both M83 and M84. In both samples, after the creep burst, the sample continued to deform until macroscopic failure.
The development of faults in laboratory creep experiments has been observed in granite (Lei et al., 2000) and basalt (Heap et al., 2011). In amphibolite (Satoh et al., 1996), a fault developed along a preexisting joint that acted as a nucleation site for localized deformation. In these experiments on granite, basalt, and amphibolite, an increase of the number and the energy of acoustic emissions occurred when approaching failure. This increase was interpreted as the development of microfractures at the grain scale. In Experiment M83 (dry), we directly observe that fault growth occurs slowly during a creep burst, with varying rates of creep, rather than a faster evolution toward catastrophic failure. Moreover, once the fault network spans the system, slip on it decreased, and creep deformation was accommodated in the lower part of the sample by the opening of new microfractures, rather than through accelerating slip along this throughgoing fault network. Conversely, in Sample M84 (wet), the creep burst precedes macroscopic failure, and occurred during the accelerating creep phase.
An important difference between our two experiments is that the macroscopic strain rate during the creep burst is around ten times faster in Experiment M84 (wet) (Figure 4), than Experiment M83 (dry) (Figure 3b). Two effects could explain the role of water in changing the deformation rate: the influence of pore pressure or chemical weakening. On one hand, the low permeability of the Carrara marble likely causes the experiment to experience undrained conditions throughout loading. Therefore, when microfractures open during the creep burst, the pore fluid pressure might not have time to equilibrate in the sample. Locally, this effect would reduce the pore fluid pressure inside the sample and therefore increase the confining pressure locally. Increasing confining pressure should stabilize fault development (e.g., Brantut et al., 2013), bringing the sample closer to the ductile transition. Moreover, the hardening effect expected during dilatant deformation of wet rocks (e.g., Rice, 1975) should also stabilize fault development. Therefore, the pore pressure effect due to the undrained conditions is expected to stabilize deformation and not accelerate it as observed in Experiment M84.
On the other hand, the influence of water on accelerated creep has been observed in granite (Kranz, et al., 1982), sandstone (Baud et al., 2000), and single calcite crystals (Røyne et al., 2011). The presence of water lowers the surface energy between calcite and water during fracture propagation, producing faster rates of fracture growth, as observed in Experiment M84. Another explanation besides the reduction of surface energy is that water activates pressure solution creep (e.g., Gratier et al., 2013). However, due to the limited spatial sampling size, microstructures responsible of such effect would not be detectable in our tomography data.

Creep Burst and Fault Locking
Creep experiments in Darley Dale sandstone (Figure 1 in Heap et al., 2009)  damage and microfracture evolution was not possible. Therefore, the deformation mechanisms were identified by imaging thin sections of the samples after unloading and removing the samples from the deformation rigs.
In our two experiments, the creep burst corresponds to an acceleration of creep coincident with the slow nucleation and growth of a system spanning fault network that we can image using dynamic in situ X-ray microtomography. In Experiment M83, after the fault network grew (Scans 141-150), shear along the fault stopped and creep was accommodated by dilation in the bottom part of the sample (Scans 150-162), see Figure 11. This observation indicates that the fault network locked at the end of the creep burst, and volumetric strain concentrated below this fault network, producing a lower creep rate. In a series of XCT experiments on granite saw-cut samples where the initial fault surface roughness was varied, our group has observed that rougher faults could lock and produce off-fault damage . The interpretation was that asperities along a fault surface with a larger roughness amplitude would lock, allowing larger stress concentrations to propagate in the off-fault volume and creating damage there.
Therefore, the mechanism of fault locking could be related to an interplay between the interlocking of stronger asperities and off-fault deformation, a process proposed to explain how the slip rate along active faults may vary with time and space (e.g., Herbert et al., 2014). A switch between fault reactivation (or locking) and off-fault ductile deformation has been observed in laboratory experiments where saw-cut samples of Carrara marble were deformed under different confining pressures (Meyer et al., 2019). In our experiments, the confining pressure was constant during the creep bursts and fault locking was still observed. This may indicate that, along the damaged fault zone, material properties have evolved such that the rock became stronger than the surrounding host rock where deformation continued.

Geophysical Implications and Conclusion
The experiments reported in the present study indicate that, under constant stress conditions at room temperature and confining pressure below 30 MPa, cores of Carrara marble deform macroscopically in creep due mainly to the growth of microfractures. Both quasi-linear and accelerating creep occur in the experiments (Figures 1, 3, and 4). The stage of accelerating creep includes creep bursts: transient accelerations of creep, similar to slow slip events measured on continental faults that have been interpreted as an acceleration of slip along a preexisting fault plane (Jolivet et al., 2013(Jolivet et al., , 2015Linde et al., 1996;Rousset et al., 2016). In our experiment, the creep bursts did not occur on a preexisting fault but coincided with the slow development of a system-spanning connected fault network. Our data therefore demonstrate that a creep burst, in laboratory experiments, does correspond not only to a transient slip acceleration on a preexisting interface, as often interpreted in field studies of active faults, but also to the slow nucleation and growth of new faults. At the microscale, the creep burst coincides in time with the (1) localization and coalescence of microfractures along a quasi-planar structure, (2) increase of the aperture of the microfractures, (3) rotation of the orientation of the microfractures to 30°to σ 1 , more favorable for shear deformation, and (4) localization of the highest magnitudes of the local dilation and shear strain into a subplanar inclined structure. Once formed, slip on the fault network decreased, the high magnitudes of strain delocalized from this fault zone, the rate of fault propagation slowed, and creep continued in the volume around the fault. Our results therefore demonstrate that creep bursts in laboratory experiments may indicate the birth of a new fault and not only an acceleration of transient slip on a preexisting fault, which is the common explanation for slow slip events observed in active faults.