Fiber Optic Sensing for Geomechanical Monitoring: (2)- Distributed Strain Measurements at a Pumping Test and Geomechanical Modeling of Deformation of Reservoir Rocks

In this study distributed fiber optic sensing has been used to measure strain along a vertical well of a depth of 300 m during a pumping test. The observed strain data has been used in geomechanical simulation, in which a combined analytical and numerical approach was applied in providing scaled-up formation properties. The outcomes of the field test have demonstrated the practical use of distributed fiber optic strain sensing for monitoring reservoir formation responses at different regions of sandstone–mudstone alternations along a continuous trajectory. It also demonstrated that sensitive and scaled rock properties, including the equivalent permeability and pore compressibility, can be well constrained by the combined use of water head and distributed strain data. In comparison with the conventional methods, fiber optic strain monitoring enables a lower number of short-term tests to be designed to calibrate the parameters used to model the rock properties. The obtained parameters can be directly used in long-term geomechanical simulation of deformation of reservoir rocks due to fluid injection or production at the CO2 storage and oil and gas fields.


Introduction
Injection of fluids into underground formation, as in the case of geological sequestration of greenhouse gases, can alter the stress state in the porous skeleton of the medium as well as the formation pressure.Stress redistribution within and around the reservoir and caprock system may lead to geophysical changes, ground surface deformation [1,2], microseismicity, fault reactivation, and may even induce damaging earthquakes [3][4][5].In recent years, following the rapid increase in applications, such as injection of fluids into deep formations of the earth's crust, including geological sequestering of CO2, enhanced geothermal systems (EGS), the extraction of unconventional oil and gas, and wastewater disposal, injection-induced earthquakes and other risks related to injectioninduced surface deformation and fault reactivation have attracted increased attention from scientists and the public alike [5][6][7][8][9][10].
Geomechanical effects induced by reservoir production can be particularly pronounced in stress-sensitive reservoirs, such as poorly compacted or fractured reservoirs [11].For both monitoring and risk assessment purposes, it is important to assess the influence of injection operations so that the geomechanical implications, particularly in terms of the potential for fault reactivation, can be adequately assessed [12].How to predict and control the risks associated with fault reactivation, including fluid leakage and damaging earthquakes, are key issues to allow fluid injections to be carried out safely and accepted by the public.In these regards, geomechanics plays a role in bridging the gap from geophysics to engineering and is also a key issue in site selection, injection operation, and post-injection management in CO2 storage and other applications.Given that changes of deformation in a reservoir-caprock system are directly linked to changes in pore pressure and stress state, coupled fluid flow and geomechanical modeling is potentially a powerful tool to understand and predict the mechanical behaviors of the reservoir-caprock system.It can contribute to the characterization, prediction, and optimization of injection operations.In recent years, significant progresses in numerical modeling of reservoir scale have been made for a number of purposes such as quantification of challenge in geothermal development [13], analysis of induced seismicity by shale fracking [5], stimulation of heavy oil reservoir [14], and CO2 geological storage [15].
In order to make accurate predictions, a geomechanical model needs to be first calibrated against geophysical and geomechanical data obtained from monitoring.As a new technology, continuous and distributed fiber optic sensing (DFOS) technology based on Brillouin and/or Rayleigh scattering has great potential in geotechnical monitoring at both laboratory and field scales [16][17][18][19][20]. Optical fiber sensors have advantages such as their immunity against electromagnetic interference, low weight, small size, high sensitivity, and large bandwidth.By grouting fiber optic cables into multiple boreholes, a large amount of accurate, spatially resolved data can be obtained.Xue et al. [21] carried out a fundamental study on distributed strain measurements based on the Rayleigh and Brillouin frequency shifts of two different sandstone samples (Berea sandstone and Tako sandstone) under controlled hydrostatic confining and pore pressure in the laboratory.By comparing with strain data obtained from conventional strain gauges, they showed the potential use of deformation monitoring with the DFOS technology in complex lithology reservoirs.
As a sequel to the study of fiber optic sensing for geomechanical monitoring by Xue et al. [21], this paper presents a field study of distributed strain measurement in a pumping test and also discusses the roles of distributed strain data and a combined analytical and numerical approach in providing the scaled-up formation properties for geomechanical modeling and analysis.Section 2 introduces a general modeling framework and the numerical simulation method used in the coupled thermal, hydrological, and mechanical (THM) analysis.Section 3 describes the field case study results from the pumping test.Finally, Section 4 summarizes the future practical uses of fiber optic sensing in geomechanical monitoring at the CO2 storage and oil and gas fields.

