Heidi Kjønsberg

Seismic data acquired at different times over the same area can provide insight into changes in an oil/gas reservoir. Probabilities for pore fluid will typically change, whereas the lithology remains stable over time. This implies significant correlations across the vintages. We develop a methodology for the Bayesian prediction of joint probabilities for discrete lithology-fluid classes (LFCs) for two vintages, simultaneously considering the seismic amplitude-variation-with-offset data of both vintages. By taking into account the cross-vintage correlations of elastic and seismic properties, the simultaneous inversion ensures that the individual results of both vintages, as well as their differences, are consistent and constrained by the seismic data of both vintages. The method relies on prior geologic knowledge of stratigraphic layering, the possible lithologies and fluids within each layer, and the possible cross-vintage changes in lithology and pore fluid. Multiple LFCs can be used to represent different strengths of dynamic cross-vintage changes. We test the algorithm on a synthetic data set and data from the Edvard Grieg field in the central North Sea. Synthetic results demonstrate that the algorithm is able to use dual-vintage data together with a prior model specifying their correlations to calculate joint LFC posterior probabilities for both vintages with a lower degree of uncertainty than independent single-vintage inversions. The Edvard Grieg results indicate that the underlying model is sufficiently general to explain 4D variations in seismic data using a reasonably simple prior model of 4D LFC changes.

Seismic data acquired at different times can provide insight into changes in an oil/gas reservoir. We apply time-lapse AVO inversion methodology to the Edvard Grieg field in the North Sea. The methodology we use inverts for discrete lithology-fluid classes (LFCs), where different fluid fillings and/or pressure effects in the monitor are represented by different LFCs. Our focus here is on fluid substitution. We only look at amplitude effects, and use monitor time-lapse data that is time aligned with the base seismic.

The GIG, Geophysical Inversion to Geology, consortium is a research consortium run by Norwegian Computing Center, with the aim of developing new understanding, new methods and software for obtaining reservoir properties from geophysical measurements This note gives a summary of the work in 2021 and the suggested work plan for 2022. The final work plan for 2022 will be decided on the annual meeting in Jan 2022.

We have estimated the reservoir sand thickness and internal layering in the Aurora area, a planned geological CO2 storage site in the northern North Sea. The results are obtained by stochastic Markov chain Monte Carlo (MCMC) simulations on probabilistic lithology and fluid distributions in the subsurface. The probabilistic distributions are obtained by inverting the seismic data in a Bayesian framework. The inversion is using angle-dependent pre-stack seismic data and the linearized seismic forward model of Aki and Richards to estimate the posterior probabilities of lithology and fluid classes (facies) in the subsurface. The facies are defined from well log data by identifying depth intervals with distinct elastic responses to seismic waves. The inversion methodology and the MCMC simulations are developed and implemented by the Norwegian Computing Center. The planned CO2 storage reservoir comprises the Johansen and Cook sandstones belonging to the Early Jurassic Dunlin Group. The injection site (well 31/5-7) is located south and west of the Troll field in the northern North Sea and is currently being developed by Equinor in partnership with Total and Shell as part of the Northern Light project. The seismic inversion is mapping structural details like faults and internal layering of the sandstones, and the MCMC simulations estimate the probability of sand thickness and expected layering. Results show that the Johansen Formation sandstone has a tendency of layering towards the west and largest thickness to the east in the inversion area. The sandstone in the Cook Formation is generally thinner, and probability maps indicate that it pinches out to the east. Cumulative thickness distributions provide low (P90), median (P50), and high (P10) thickness maps of both Johansen and Cook Formation sandstones. The presented methodology defines a functional workflow for quantifying the thickness uncertainty and possibly internal layering of the reservoir sands. Such results may provide important input for future field development strategies and CO2 migration predictions.

The GIG, Geophysical Inversion to Geology, consortium is a research consortium run by Norwegian Computing Center, with the aim of developing new understanding, new methods and software for obtaining reservoir properties from geophysical measurements This note gives a status of the GIG consortium after its fourth year of operation, and also gives the suggested work plan for 2021.

Lithology and fluid prediction from seismic data is traditionally done in two steps; first an inversion of the seismic data to elastic parameters, and subsequently a prediction of lithology and fluid based on a rock physics model linking the elastic parameters to individual lithology and fluid combinations. Recently, a number of inversion algorithms have been developed that, based on Bayesian statistical methodology, estimate the probability of lithology and fluid directly from seismic data. In this paper we compare the performance of two state-of-the-art Bayesian inversion algorithms on a real data set from the Volund field in the North Sea. The first algorithm follows the traditional two-step approach and cannot take into account the stratigraphic ordering of lithology and fluid. The second algorithm, referred to as one-step, evaluates possible lithology and fluid combinations within a vertical window around each inversion point enabling correct stratigraphic ordering. We find that the one-step inversion resolves more details and honours the data more strongly than the two-step approach. The latter is more prone to return the prior model if information in the seismic data is not sufficiently strong. Both models detect hydrocarbon filled sand injectites that are typical for the field.

