Stochastic Inversion Approach Research Articles

Stochastic inversion for equivalent hydraulic fracture characterization using low-frequency distributed acoustic sensing data

Distributed fiber-optic sensing, a game-changing technology, has been successfully applied in monitoring hydraulic fracturing treatments. The cross-well fiber-optic strain/strain-rate responses are directly related to the fracture-propagation process and dynamic fracture geometry. Extracting the insightful information regarding hydraulic fractures requires robust interpretation workflows. In this paper, we aim to develop a stochastic inversion approach to characterize the dynamic equivalent fracture geometry and quantify the associated uncertainties using cross-well low-frequency distributed acoustic sensing (LF-DAS) measurements. A fast analytical forward solver based on the 3D displacement discontinuity method is developed to model the strain. The objective function is the sum-of-squares error between the modeled strain and the measured strain. We focus on an equivalent fracture geometry with a predefined length centered around the monitoring well, which is a reasonable assumption considering the spatial sensitivity of LF-DAS. Therefore, the model parameters reduce to fracture width and height of the equivalent fracture. We adopted the delayed rejection adaptive metropolis (DRAM) algorithm, an efficient Markov-chain Monte Carlo (MCMC)-based method, to perform the optimization. The efficiency and accuracy of the algorithm are tested using a synthetic case that has three monitoring wells at different locations. The mean values of fracture width and height are both close to the true values with negligible standard deviations for all the three monitoring wells. The algorithm performance is not sensitive to the initial guess and the predefined length needs to be large enough. Then, a field case with a single cluster from the HFTS-2 project is presented. The fracture width increases from about 0.5 mm to about 0.9 mm during the pumping stage and gradually decreases during the shut-in period. The fracture height changes during the pumping stage and stabilizes at about 150 m during the shut-in period, the trend of which is consistent with the LF-DAS signals from the offset vertical monitoring well. To the best of our knowledge, this is the first algorithm that can efficiently quantify the fracture width and height simultaneously. The results could provide novel insights on the characteristics of fracture geometry and significantly improve the understanding of complex fracture propagation processes in subsurface reservoirs.

Multivariate stochastic seismic inversion with adaptive sampling

We have developed a fully coupled categorical-multivariate continuous stochastic inversion with a combined petro-elastic model and convolution. The new multivariate stochastic seismic inversion approach simulates multiple reservoir properties simultaneously and conditions them to the well and seismic data at the same time through the close integration of multivariate geostatistical modeling and stochastic inversion. This approach combines a trace-by-trace (column-wise) adaptive sampling algorithm with multivariate geostatistical techniques to select reservoir properties that match the seismic data. The adaptive sampling method uses an acceptance-rejection approach to condition geostatistical models to the well and seismic data. The adaptive sampling algorithm defines a practical stopping criteria based on the inherent uncertainty due to modeling assumptions and the size of the uncertainty space. This technique samples the realizations inside the space of uncertainty; the number of realizations attempted increases with the size of the space of uncertainty. Characterizing multiple reservoir properties simultaneously through the close integration of seismic inversion and multivariate geostatistical techniques leads to improved high-resolution reservoir property models that reproduce the original seismic data. A case study is considered to compare the proposed stochastic inversion approach with the conventional methods. The case study represents multivariate stochastic inversion provides high-resolution facies and reservoir physical properties simultaneously that reproduce the original seismic data within quality of data better than the other approaches.

Impact of near-surface fault geometry on secular slip rate assessment derived from uplifted river terraces: implications for convergence accommodation across the frontal thrust in southern Central Bhutan