General Framework of Geomechanical Modeling
The governing equations for heat and fluid flows, and mechanics are equations for mass and energy conservation, and momentum balance.The mass and energy conservation equations are of a similar form [22]: In the mass-conservation Equation (1), the subscript j denotes a particular fluid phase, mj (kg/m 3 ) is the unit fluid mass, wj (kg/m 2 /s) is the unit mass-flux, and qj (kg/m 3 /s) is a unit source term of fluid phase j.In the energy-conservation Equation (2), the subscript θ denotes the heat component, mθ (J/m 3 /s), fθ (J/m 3 /s), and qθ (J/m 3 /s) are unit heat, heat flux, and heat source, respectively.The volumetric flux of fluid phase j,  (/) =   ⁄ , is assumed to follow the multiphase Darcy's law: where g is the vector of gravitational acceleration, k (m 2 ) is the absolute permeability, krj, μj (Pa•s), and ρj (kg/m 3 ) are the relative permeability, viscosity, and density of fluid phase j, respectively.The heat and heat flux are given by where ∅,  ,  (J/kg/K), and  (J/m 2 /s/K) are the porosity, density, heat capacity, and thermal conductivity of the porous media, respectively. , and ℎ (J/kg), are the saturation and specific enthalpy for fluid phase j, respectively.The linear momentum balance equation is given by where, σ is the Cauchy total-stress tensor, and ρb is the bulk density.Based on linear poroelasticity theory [23], the coupled governing equations for pore pressure (here, the pore is saturated with a single "effective fluid") and solid mechanics can be written as [24] where u (m) and p (Pa) are the displacement vector and excess pore pressure, respectively.λ and μ are Lamé coefficients, α is the dimensionless coefficient representing the change in pore pressure per unit change in bulk volume under drained conditions, Q −1 (Pa -1 ) is the bulk compressibility, D is Darcy conductivity.The f (Nf/m 2 ) and q (m 3 /s) are body force and fluid source, respectively.By assuming the infinitesimal transformation, the strain tensor (ε) is given by the symmetric gradient of the displacement vector (u): As seen from these equations, fluid flow, heat transfer, and the mechanic responses of the reservoir systems are coupled with each other.Its solution usually requires a coupled modeling approach, thus, coupled THM analysis is important in geomechanical modeling [25].
Figure 1 shows a schematic of the geomechanical modeling framework, which contains a coupled THM simulation and history matching for calibration of the numerical model.Although full coupling schemes are more accurate, they are also more complex.On the other hand, sequential coupling schemes, in which fluid and heat flows and mechanics are explicitly coupled, may be sufficiently accurate for problems in which coupling is relatively weak.The general steps of developing a numerical model for coupled THM simulation start with building a conceptual geological model with existing geological and geophysical data and then creating a numerical model with proper constitutive relations and scaled rock properties obtained through analytical approaches.
To calibrate the numerical model it needs to run the coupled THM simulation for history matching of the monitoring data by adjusting the input parameters, such as pore pressure (Pp), temperature (T), water saturation (Sw), and rock property.An accurate numerical model can be useful for long-term prediction of geomechanical integrity of a reservoir-caprock system.
Geological formations contain heterogeneities at scales down to the pore and grain scale; however, a practically useful model cannot fully account for this level of detail.There are practical needs to deal with a scaled-up model that accounts for most of the important physics in subsurface geoscience.In a scaled-up model, an element represents a volume containing heterogeneities and fractures at scales smaller than the scale of the element, down to the pore and grain scale, and these features may distribute uniformly or in a patchy fashion.Because there are uncertainties in scaling up such models, probability-based uncertainty analysis and prediction are also necessary.
Within the geomechanical modeling framework, the coupled THM simulation, which predicts injection-or production-induced changes in the rock properties, rock mass deformation, stress distribution, and the fracture or fault stability, is the key element to ensure safety of the operations.There are a number of choices of research-oriented or commercial software for reservoir simulation and/or stress analysis.In this study, TOUGH2 and FLAC3D were selected for these purposes.TOUGH2 is a package of programs for multiphase reservoir simulation, and FLAC3D [26] is a commercial software for geotechnical analysis.With sequential couplers [27], this combined TOUGH-FLAC approach has proven useful in the analysis of deformation coupled with fluid flows within hard and soft rocks in injection-induced seismicity in shale gas production [5], geological CO2 storage [28][29][30], geothermal and natural analog studies [31,32].
The permeability of rock mass may change greatly as a result of deformation and fracturing [12,33].In previous works, the permeability has been expressed as a function of volumetric strain εv or shear strain εs [32,34,35] as where φ is porosity, k (m 2 ) is permeability,  and  are their respective initial value, and n is constant integer.