A multigrid Markov mesh model for geological facies is formulated by defining a hierarchy of nested grids and defining a Markov mesh model for each of these grids. The facies probabilities in the Markov mesh models are formulated as generalized linear models that combine functions of the grid values in a sequential neighborhood. The parameters in the generalized linear model for each grid are estimated from the training image. During simulation, the coarse patterns are first laid out, and by simulating increasingly finer grids we are able to recreate patterns at different scales. The method is applied to several tests cases and results are compared to the training image and the results of a commercially available snesim algorithm. In each test case, simulation results are compared qualitatively by visual inspection, and quantitatively by using volume fractions, and an upscaled permeability tensor. When compared to the training image, the method produces results that only have a few percent deviation from the values of the training image. When compared with the snesim algorithm the results in general have the same quality. The largest computational cost in the multigrid Markov mesh is the estimation of model parameters from the training image. This is of comparable CPU time to that of creating one snesim realization. The simulation of one realization is typically ten times faster than the estimation.

Traditionally fault seal calculations take place directly within the simulation grid. This approach works well for grids where all the faults are aligned along the grid pillars, but implementing an algorithm that works with stair-stepped representation of the faults has proven to be very difficult. Especially the calculation of the displacement field used both indirectly in the fault seal parameter calculation and directly in the calculation of fault zone permeability is challenging. It is hard to find where the different grid layers intersect the fault trace, and the layers are not always completely represented on both sides of the fault. We present a novel algorithm where the calculation of the fault zone permeability is carried out on a 2D plane representing the fault surface. The input parameters needed for calculating the fault zone permeability are resampled from the simulation grid onto the 2D plane, while the resulting fault zone permeability is resampled back into the simulation grid, prior to calculation of the fault transmissibility. The new approach is shown to generate good results both for pillar-faulted grids, and for grids with stairstepped faults, and also works well near complex truncations.

Traditionally fault seal calculations take place directly within the simulation grid. This approach works well for grids where all the faults are aligned along the grid pillars, but implementing an algorithm that works with stair-stepped representation of the faults has proven to be very difficult. Especially the calculation of the displacement field used both indirectly in the fault seal parameter calculation and directly in the calculation of fault zone permeability is challenging. It is hard to find where the different grid layers intersect the fault trace, and the layers are not always completely represented on both sides of the fault. We present a novel algorithm where the calculation of the fault zone permeability is carried out on a 2D plane representing the fault surface. The input parameters needed for calculating the fault zone permeability are resampled from the simulation grid onto the 2D plane, while the resulting fault zone permeability is resampled back into the simulation grid, prior to calculation of the fault transmissibility. The new approach is shown to generate good results both for pillar-faulted grids, and for grids with stair-stepped faults, and also works well near complex truncations.

Traditionally fault seal calculations take place directly within the simulation grid. This approach works well for grids where all the faults are aligned along the grid pillars, but implementing an algorithm that works with stair-stepped representation of the faults has proven to be very difficult. Especially the calculation of the displacement field used both indirectly in the fault seal parameter calculation and directly in the calculation of fault zone permeability is challenging. It is hard to find where the different grid layers intersect the fault trace, and the layers are not always completely represented on both sides of the fault. We present a novel algorithm where the calculation of the fault zone permeability is carried out on a 2D plane representing the fault surface. The input parameters needed for calculating the fault zone permeability are resampled from the simulation grid onto the 2D plane, while the resulting fault zone permeability is resampled back into the simulation grid, prior to calculation of the fault transmissibility. The new approach is shown to generate good results both for pillar-faulted grids, and for grids with stair-stepped faults, and also works well near complex truncations.

This note describes how to find a prior 4D model for seismic parameters from an underlying rock physics model. The main ideas are expressed in terms of the Markov property for the dependency between different time steps, time locality in the sense that the seismic parameters at any time instance are found from the rock physics parameters as they are at that same time instance, and the use of Gaussian-linear models.

This note discusses the general framework of 4D inversion with multiple data-sources, phrased as a filtering and smoothing problem. Statistically this is formulated as a Markov process prior distribution, with multi-dimensional state variables. This is a well-established field, and the purpose of the current paper is twofold; to recapture the general equations and relations of the filtering and smoothing theory for Markov processes; and to outline the general framework we will work under in the NFR funded MonCO2 project.

The aim of this project is to develop new and improved methods for modelling geological facies by combining the efficiency of multipoint methods with the robustness and consistency of Markov random field methods. It is a goal to develop software tools and test new methods on real cases and data. The project started in 2006 and is planned to finish in December 2008. The project is a mutually beneficial collaboration between the Norwegian Computing Center, the Norwegian University of Science and Technology and Stanford University. It is sponsored by the Research Council of Norway, ENI and StatoilHydro. The industry partners have contributed with valuable input for the discussions as well as useful guidelines for the choice of research topics. In this document we first provide a short summary of the main achievements in the project. This is followed by an introduction to multiple point methods and a somewhat more detailed account of each theme of the project.

We consider the problem of formulating Markov Mesh models as Markov Random field. Since Markov Mesh models are a subclass of Markov Random Fields, this can in principle always be done. In these notes we explore the details of the parameter translation for the case of a 1D stationary Markov Mesh model consisting of a homogeneous external field and general two-particle interactions. We show theoretically that the corresponding Markov Random Field includes also higher order interactions, the complexity of the interactions being limited by the neighbourhood structure. Explicite formulas and recursive algorithms expressing the MRF parameters in terms of the independent Markov Mesh parameters are provided. A matlab implementation of the parameter translation is also described. Using simulations we explore to which extent the higher order interactions in the MRF formulation are necessary in order to reproduce the statistics of the considered 2nd order Markov Mesh models. We find that in general the higher order interactions have significant impact.

The main purpose of these notes is to provide documentation of version 0.1 of the Markov random field implementation done for the Multipoint project. The report is mainly meant for internal use in the project, serving as a user manual for running the implemented program and as a guide for programmers wanting to modify it.