Vertical velocities obtained from uplifted river terrace dating near mountain fronts are commonly converted into overthrusting slip rates assuming simple geometry of the fault at depth. However, the lack of information on the dip angle of these shallow structures can lead to misinterpretation in the accommodation of convergence, and thus to erroneous conclusions on the transfer of shortening to the emergent thrust faults. Here, to assess the impact of fault geometry, we focus on the eastern Himalayan region in the south Central Bhutan, where the topographic frontal thrust (TFT) has been already documented by GPS, paleoseismic, geomorphic and geological studies. This study is based on high-resolution near-surface geophysical investigations, including electrical resistivity, seismic and gravity measurements. Using a similar stochastic inversion approach for all data sets, new quantitative constraints on both fault geometry and petrophysical parameters are obtained to image shallow depths, in the upper ca. 80 m. The combined results from both surface observations and geophysical measurement provide a TFT geometry that is dipping northwards with a shallow angle at the top (0–5 m), steeply dipping in the middle (5–40 m) and flattening at deeper depths (>40 m). Together, our new constraints on the fault geometry allow us to estimate a minimum overthrusting slip rate of 10 ± 2 mm yr−1, which is only a part of the ca. 17 mm yr−1 GPS convergence. This suggests that, in the study area, significant deformation partitioning on several faults including TFT and the Main Boundary Thrust cannot be ruled out. More importantly, assuming constant slip rate, the obtained dip angle variations lead to uplift rate changes with distance to the TFT. This underlines that taking into account uplift rate from terrace dating only at the front location and assuming a constant dip angle fault geometry based on surface observations may significantly bias the slip rate estimates.

Coupled land surface–subsurface hydrogeophysical inverse modeling to estimate soil organic carbon content and explore associated hydrological and thermal dynamics in the Arctic tundra

Abstract. Quantitative characterization of soil organic carbon (OC) content is essential due to its significant impacts on surface–subsurface hydrological–thermal processes and microbial decomposition of OC, which both in turn are important for predicting carbon–climate feedbacks. While such quantification is particularly important in the vulnerable organic-rich Arctic region, it is challenging to achieve due to the general limitations of conventional core sampling and analysis methods, and to the extremely dynamic nature of hydrological–thermal processes associated with annual freeze–thaw events. In this study, we develop and test an inversion scheme that can flexibly use single or multiple datasets – including soil liquid water content, temperature and electrical resistivity tomography (ERT) data – to estimate the vertical distribution of OC content. Our approach relies on the fact that OC content strongly influences soil hydrological–thermal parameters and, therefore, indirectly controls the spatiotemporal dynamics of soil liquid water content, temperature and their correlated electrical resistivity. We employ the Community Land Model to simulate nonisothermal surface–subsurface hydrological dynamics from the bedrock to the top of canopy, with consideration of land surface processes (e.g., solar radiation balance, evapotranspiration, snow accumulation and melting) and ice–liquid water phase transitions. For inversion, we combine a deterministic and an adaptive Markov chain Monte Carlo (MCMC) optimization algorithm to estimate a posteriori distributions of desired model parameters. For hydrological–thermal-to-geophysical variable transformation, the simulated subsurface temperature, liquid water content and ice content are explicitly linked to soil electrical resistivity via petrophysical and geophysical models. We validate the developed scheme using different numerical experiments and evaluate the influence of measurement errors and benefit of joint inversion on the estimation of OC and other parameters. We also quantify the propagation of uncertainty from the estimated parameters to prediction of hydrological–thermal responses. We find that, compared to inversion of single dataset (temperature, liquid water content or apparent resistivity), joint inversion of these datasets significantly reduces parameter uncertainty. We find that the joint inversion approach is able to estimate OC and sand content within the shallow active layer (top 0.3 m of soil) with high reliability. Due to the small variations of temperature and moisture within the shallow permafrost (here at about 0.6 m depth), the approach is unable to estimate OC with confidence. However, if the soil porosity is functionally related to the OC and mineral content, which is often observed in organic-rich Arctic soil, the uncertainty of OC estimate at this depth remarkably decreases. Our study documents the value of the new surface–subsurface, deterministic–stochastic inversion approach, as well as the benefit of including multiple types of data to estimate OC and associated hydrological–thermal dynamics.

A stochastic inversion workflow for monitoring the distribution of CO2 injected into deep saline aquifers