Field Case Study of Pumping Test
In Japan, the major potential sites for CO2 storage are characterized as soft rock formations.Thus, the efficiency of the TOUGH-FLAC approach should be examined for such sites.To this end, a pumping test with distributed strain monitoring by fiber optic sensors [20] was conducted.In this study, the water level and strain data obtained from the pumping test were used for geomechanical modeling and analysis.

Mobara Test Site
The pumping test site selected is located at Mobara city (Chiba Prefecture), north Boso Peninsula of the Kanto region of Japan, where there are many wells for agricultural use with depths of up to 200 m.The test site consisted of two pumping wells and two monitoring wells (Figure 2).A singlemode optical fiber cable was installed behind one of the monitoring wells (Well-a) and it was cemented in the annulus between well casing and the geological formation.The two pumping wells (Well-3 and Well-4) are located at a distance of 175 m and 280 m from the Well-a respectively.The other monitoring well (Well-b), located at a distance of 5.5 m from Well-a, was used to monitor the water level during the pumping test.Figure 3 shows a schematic view of geological cross-section between the wells.Sediments at the test site, from the surface downward, are of the Kasamori, Chonan, and Kakinodai Formations, all of which belong to the upper Kazusa Group [36].Kakinokidai Formation: Within the site, the Kakinokidai Formation is formed of sandy siltstone.Indistinct stratification is present in the upper part of the formation.Over a wider region, the formation coarsens upward and fines eastward.
Chonan Formation (approximately 120 m thick): The Chonan Formation is subdivided into two parts: the lower and upper parts.The lower part (approximately 40 m thick) is formed mostly of siltstone-dominated alternations with thin slumped beds.The upper part (approximately 80 m thick) is composed mostly of sand-dominated alternations with numerous slumped beds.
Kasamori Formation (154 m thick): The Kasamori Formation is composed mostly of highly bioturbated deposits such as sandy siltstones and silty sandstones.The Kasamori Formation is characterized by low permeability and can be considered as a caprock layer in the site.

Optic Fiber Measurement
The fiber cable was custom-made for this study to improve the coupling strength between casing and cement.The frequency shifts were taken using a Neubrescope NBX-8000 device (2012, manufacturer, city, country), which can measure both Rayleigh and Brillouin as in the experimental study [21].Unfortunately, only the Rayleigh shifts were thoroughly recorded during the pumping test because the signal to noise ratio of Brillouin shift was too low for measuring the strain data.During water pumping the Rayleigh frequency shift was recorded every 5 cm along the fiber cable and the frequency shift measurement was repeated in a time interval of 5 minutes.Deformation of the target reservoir formation caused by water pumping was estimated from the Rayleigh frequency shift following the same approach demonstrated in the experimental study [21].To examine the difference in hydraulic communication between Well-a (fiber sensor well) and the two pumping wells, pumping at Well-4 was paused at 14:00 on 3 November.Then pumping was restarted at 10:00 on 4 November concurrently with stopping pumping at Well-3.Compressive strain increased to 40 µε at the depth from 150 to 220 m when pumping at Well-4 continued from 10:00 on 4 November through 11:00 on 11 November.The total volumes of water pumped out from Well-3 and Well-4 were 1,382.4 and 4,233.6 m 3 , respectively.After the termination of pumping at Well-4 the compressive strain decreased gradually in the measurement time period.Later in 2015, additional pumping tests were conducted separately at Well-3 and other wells to investigate the heterogeneity of the Chonan Formation.

Data Preprocessing
The spatial and temporal sampling resolution of strain measurement during the pumping test was 5 cm and 5 min respectively and both of them were considered to be too fine for calibration of numerical model.Therefore, the Bayesian method [3837] was applied to fit the raw strain data.Splines with a fixed order of 3 were used, and the dimension number dim was determined by minimizing the Bayesian information criterion (BIC), defined by where N is the number of data points and Q is the weighted error sum of squares from a least-squares spline approximation [38].
As demonstrated by the examples shown in Figure 5, the short-wavelength noise and signals were removed by the Bayesian method.From the resulting spline functions, the strain data can be resampled with any sampling interval.

Geological and Numerical Models
As mentioned previously, formations in the test site are alternations of mudstone, silt, and sandstone.Each individual layer is typically thin as characterized by physical logging.From the observed strain data, it was decided to divide the geology into four major types: the Kasamori, upper and lower Chonan, and Kakinodai Formations (Figure 3).For this type of geology mixed with different rocks, the scaling up of rock properties is a key issue in geomechanical modeling and numerical simulations.

Two-dimensional Radial Flowing Model
The sedimentary layers in this region are quite stable and almost flat over a few tens of kilometers.Thus, it is valuable to make a rough estimate of some key parameters using simple models for which theoretical solutions are available.For this end, a two-dimensional (2D) radial flowing model of an infinite and isotropic horizontal layer was tested.Under the assumption that Darcy's law holds, the theoretical solution of the pore pressure change is given by [39] In Equations ( 12) and ( 13), t (s) is time since the start of pumping, r (m) is the distance from pumping well.Q(m 3 /s) is the pumping rate, H (m) and K (m 2 ) are thickness and permeability of the layer, D (m 2 /s) is hydraulic diffusivity, η (Pa⋅s) is dynamic viscosity of water, Sa (Pa −1 ) is unconstrained specific storage coefficient, βfl (Pa −1 ) and βpv (Pa −1 ) are compressibility of the fluid and pores, respectively.
In the present problem, all parameters except K and βpv can be estimated based on the borehole data.In this study history matching approach was applied to estimate the values of K and βpv.First, sensitivity tests were performed to determine if K and βpv are in a tradeoff relationship.As shown in Figure 6, the predicted water level response during pumping and after shutdown is sensitive to both parameters.Thus, the permeability of the layer and the compressibility of the pores, which are poorly known and show large uncertainty in many operation sites, can be determined from the observed water level data.The parameters K and βpv can be determined by fitting the model predictions to the observed data using a grid-searching algorithm.
The water level change in the upper Chonan formation during the initial pumping period from Well-3 and Well-4 were then considered.It was found that, for a single well pumping test, the water level change observed at Well-2, located 630 m away from Well-3 (Figure 2), can be represented fairly well using a single-layer 2D radial flowing model (Figure 7a).For multiple tests, when a single set of K and βpv values were used, the observed data were not well represented (Figure 7b).However, when using different sets of K and βpv values, the fitting results greatly improved (Figure 7c), indicating lateral inhomogeneity in the physical properties.K and βpv values estimated in such a manner may be considered as a type of "mean" or "equivalent" value and can be used as initial values in the coupled model.
A rough estimate of the bulk modulus of the Chonan formation can also be obtained from the maximum change in the water level and the maximum fiber strain observed during the pumping test.Under drained conditions, the volumetric strain εv caused by a change of pore pressure ΔP is correlated with the bulk modulus as where c is a constant ranging from 1 to 3.