We present a two-step stochastic inversion approach for monitoring the distribution of CO2 injected into deep saline aquifers for the typical scenario of one single injection well and a database comprising a common suite of well logs as well as time-lapse vertical seismic profiling (VSP) data. In the first step, we compute several sets of stochastic models of the elastic properties using conventional sequential Gaussian co-simulations (SGCS) representing the considered reservoir before CO2 injection. All realizations within a set of models are then iteratively combined using a modified gradual deformation algorithm aiming at reducing the mismatch between the observed and simulated VSP data. In the second step, these optimal static models then serve as input for a history matching approach using the same modified gradual deformation algorithm for minimizing the mismatch between the observed and simulated VSP data following the injection of CO2. At each gradual deformation step, the injection and migration of CO2 is simulated and the corresponding seismic traces are computed and compared with the observed ones. The proposed stochastic inversion approach has been tested for a realistic, and arguably particularly challenging, synthetic case study mimicking the geological environment of a potential CO2 injection site in the Cambrian-Ordivician sedimentary sequence of the St. Lawrence platform in Southern Quebec. The results demonstrate that the proposed two-step reservoir characterization approach is capable of adequately resolving and monitoring the distribution of the injected CO2. This finds its expression in optimized models of P- and S-wave velocities, density, and porosity, which, compared to conventional stochastic reservoir models, exhibit a significantly improved structural similarity with regard to the corresponding reference models. The proposed approach is therefore expected to allow for an optimal injection forecast by using a quantitative assimilation of all available data from the appraisal stage of a CO2 injection site.

Large‐scale hydraulic tomography and joint inversion of head and tracer data using the Principal Component Geostatistical Approach (PCGA)

AbstractThe stochastic geostatistical inversion approach is widely used in subsurface inverse problems to estimate unknown parameter fields and corresponding uncertainty from noisy observations. However, the approach requires a large number of forward model runs to determine the Jacobian or sensitivity matrix, thus the computational and storage costs become prohibitive when the number of unknowns, m, and the number of observations, n increase. To overcome this challenge in large‐scale geostatistical inversion, the Principal Component Geostatistical Approach (PCGA) has recently been developed as a “matrix‐free” geostatistical inversion strategy that avoids the direct evaluation of the Jacobian matrix through the principal components (low‐rank approximation) of the prior covariance and the drift matrix with a finite difference approximation. As a result, the proposed method requires about K runs of the forward problem in each iteration independently of m and n, where K is the number of principal components and can be much less than m and n for large‐scale inverse problems. Furthermore, the PCGA is easily adaptable to different forward simulation models and various data types for which the adjoint‐state method may not be implemented suitably. In this paper, we apply the PCGA to representative subsurface inverse problems to illustrate its efficiency and scalability. The low‐rank approximation of the large‐dimensional dense prior covariance matrix is computed through a randomized eigen decomposition. A hydraulic tomography problem in which the number of observations is typically large is investigated first to validate the accuracy of the PCGA compared with the conventional geostatistical approach. Then the method is applied to a large‐scale hydraulic tomography with 3 million unknowns and it is shown that underlying subsurface structures are characterized successfully through an inversion that involves an affordable number of forward simulation runs. Lastly, we present a joint inversion of head and tracer test data using MODFLOW and MT3DMS as coupled black‐box forward simulation solvers. These applications demonstrate the advantages of the PCGA, i.e., the scalability to high‐dimensional inverse problems and the ability to utilize multiple forward models as black boxes.

Fining up the seismic stochastically for reservoir characterization, offshore North West Australia

SUMMARY The vertical resolution of conventional seismic data limits the illumination and definition of thin sand reservoirs offshore North West Australia. With the availability of more wells within a field it is justified to use a well-centric stochastic seismic inversion approach to characterize the reservoirs. The stochastic AVO seismic inversion yields multiple triplet realizations of acoustic impedance, Poisson’s ratio and density at the fine enough vertical sampling for a more relevant reservoir modelling procedure. The paper will discuss the stochastic seismic inversion approach, followed by a structured analysis workflow which condenses the generated multiple realizations to more interpretable data volumes.

Geostatistical Data Integration Model for Contamination Assessment

Soil contamination assessments can be improved with new methods aimed at the accurate estimation of the volume and extension of contaminated soil to be remediated. Geostatistical models that use secondary information to characterize soil contamination are incorporated into a new integration model to provide accurate three-dimensional maps. The proposed integration model is based on a stochastic inversion approach and uses sequential indicator simulation. A two-dimensional reference image representing the areal extension of the contamination is combined with local measurements of contamination in the vertical direction, to render a three-dimensional contamination map. To demonstrate how well the integration model performs, the case study presented focuses on geophysical data and how it can be integrated with soil contamination measurements to improve the characterization of a contaminated site. The results show that the model reproduces successfully the reference image thus providing an accurate three-dimensional contamination map.

A Stochastic Analysis of Cross‐Hole Ground‐Penetrating Radar Zero‐Offset Profiles for Subsurface Characterization

Cross‐hole ground‐penetrating radar (GPR) zero‐offset profiling (ZOP) is a well‐established technique for the measurement of the one‐dimensional (1‐D) vertical distribution of soil dielectric permittivity and has often been used in time‐lapse mode for the monitoring of natural or man‐made infiltration and soil moisture redistribution processes in the vadose zone. However, in spite of its widespread use, the quantitative interpretation of ZOPs in terms of dielectric permittivity profiles is known to be fraught with difficulties. Often a simplified approach is adopted that translates directly ZOP travel times into electromagnetic (EM) velocities and these, in turn, into dielectric permittivities and soil volumetric moisture content. This approach is known to lead to over‐smoothed moisture content profiles, which are a consequence of the ZOP measurement scale, controlled by critical refractions along fast layers and averaging effects within the first Fresnel zone. Such smooth profiles are not necessarily compatible with the true soil moisture content that should be reproduced by water flow models. In this paper we present a stochastic inversion approach that aims at solving these issues. The approach is based on a forward stochastic simulation that generates the expected travel times by means of a two‐dimensional (2‐D) full waveform modeling, thus reproducing all physical processes that contribute to profile smoothing. We first assess the robustness of this approach on synthetic data, where the true dielectric permittivity profile is known. Then we apply the technique to two different case studies, where the results of the proposed technique show supporting evidence from independent information on the sites' stratigraphy. The proposed approach proves to be capable of reconstructing sharp dielectric profiles, in addition to assigning relevant uncertainty bounds derived from the expected errors in travel time picking. These results may prove very useful in a correct assessment of the hydrological implications of the measured dielectric permittivity, and thus moisture content, profile.

Shallow Seismic Velocity Structure of the Betic Cordillera (Southern Spain) from Modelling of Rayleigh Wave Dispersion

A detailed dispersion analysis of Rayleigh waves generated by local earthquakes and occasionally by blasts that occurred in southern Spain, was undertaken to obtain the shear-wave velocity structure of the region at shallow depth. Our database includes seismograms generated by 35 seismic events that were recorded by 15 single-component short-period stations from 1990 to 1995. All these events have focal depths less than 10 km and body-wave magnitudes between 3.0 and 4.0, and they were all recorded at distances between 40 and 300 km from the epicentre. We analysed a total of 90 source-station Rayleigh-wave paths. The collected data were processed by standard digital filtering techniques to obtain Rayleigh-wave group-velocity dispersion measurements. The path-averaged group velocities vary from 1.12 to 2.25 km/s within the 1.0-6.0 s period interval. Then, using a stochastic inversion approach we obtained 1-D shear-wave velocity–depth models across the study area, which were resolved to a depth of circa 5 km. The inverted shear-wave velocities range approximately between 1.0 and 3.8 km/s with a standard deviation range of 0.05–0.16 km/s, and show significant variations from region to region. These results were combined to produce 3-D images via volumetric modelling and data visualization. We present images that show different shear velocity patterns for the Betic Cordillera. Looking at the velocity distribution at various depths and at vertical sections, we discuss of the study area in terms of subsurface structure and S-wave velocity distribution (low velocity channels, basement depth, etc.) at very shallow depths (0–5 km). Our results characterize the region sufficiently and lead to a correlation of shear-wave velocity with the different geological units features.