Three-dimensional Model for THM-coupled Simulation
The geometry of the numerical model for THM-coupled simulation is shown in Figure 8.The far-field boundaries were placed at a distance of 5 km to approximate infinite boundaries.The in-site stress was applied in all zones and as loads acting on the far-field boundaries.The model has dimensions of 10,000 m × 10,000 m × 300 m centered at the fiber well.Roller boundary conditions were imposed on the four sides and bottom of the model.This model is divided into grid elements with dimensions in the horizontal plane ranging from 20 to 2500 m according to the distance from the pumping well and depth dimensions of 10 or 20 m.The total number of active elements (zones in FLAC3D) is 23,896.To count for the inhomogeneity in the Chonan Formation, two sub-models, model-A and model-B, were examined.In model-A, the Chonan Formation is divided into two parts, L-2s and L-2n, with different permeability.In model-B, the permeability was assumed to linearly decrease along the dip direction of the layers (Figure 8).

History Matching
In the initial model used in the simulations, the permeability and pore compressibility of the aquifer (L-2) were obtained from the 2D radial flow solutions.The other parameters of the aquifer and the other layers were determined based on experimental results.Table 1 lists the optimal hydraulic and mechanical properties determined using the history matching approach.Aided by the aforementioned 2D radial flow analysis results, the determination of the optimal parameters required only a few tuning steps.,b show the water level and strain, respectively, from the observed data and simulated results.The water level data were represented very well by the numerical results.The observed strain data showed small-scale fluctuations, corresponding to the sand-silt alternation structures of the Chonan Formation.The mean strain data were also represented fairly well by the numerical results.

Concluding Remarks
The key topics in geomechanical monitoring at the oil and gas fields and geophysical CO2 storage sites are (1) what depth and how much strain occurs in the subsurface, due to the fluid production or injection?(2) How to interpret ground surface uplift or subsidence as a result of deformation in the subsurface?In this study the distributed fiber optic strain measurement was carried out along the wellbore in the pumping test.The optic fiber installed behind well casing clearly detected the deformation resulting from pumping water from the target formation.Compressive strains responded well to the reservoir rock property.Higher strain magnitude was observed in the sanddominated sediment layer during pumping.This field application demonstrated the practical uses of the distributed strain sensing technology for geomechanical monitoring, which was proposed by Xue et al. [21] based on the laboratory experiment results of core samples.The distributed strain sensing technology not only helps us to understand the development of subsurface deformation caused by fluid injection or production, it also enables us to relate the ground surface subsidence or uplift to the subsurface deformation, as well as interpret how the subsurface deformation migrates to the ground surface.
In this study the observed distributed strain data, in combination with water level data, was used to estimate the scaled-up hydraulic and mechanical properties required in geomechanical modeling and analysis.Generally geological materials contain inhomogeneous structures over a wide scale range from the scale of grains and pores to that of reservoirs; the scaling up of permeability and elastic modulus is commonly required in reservoir simulations and rock mechanics.Numerical simulations are widely conducted to obtain the scaled-up elastic stiffness tensors [40,41] and permeability [42].The scaled-up parameters are expected to be equivalent to the elastic stiffness tensors that include the fine-scale contribution to coarsen the geomechanical mesh.However, such numerical simulations for scaling up rely on some unknown parameter, such as a statistical model of heterogeneities and fractures, and thus must be examined with observed data taken into consideration.The present study demonstrates that the distributed strain data are straightforward and very useful for this purpose.
When processing the fine-resolution fiber optic strain data obtained from the sandstonemudstone alternation along a continuous trajectory, we confirmed that the Bayesian curve-fitting is an effective approach for modeling purpose.In the first-order approximation, the simplified analytical 2D radial flow model of an infinite and isotropic horizontal layer provide an effective and quick way to roughly estimate the permeability and thickness of the target formation by history matching water level change at the monitoring well during a single well pumping.Multi-well pumping test results can be used to quantify permeability anisotropy in horizontal planes within a formation.It has also been demonstrated that the scaled-up sensitive rock properties, including the equivalent permeability and pore compressibility of the reservoir, can be well constrained by the combined use of water level and distributed strain data.More generally, this study illustrates that the combination of analytical and numerical methods provides a promising approach for the efficient optimization of scaled-up formation properties.In comparison with conventional methods, distributed fiber optic strain monitoring allows a lower number of short-term tests to be designed for the calibration of parameters describing the rock properties, which can be used for long-term geomechanical simulations.

Figure 1 .
Figure 1.Framework of geomechanical modeling with coupled thermal, hydrological, and mechanical (THM) simulation and history matching using distributed strain data, for which no geophysical post processors are required.

Figure 2 .
Figure 2. Map of the test site and configuration of wells for pumping, fiber sensor monitoring, and water head measurement.

Figure 3 .
Figure 3. Profile showing well locations and simplified geological columnar sections.

Figure 4
Figure 4 shows the water pumping histories from Well-3 and Well-4, the change of water level measured at Well-b, and the strain estimated from the fiber optic sensing at Well-a.Pumping started simultaneously at Well-3 and Well-4 at 10:00 on 2 November 2015, at a pumping rate of 480 L/min.Several hours after starting pumping, warm colors appeared in the upper layer of Chonan formation.The yellow and red colors indicate compressive deformation of the reservoir rocks.Variations in the magnitude and the times of compressive strain with depth in the Chonan formation suggested the different responses of the sandstone-mudstone alternation.

Figure 4 .
Figure 4. Pumping history and vertical strain observed by the fiber sensor.The spatial and temporal sampling intervals are 5 cm and 5 min, respectively.

Figure 5 .
Figure 5. Example of smoothed profiles of strain data obtained using a Bayesian method of spline fitting.The order of the spline functions was fixed at 3, and the dimension was determined by minimizing the Bayesian information criterion (BIC).

Figure 6 .
Figure 6.Water head changes estimated from the theoretical solution of a 2D radial flow in a homogeneous and isotropic infinite horizontal layer for pumping at a rate of Q = 480 L/min calculated at a distant of r = 175 m.(a) Results for different permeability (K) and hydraulic diffusivity (D).(b) Results for different permeability (K) and layer thickness (H).

Figure 7 .
Figure 7. Water head changes estimated from the theoretical solution of a 2D radial flow in a homogeneous and isotropic infinite horizontal layer.(a) Results of a single test.(b) Results of multiple tests in two wells, water head data were fitted by a single model.(c) Same to (b), but water head data were fitted by a combination of two single models to count the effect of spatial inhomogeneities.The observed data (dashed line) were represented very well by the theoretical solution with the optimal values of the permeability K of the layer and the compressibility βpv of the pores.

Figure 8 .
Figure 8. Geological and numerical models for simulation of the pumping tests at Well-3 and Well-4.

Figure 9 .
Figure 9. Observed and simulated (a) water head and (b) vertical strain (εzz).Both the water head and strain data were represented fairly well by the numerical results.

Table 1 .
Major hydraulic and mechanical properties.