An effective field scale reservoir model for history matching and reservoir prediction. Part 1. Model architecture and specification

Madatov A.G.

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

УДК [622.276 : 622.031] : 004.9(045)

An effective field scale reservoir model for history matching

and reservoir prediction

Part 1. Model architecture and specification

A.G. Madatov

Polytechnic Faculty of Murmansk State Technical University, Higher Mathematics Chair

Abstract. An effective reservoir model concept has been introduced to ensure a practical solution to the history matching problem for a multi-well data set. A deterministic process-based approach to its architecture and specification has been implemented. The common hydraulic element concept is combined with generalised empirically grounded permeability functions for mixed clastic lithology to get a 3D grid containing a reservoir transport property. The clay content constant is used as a universally sensitive control model parameter for specification of absolute and relative permeability distribution for both laterally continuous and discontinuous cases. The power trend functions are used for approximation of inherent lateral distribution of clay content in field scale. Being in close relation with relevant depositional environment they are considered as 'soft' model constraints. The relative sensitivity of model vector components is analysed based on plausible multi-well model examples.

Аннотация. Концепция эффективной модели залежи углеводородов вводится для обеспечения практического решения задачи инверсии промысловых данных, накопленных по системе добывающих скважин в процессе эксплуатации месторождения. Используется детерминистический, депозиционный подход к построению модели и ее спецификации. Концепция обобщенного потокового элемента комбинируется с генерализованным представлением об эмпирических функциях проницаемости для терригенных осадочных пород, слагающих коллектор, для получения трехмерной сеточной модели транспортных свойств нефтяного пласта. Г линистость используется как универсальный чувствительный параметр модели при спецификации свойств абсолютной и относительной проницаемости. Рассматриваются непрерывный и разрывный (тектонически блокированный) варианты представления модели пласта. Степенные функции используются для описания латерального тренда глинистости пласта, присущего большинству условий осадконакопления. Свойства трендов выступают в роли априорных ограничений, накладываемых на подбираемую эффективную модель. Относительная чувствительность компонент модельного вектора параметров анализируется на реалистичных примерах разработки нефтяного месторождения.

Ключевые слова: спецификация модели, проницаемость, пористость, инверсия данных, анализ чувствительности Key words: model specification; permeability; porosity; data inversion; sensitivity analysis

1. Preface

Management of oil-gas field production is traditionally based on construction and maintenance of the relevant field scale Reservoir Model (RM) suitable for simulation of a field performance (Schelkachev, 2002; Zeltov, 1998).

Historically more attention was paid to the generation of accurate and sophisticated solutions of the forward reservoir modelling problem, whereas the issue of matching the model parameters to the available dynamic data had secondary priority. Until recently the history matching problem needed advanced reservoir engineering competence. It was resolved by expert-based manual adjustment of the control RM parameters that were least certain and most sensitive to the result. In the last decade the problem of history matching in reservoir engineering or establishment of a data-constrained reservoir model was re-addressed throughout the process from human expertise to the data inversion problem (Christie et al., 2002). This is because of sharp jump in data amount and quality (time-lapse seismic) and also due to great advance in computer power in particularly in connection with the inverse problem solutions (Tarantola, 1987). Additional impetus in this data inversion direction was generated by requirements of real time reservoir management such as "Closed loop strategy" (Sarma et al., 2005).

Thus present day viewing of an oil/gas field production management highlights development of a new generation of reservoir models capable of fast and adequate adjustments to the multi source static and dynamic data (Aziz, 1989). Namely: the real time management of reservoir performance requires new approaches to the discrete representation of a reservoir model, its parameterisation, upscaling and calibration in order to integrate multi-discipline input and facilitate monitoring for reservoir engineers.

This can be considered as a main motivation of the research study, results of which are presented here.

457

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

2. Introduction

Successful reservoir management within the range of a single field can be understood in terms of optimal production planning (Calvert, 2005) as following:

The general goal of reservoir engineering is harmonisation between collected knowledge-data and current reservoir model conception to provide sufficiently good input to the field drilling and production programs that ensures most effective extraction of the remaining reserves.

This goal could be reached in three different ways:

1. Optimisation of well placement for both injector and producer.

2. Maximising the production potential (production - injector rates).

3. Minimising remaining reserves, including hidden or poorly recoverable reserves.

Generally some production cost function can be constructed as a combination of the above-listed criteria and fixed as a specific Target Performance Criteria (TPC). Regardless to this criterion, there is some sequence of interconnected data processes, which can be represented as a general workflow diagram (see fig. 1).

According to this scheme, the managing process begins from initialization of a Reservoir Model (RM) in terms of its seismically visible sub-surface structural elements (fixed framework part) and textural seismically invisible sub-surface elements (tuneable textural part). Depending on static well data (wire logs, core measurements) and 3D seismic available, the initial oil/gas-saturation, geopressure and porosity distributions within the field scale could be referred to the fixed or tuneable parts of the RM.

Based on this setting, the detailed knowledge about reservoir structure and assumptions about possible texture variations from the set of available data then can be represented in a finely gridded 3D cellular model.

The architecture of such multiscale and multidisciplinary models then should be adjusted for data fitting purposes (upscaled). The RM adjustment in this context implies simplification of the RM architecture, optimization of grids and reduction of the number of tuneable model parameters.

Properly upscaled RM is then matched to the available production data records. Note that the level of RM details remaining after upscaling for subsequent data inversion depends greatly on the stage of production history and content of production history record.

The set of matched reservoir models then is entered into the RM validity analysis routine, where calibrated models are ranked in agreement with the fitness to the Target Performance Criteria.

Both data fitting and model ranking processes are represented by recursive loops with active user interface.

Finally, the reservoir management system keeps an upper rank calibrated reservoir model updatable in accordance with updating of the production history record. The calibrated RM contains all necessary input for an industry available reservoir simulator to produce synthetic dynamic data. As the model framework is consistent with available static data and the core textural part of the model is fitted against available dynamic data, the relevant synthetic production history ahead of the present time could serve as a prediction scenario.

This scheme seems to be common for reservoir management regardless of the stage of the production history and class of the model used for history matching (Talukdar, Brusdal, 2005; Tipping et al., 2006). At this the RM architecture, the way used for its specification, simulation, upscaling and fitting to the real production data should be adjusted to the current stage of production history, amount of available data and purposes of RM identification.

Fig. 1. Upper level workflow diagram of reservoir management system

458

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Note that this scheme essentially contains several reservoir models: initial, upscaled and matched to the production history. The comprehensiveness of every model should evidently be quite different. The initial RM can be considered as a static reference one, where all available static data are accommodated. The RM that is simplified during upscaling is assumed to keep level of details suitable for fitting, which should be well balanced with the stage of production history. The history matched (calibrated) RM can be considered as a set of updatable versions of the upscaled RM. This is the key element of the scheme since it is the crossing point of two main processes: history matching and uncertainty analysis.

Let us distinguish three consequent and unavoidable stages of oil/gas production history (Krylov et al., 1948; Ter-Sarkisov, 1999; Kanevsaya, 2002):

• initial stage, when a field scale reservoir model is being built and oil/gas reserves in place being defined;

• active production stage, during which most of the oil/gas reserves are being extracted;

• mature stage, when long production history record is collected, oil/gas production rate becomes unstable and water flooding risk is high.

This subdivision is of course rather conditional. The duration of each stage very much depends on specifics of the given field subsurface architecture and production-injection well system implemented for the reservoir engineering. Still, it is certainly true that the level of knowledge collected about reservoir increases directly proportional to the depletion of oil/gas reserves. Thus, the comprehensiveness of the reservoir model should be increasing consistently in line with the amount and quality of collected dynamic well data, which are increasing during production history.

In practice, however the logic "the more detailed reservoir model the better the management" remains as the most popular, despite the fact that the model details and ways of specification used for data simulation and data inversion are quite different. Adherence to this logic often leads to generation of an overcomplicated reservoir model having billions of calls to be specified and thousands of model parameters to be tuned. Evidently, history matching as an inverse problem by definition cannot be successfully resolved based on such a heavy model.

We claim that the reservoir model that is suitable for prediction of reservoir performance and updatable ahead of the current edge of production history should belong to an effective model class. Here the level of details and dimension of relevant forward problem are in close balance with the amount and quality of data available for model identification, which have in our particular case a close link with the stage of production history.

We propose a concept of an effective reservoir model generation - upscaling and calibration for quick and optimal reservoir management. It does not imply any new forward problem solutions, but is based on implementation of industry available reservoir simulators for production data inversion.

The main efforts here are focused on preparing of a properly organized RM (model specification), adopting this model to the data fitting process and producing the history matched reservoir model set. Below we will discuss the achieved solutions in the context of industry available approaches and strategies. This allows better recognition of new features of the presented concept.

Despite the fact that a data inversion definition and relevant requirements to the model parameter space and forward modelling operator served as guide lines through all of the stages of the RM creating from initialization of a reference model up to its parametrical upscaling and calibration, we will present the results in their natural order: part 1 - model architecture and specification; then part 2 - its parametrical upscaling and calibration. Both parts are presented separately only for reasons of better paper composition and easier understanding by the reader. Essentially they are considered by the author as inseparable components of the entire concept.

3. Earth model architectures

Understanding of reservoir performance during production life depends on links established between RM elements and observed field dynamics. The concept behind this link could be treated generally as a reservoir model architecture.

As the production life becomes longer and relevant dynamic data set includes more and more real information, the relevant RM architecture is evolving. In particular, a pure formal RM architecture like Response Surface Model (Carrears, Turner, 2006) appropriate for initial production stage is transformed into the more and more detailed earth model architecture (Donselaar, Geel, 2003; Joseph et al., 2004) during later production stage.

The formal model parameters are normally associated with sensitive factor descriptors (Joreskog et al., 1976) or coefficients used in relevant surrogate equations, regressors in case of linear surrogate equations (Myers, Montgomery, 1995). Such kind of model specification can be considered as built for particular fitting purposes. The relevant formal RM architectures are empirically driven, where there is no room for prior model constraints and a posterior interpretations connected with a local geology setting. The advantages of the formal RM architectures are related first of all to fast and unique model identification. This feature is quite convenient at the initial stage of production history, when the production record is short and represented mainly by the set of

459

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

production indicators (Cheong, Gupta, 2005). Applicability of this approach to the reservoir understanding and the prediction of a later stage of production history gradually decreases (this topic will be considered in more details in part 2).

In contrast, a reservoir model in the later stages of production history contains more and more details of the earth model architecture, where flow units and flow barriers become its inherent elements. Both these elements associated with non-formal (later stage) RM can be described at the seismic or sub-seismic scale depending on details required for reproducing reservoir heterogeneity.

It is convenient to distinguish between framework (structural) and texture (core) parts of any geologically grounded reservoir model. Let us define these two constitutive components of effective reservoir model architecture as follows:

A framework component is the a priori defined part of RM, which describes the structure of reservoir subsurface geometry including its tectonic features according to seismic scale resolution.

A textural component is a priori unfixed part of RM, which describes the reservoir heterogeneity at sub-seismic scale in the form of multi-layered and/or gradient litho-facies system.

Note, that both framework and textural components could also be treated as deterministic and stochastic components of RM architecture respectively.

The main advantage of such kind geologically grounded RM is the ability for reservoir engineers to account for a priori available geology knowledge, which could include some input from litho-facies analysis, sequence stratigraphy, simulation of sedimentation processes etc (Miall, 2000; Franseen et al., 1991; Merriam, Davis, 2001).

Traditionally the workflow of an earth model setting includes the following three sequential steps (Baltrusch et al., 2004):

1. reconstruction of the seismic scale geometry features, during which the seismic interpreter picks tops and faults to generate the framework part of reservoir model;

2. identification of the flow unit and flow barrier elements in agreement with hypothesis about RM texture during which reservoir geologist produces reservoir layering controlled by relevant depositional system hierarchy;

3. aggregation of framework and textural components in form of cellular earth model.

Practically this process implies multiple looping and reconstruction of every sub-stage based on currently obtained results in agreement with the Shared Earth Model concept (Collazos et al., 2003).

For example, inter-well correlation can first be done formally based on Kriging variance concept (Wong, Boerner, 2004), or Neural Network method (Aminzadeh et al., 2000), or combined 3D seismic - well-log correlation principles (Sandjivy et al., 2004). After reinterpretation of achieved 3D formal attributes' model against available sedimentological concept (Miall, 2000) or process based depositional model (Franseen et al., 1991; Merriam, Davis, 2001) iterative refinements of textural component could be accomplished (Al-Deeb et al., 2002).

Below each step of the earth model setting is considered in more detail.

Fig. 2. Formal stream-tube based (A) and non-formal geology driven (B) earth model architecture established for

the same speculative reservoir model

a - the gray mesh represents stream-tube position with stream barrier indicated by dotted lines; b - the flow units and barriers given by stacked discontinuous sub-surfaces with facies lithology variation (clay content) attached in form of contour lines. Red sticks indicate well positions which are equivalent for (a) and (b).

460

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

3.1. Reconstruction of the seismic scale geometry features

This part of earth model setting includes all stages of standard 3D seismic data processing aimed at full scale reconstruction of seismically visible litho-stratigraphy/ tectonic features of subsurface (Moushin et al., 1990). Major factors affecting uncertainty in a structural framework are essentially also those that control unique and accurate TWT-Depth conversion of seismic horizons and fault recognition. Namely: inherent seismic scale resolution; finite density of an observation system; signal/noise ratio; acoustic contrasts and so on (Averbukh, 1982). An additional component of the framework uncertainty is inevitably introduced by human interpretation of the seismic wave field at all stages of processing.

Still the seismic data within the field area normally have good quality, dense coverage and firm link with well system. Thus, a joint interpretation of 3D migrated seismic data versus petrophysical network ensures a framework model with sufficient accuracy for reservoir simulation purposes (Calvert, 2005).

3.2. Identification of the Flow Unit / Flow Barrier elements

The earth model of reservoir at field scale includes at least 3 hierarchically embedded structure levels corresponding to 4-6 orders of cycle associations represented in Table 1 (Liu et al., 2002).

The levels of 4th and higher order of cyclicity (accumulation time <100,000years) represent level of fluvial channel-belt deposits; delta complex and alluvial fans, shelf sand-ridge fields, major lobes of submarine fans, etc. In most of the cases these textural elements are recognizable from well data analysis. But the inter-well correlation of relevant boundaries by using seismic data is very uncertain (Bouvier et al., 1989; Eide et al., 2002). This is because such high levels of sediment aggregation into the thin layer/ irregular channel system generally imply a non-steady state sedimentation rate model for a deposition process. The majority of bounding surfaces within this group is much less regular than at entire reservoir level since they might be caused by autogenic mechanisms: channel of migration; deposition controlled lithology substitution etc. It is important also to stress that the coverage of the well data is rather limited and irregular.

The common approach at this stage of earth model establishment is usage of pure stochastic model concept (Donselaar, Geel, 2003; Sandjivy et al., 2004; Raisson, Temple, 2004) or sedimentologically controlled model concept (Franseen et al., 1991; Merriam, Davis, 2001).

Both options include recognition of main irreducible architectural elements in form of "elementary flow channel" and/or "elementary flow barrier", specifying their geometry characteristic (length, width, thickness, shape, margin collapse, etc.) and generation of an envelope of architecture scenarios (sedimentological images) with the given probabilities predefined for specified characteristics (Labourdette et al., 2003).

At this the geometry characterisation of architectural elements in the purely stochastic approach follows the normal or log-normal distributions with Probability Distribution Functions (PDF) deduced from geostatistical variogram analysis for relevant spatial variabilities.

Table 1. The generalised hierarchies of architectural units after (Miall, 2000)

Order of cycle Aver. Thick ness [m] Tectonic process Stratigraphy signature Sediment rate [m/My] Litho- stratigraphic units Tectono-stratigraphic cycle

1-2 102- 103 Regional flexural loading, imbricate stacking Basin fill complex or depo-system 100-1 Group Global super-continent cycle

3 101 -103 Effects of basement heterogeneities during crustal shortening Para-sequence complex 101-3 Formation (multi- layer association) Eustatic sea level oscillation and regional basement movement

4 1 О <N О О Fault-propagation anticline and foreland syncline Channel belt, shelf sheet 102-4 Member (Single layer) Local basement movement & changes in intraplate stress regime

5 10-1 -101 Thrust overstep branches developing inside fault-propagation anticline Macro-form, Fan, Facies package, Major-to-minor lobe 102-5 Facies Major orbital forcing cycles in glacio-eustasy, productivity etc.

6 10-3 -10-1 Movement of individual thrust plates, normal listric faults, minor folds Mesoform, sequence, single facie 103-6 Lamina Minor orbital forcing cycles in glacio-eustasy, productivity

461

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

In the sedimentological approaches the geometry characterisation is controlled by a combination of depositional processes, which produce distinct heterogeneous sediment sequences with unique scaling characteristics (Miall, 2000). A hierarchy of scales for a number of the depositional systems including fluvial-deltaic, coastal-estuarine, shelf, turbidite and eolian have been proposed here (Van Wagoner et al., 1990).

3.3. Aggregation of deterministic and stochastic components

The fine-scale reservoir heterogeneities have a non-random nature but are closely controlled by facies type of the considered field. Thus correct choice of the litho-facies origin of the given reservoir could be considered as an additional (sedimentological) component to the structural component of relevant earth model architecture. Moreover, it could be used as a guideline for population of static petrophysical data at the RM initialization (see Appendix 1 for more details).

The stochastically driven texture component of earth model architecture can be grounded on statistical integration of similar geo-conditions or process-oriented deposition modelling (Fractional Kriging, Sequential Gauss Simulation; Sedimentary Bedding (SBED) method (Wen et al., 1998) etc.) capable of producing textural templates for established and fixed framework.

In general, the multi-level RM architecture design could be summarised as following hierarchically built consequence of embedded models:

1. micro level: appropriate rock model, which provides a process-dependent pore-(fracture)-scale flow parameter set through Network methods (0ren, Bakke, 2003; Lerdahl et al., 2005);

2. local level: litho-facies scale model, which comprise several rock models, and flow parameters in common flow unit parameterisation (Al-Deeb et al., 2002);

3. macro level: fine scale geological model which incorporates seismic scale input (Framework part) with sub-seismic scale information (Textural part) in a common fine scale earth model to be built on 106-9 cells (Labourdette et al., 2002).

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

The first two levels of RM architecture define its textural component, which belongs to the sub-seismic resolution interval of geological space. They are described stochastically. These levels of RM architecture are obviously justified by simulation of flow phenomena in near well locality.

Lower level incorporates both deterministic (framework) and stochastic (texture) components of RM architecture. This level of RM architecture is appropriate at simulation of flow phenomena on filed block - entire field scale.

4. Initialization of a non-formal reservoir model

Getting the framework and textural components of RM architecture can be treated as an obligatory preprocessing part of the initialization of a non-formal reservoir model. In order to prepare cellular earth model for reservoir simulation one should specify rock matrix properties (porosity, permeability, capillary barriers) and pore fluid properties (initial oil/gas concentration) in agreement with established RM architecture. The PVT fluid properties (density, viscosity, compressibility, wettability) are normally described as physical constants independent from RM architecture. The relative permeability and fault transmissibility multipliers could be specified explicitly in table form or implicitly by using empirical relationships (see below in this paper). In addition, a production-injection well system should be specified in terms of position, schedule and regime control (Aziz, Settari, 1979; King, Datta-Gupta, 1998) (see also Appendix 2).

The specification of pore fluid properties normally does not create serious problems for geologist, petrophysicist and reservoir engineers. Indeed, the 3D distribution of initial oil/gas saturation normally comes from water-oil/gas-water/gas-oil contact positioning (Dake, 1978). The PVT fluid properties are easily obtained from lab measurements (Denesh, 1998). Also capillary barriers could be estimated based on lab measurements (mercury injection tests) and linked to the reservoir lithology type without big risk for simulation accuracy (Calhoun Jr., 1982).

In contrast, the specification of the rock matrix properties inside a 3D cellular reservoir model requires resolving of scale transfer and inter-well correlation problems.

The strategy here could be based on direct transformation algorithms, or on inter well population algorithms, or on their combination. Significant progress in the 3D simulation of sedimentation processes (Wen et al., 1998; Merriam, Davis, 2001; Labourdette et al., 2003) allows also usage of process based algorithms, where synthetic prototypes of relevant petrophysical data are populated in agreement with predefined reservoir texture (Norris, Lewis, 1991; Ringrose et al., 2003; Lopez et al., 2004).

Note that the in situ measurement of absolute permeability of porous rock is universally recognized as a one order more complex petrophysical problem than that for porosity. The main reason for this is the high and nonlinear scale dependence of any permeability measurements and extremely low rate of background filtration process (Bear, Bachmat, 1991). It could be estimated directly only based on core measurements, or indirectly

462

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

based on Repeated Formation Tests (RFT) or Sweep Rate Tests (SRT). At this some scale transfer problems exist, which could be summarised as follows:

1. Permeability core measurements are representative at pore scale [mc-mm] but not at facies scale [cm-m] and even less at field scale [km] (Nordahl et al., 2003).

2. Both RFT and SRT data represents draining volume scale [m] at near well locality (Idrobo et al., 2000), which is often damaged by drilling process and generally cannot be representative for in situ reservoir heterogeneity.

Thus, it is appropriate to talk about permeability prognosis (not estimation) based on petrophysical data and litho-facies input in the considered context. Note, that usage of fuzzy logic (Cuddy, 2000) or artificial neural network methods (Bhatt, Helle, 2002) for predicting of absolute permeability based on petrophysical data so far has restricted validity. The most common way consists in usage of porosity data as an input for further transformation into permeability scale by using empirically derived relationships. The available approaches are considered in more detail in section 6. Below we consider approaches to the specification of the rock matrix properties in connection with porosity.

The direct transformation algorithms normally imply transforming of seismic velocity data (Berryman et al., 2002), or another porosity-sensitive seismic attributes (Kalkomey, 1996; Knoth, Haller, 2004) available in form of data cube, into a porosity cube associated with an RM grid. Note that the above-mentioned scale transfer problem exists also here because of inherent low resolution of seismic method, that could be considered as only representative at a field scale (102-3 m), whereas reservoir heterogeneity requires 1-2 orders less scale (Averbukh, 1982; Furre et al., 2003). In addition, a seismically driven response expected from theory, even for shallow reservoir, appeared to be too weak and mixed with all kinds of geophysical artifacts (Madatov, 2005). Thus the reliable conversion of seismic attributes into porosity scales often becomes unfeasible (Calvert, 2005).

The conversions of petrophysical data and in situ core measurements into porosity scale are generally free from the scale transfer problem. This approach is commonly recognized as most reliable for RM initialization purposes. In particular, transformation of ultra-sonic, density or neutron logs into porosity scale can be considered as most common and stable (see Appendix 1 for more details). In connection with this kind of petrophysical transformation results for population of porosity values inside a cellular 3D RM it is important to remember the inter-well correlation problem caused by very sparse and irregular well grid within reservoir volume involved in model specification. According to (Alsos et al., 2002) the petrophysical data typically represents just a tiny part (about 10-13 ) of the total reservoir volume. Although the coverage of a field area with good quality 3D seismic data could be several orders higher, it does not help in direct inter-well correlation of micro-to-local level reservoir heterogeneity (sub-seismic heterogeneity). It was shown that up to 70 % of the vertical sand layering is below of the effective seismic resolution and thus these small sand layers cannot be resolved by any kind of inversion (Sandjivy et al., 2004). So there is a great deal of uncertainty about the spatial distribution of the sub-seismic scale features like fault-fracture network, thin layering; pillars-channels heterogeneity and so on, all of which influence on the fluid flow. Thus, a robust population algorithm should provide interpolation of log-derived porosity to all grid cells within a fine scale RM whilst preserving the desired heterogeneity level.

We distinguish between two groups of existing approaches: stochastic and hybrid (seismically constrained). Note that all of them are closely related to the chosen earth model architecture. In other words they are litho-facies controlled.

Stochastic approaches such as kriging, multiple regression (Wong, Boerner, 2004; Labourdette, 2003); neural networks (Bhatt, Helle, 2002), fractals, sequential Gauss simulation (Raisson, Temple, 2004), Karhunen-Loeve expansion (Reynolds et al., 1996); Space renormalization (Christie et al., 1995) are based on statistically driven distribution laws, variograms and correlation information deduced, for each individual facies association, from well logs and core data. These methods provide background for multiple RM simulation and risk assessment. Attachment of any given population result to one particular class of stochastic model ("hydrology stochastic") assumes normal or lognormal probability distribution function and short-range correlation functions. Unfortunately, such hypothesis about variability and uncertainty distribution normally cannot be justified because of insufficient amount of data. Additionally, extremely low and high porosity values most likely are disconnected due to requirement of maximal spatial entropy (Hristopulos, Christakos, 1999).

Hybrid approaches imply implementation of seismically driven transformations, likewise "seismic porosity" or sedimentologically driven simulations as an external constraining/regularization of inter-well interpolation routine. Here co-kriging algorithms (Raisson, Temple, 2004) and Bayesian updating procedures (Christie et al., 2002) could be considered as most popular. Despite a large amount of expensive data involved they are not free from the above-mentioned drawbacks.

463

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

5. Approach to the setting of an effective reservoir model

The setting of non-formal reservoir model includes generation and specification of a cellular earth model in agreement with chosen architecture and its initialization, during which the input data necessary for reservoir simulation should be specified. The approaches to the earth model design and RM initialization that are reviewed above were developed mostly in the context of demands for accurate forward modelling results, where computing costs were presumed to be of secondary importance.

The complexity of such models is often too high to be applicable for history matching in reservoir engineering even at the latest stage of production history. This is mainly because demands to a reservoir model complexity coming from a data inverting view point are ignoring. In certain degree these demands are implicitly satisfying by implementation of different upscaling procedures, which are aimed for simplification of similar none purpose built models in terms of grid architecture and specification rules. In part 2 of this paper we will discuss available approaches to the upscaling problem in more detail (see also (Ates et al., 2003; Boassen, Palatnik, 2005; Lohne et al., 2004; Virnovsky et al., 2004). Here we just stress that most standard upscaling approaches are aimed for optimal grid coarsening. In context of stochastic models the term "upscaling" is used to denote procedures by means of which the large-scale flow and transport properties can be inferred from the statistical information characterizing the medium at smaller scales. Here some renormalization methods are used for RM simplification (Hristopulos, Christakos, 1999).

In effective model concept we proceed based on the requirements of data inversion theory from the very first stage of reservoir model creation, i.e. from designing of relevant earth model architecture and model parameter specification (Madatov, Sereda, 2005)1. These requirements briefly can be reduced toward two main features of an effective reservoir model. Namely:

1. The control model parameters used for specification of transport properties must be highly sensitive to the production history simulation and also physically and/or geologically interpretable. The last is important for ensuring of model constraints and regularization at data inversion process (Madatov, Sereda, 2000).

2. The control model parameter space X where relevant forward modelling operator M[x] is defined must have minimal dimension necessary for reproducing of production data prototype with predefined accuracy. This is important from both inverse problem solvability and computer costs points of view (Williams, Madatov, 2007).

We assume here that the structural framework of a reservoir model should be defined and fixed during joint seismic data and wireline log data interpretation. Uncertain components of this interpretation like 3D fault geometry and reservoir thickness could be finely tuned during identification of early stage reservoir model together with defining of initial geopressure regime and oil/gas saturation field (OWC/GWC positions).

We also assume that 3D porosity distribution could be built during RM initialization based on one of the available approaches reviewed above (see also Appendix 1).

Then choice of interpretable Control Model Parameters (CMP) and minimization of their total amount entirely depend on approach to absolute and relative permeability specification. We implement for this purpose a generalized permeability function approach. Here the permeability scalar/vector is defined for every model cell based on its porosity-litho-facies specification and predefined type of functional link (Corbett, Ellabad, 2003; Ringrose et al., 2003). Clearly, the simpler functional link and the fewer total amount of litho-facies elements (flow units) introduced for the given RM texture are, the shorter the dimension of the relevant model parameter space X.

Specification of litho-facies elements to be included into RM texture at a fine scale level of the relevant earth model is a matter of well data interpretation and sedimentological upscaling (Liu et al., 2002). Here a hierarchy of scales for a number of depositional systems including fluvial-deltaic, coastalestuarine, shelf, turbidite, and eolian has been proposed in close relation to cycle ordering (see also Table 1). Depositional units at each size scale are formed in response to processes occurring over a particular time scale and are physically separable from each other by hierarchy of internal bounding surfaces. The temporal and spatial scales of each representative of deposystem are closely related to its level in system hierarchy.

The conventional routine consists of ranking these surfaces and assigning different significance to them. During this process the boundaries of low level depositional hierarchies (high cycle order in Table 1) may be ignored after a few folds of iterations. But the boundaries of upper level in hierarchy may need to be retained within a RM texture (Van Wagoner et al., 1990).

1 Formal definition of a specific inverse problem for the production data and an effective reservoir model is available in part 2 of this article. Here and below we will use the relevant notification to be consistent. Namely: X defines multi dimensional vector space of control model parameters, where any of possible CMP vectors x and relevant forward modelling operator M[x] are defined.

464

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

We suggest to use this approach as a guide line at the inter-well correlation of flow units (Ellis, 1987). However, implementation of heterogeneity index and the fractional Levy model for optimizing of RM texture complexity is seemed to be part of traditional upscaling technique aimed for simplification of forward modelling (Liu et al., 2002). The parametrical upscaling of fine scale earth model against production data is considered to be more suitable for aimed data inversion purposes. This approach will be considered in part 2 of this paper.

Here we focus on the choice of interpretable model parameters and type of generalized permeability functions operating with them.

6. Effective clay controlled permeability models for clastic reservoirs

It is important to use for reservoir model specification such a model parameter set that could control its transport property for a wide range of litho-facies conditions. It is advisable that such parameter set will be clearly interpretable in order to limit in advance its possible range of variation.

We consider Rock Clay Content (RCC) as a perfect candidate to be key model parameter for clastic reservoirs in this context. First, because it has a clear link with various depositional settings; second, because RCC is a dimensionless characteristic of rock varying from 0.0 for pure sands up to 1.0 for pure mud-rocks (clay or shale); finally because RCC is empirically related with absolute and relative permeability of clastic rocks and transmissibility of fault zones.

Here we implement a Global Hydraulic Element theory (Amaefule et al., 1993) and Generalized Permeability Function approach (Ringrose et al., 2003) at specification of absolute permeability via RCC for mixed clayey-sandy flow units. We also use clay content as a key factor for specification of fault transmissibility multiplier (Walsh et al., 1999) and relative permeability multiplier (Blair, Berryman, 1992).

6.1. Specification of absolute permeability for a continuous reservoir interval

There are two main approaches to physically interpretable deterministic descriptions of absolute permeability applicable for reservoir engineering needs:

1. Global Hydraulic Element (GHE) approach, where isotropic absolute permeability is defined as a continuous function KA(x;0) of empirically driven model parameter vector (x) and porosity (ф) (Amaefule et al., 1993);

2. Generalized Permeability Functions (GPF) approach, where anisotropic permeability is defined as a continuous facies-constrained vector function KA(x;V) of empirically driven model parameter vector (x) and clay content (RCC) - Vin agreement with facies class (Ringrose et al., 2003).

First approach is developed and empirically proven for sedimentary basin scales. Here the absolute permeability function could be expressed in the following general form:

Ka (£;ф) = F(aba2,......a^S; ф), (1)

where parameters a1,a2,..... aNare the normally fixed as local empirical constants; S [1/m] is the tuneable model

parameter - specific surface, which has meaning of integral representative characteristic of pore space fabric (connectivity, tortuosity, etc); ф is the porosity.

Most popular are the following empirical relationships:

Г КА(\;ф) = 10 [А+&ф] - for sandstones (Mello, Karner, 1996);

< КА(\;ф) = ф3/^2 (1-ф)2] - for clay-shale (Magara, 1984); (2)

КА(х;ф) = 1015S10[1°^- 7 2] - for chalk (Baldwin, Butler, 1985).

A more universal formula is derived based on Global Hydraulic Element theory (Amaefule et al., 1993). Here some globally universal hydraulic element (unit) is defined as a representative elementary volume of the total reservoir rock, within which geological and petrophysical properties that affect fluid flow are internally consistent and predictably different from properties of other rock volume. The in situ quality parameter W of any hydraulic unit according to this empirically derived theory could be quantified via by the following quality estimation:

¥ = [0.0314(КА/ф)05]/[ф/(1-ф)], (3)

where common notifications are used. This formula allows reduction of wide range of possible combinations of porosity and permeability relationships into a manageable number of hydraulic units. Moreover, it was shown (Corbett at al., 2003; Cortez, Corbett, 2005; Kooistra, Potter, 2005) that universal link exists between facies and petrophysical parameters W which allows reduction in RM specification elements and rapid population of an inter-well region with КА(\;ф) values. Note that formula (3) coincides with Coseny-Garman formula in (2) after substituting (W/0.0314)2 on S-2.

465

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

It is important to mention that every clastic rock in situ is presentable as a mixture of some pure lithology components. In particularly, the following two-component porosity model is prevalent (Salvadores, 2004):

ф = фс V + ФХ1-П (4)

where V is RCC value varying from V = 1 for pure mud-rocks to V = 0 for pure sandstones ("c" and "s" subscript indexes respectively). We will use this linear form also for representing specific surface of clastic reservoir matrix, which obviously represents mixed terrigenous lithology.

Based on this and on formula (3) we suggest implementation of the following mixed lithology form for defining reservoir quality parameter W

W = 0.0314/[ScV+ Ss(1-V)], (5)

where Sc and Ss are the specific surfaces for pure clay and pure sand lithology respectively.

Finally the isotropic permeability model derived from Global Hydraulic Element theory for considered purposes can be written as follows:

Ka(Sc, Ss, V, ф) = ф3/{^+ Ss(1-V)]2(^)2}. (6)

Note, that continuous permeability function (6) depends on porosity (ф) and three independent lithology controlled constants (Sc, Ss and V), which serve as control RM parameters.

The simplicity and universality are obvious advantages of relevant class of absolute permeability models. The main drawback here is closely connected with main advantage. Due to universality the model is essentially limited by isotropic case.

The second approach to creating interpretable permeability models is more related to the field scale descriptions, where reservoir rock anisotropy is a necessary attribute. Here the clay content and absolute permeability values for pure mud-rock and pure sand-rock are main control parameters. The relevant generalized permeability models are based on textural templates associated with specific background depositional processes. These processes in turn can be attributed with specific bedding type or facies type of reservoir and simulated by corresponding 3D simulator (see SBED method (Wen et al., 1998) or GPM method (Tetzlaff, Harbaugh, 1989) respectively).

In the general case the Generalized Permeability Functions could be expressed as follows:

Ka(x; V) = F(ai ,a2,...a^K, К; V), (7)

where aba2,.......aN are the local empirical constants conditioned by the given facies template; Kc and Ks are the

absolute permeability values for pure mud-rock and pure sand-rock respectively; V is the RCC value.

The main advantages of this model class are in its ability to account for the inherent spatial anisotropy of

reservoir rocks. The main drawbacks are connected with simulation of process-driven constants aba2,....aN and

values Kc and Ks which are unknown in the general sense. In order to make model (7) more universal we combine the Global Hydraulic Element concept with the Generalized Permeability Functions approach. Since isotropic permeability model is generally valid for any pure lithology, we use construction (6) to get Kc and Ks in the form:

Г Kc = (Sc)-2 ф3/[ (1-ф)2];

L ks = (Ss)-2 ф3/[ (1-ф)2].

(8)

Then the anisotropic permeability model (Ringrose et al., 2003) for mixed sandy-clayey lithology could be represented as a RCC-conditioned continuous function of physically interpretable control parameters Sc, Ss and V as follows:

' Kxxyy = K/[(1+V)( kj Kc-1)], V V > V2

I Kx^yy = exp{ [(lnp-lnr) / (lnV2- lnVj)] [lnV- lnVJ + 7}, V2 > V V > V

I Kxxyy = aKc(KJ Kc - 1)-V/b, VV<V:

„ Kzz = [Ks(1-V)+ KcV]/c,

(9)

where empirical constants a,b,c are given by: a = 0.3, b = 0.5, c = 0.55; V1 = 0.3 and V2 = 0.65 define low and upper RCC limits for the relevant flow regimes respectively; n and у correspond to К'xx,yy functions defined for the upper RCC level (V2< V) and the lower RCC level (V<Vi) respectively; Kc and Ks are given by (8); Кxx = Kyy and Kzz are the lateral and vertical components of the vector function Ka given by the following approximation:

K A =

Kx 0 0

0 К 0

yy

0 0 Kzz

(10)

466

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Fig. 3. Permeability for isotropic model computed for pure sandstone (red curve), pure clay (black curve) and their 50:50 mixture (blue dotted curve) in agreement with GHU concept adopted for mixed lithology (formula 6) and Athy's porosity-depth relationship (formula 11)

Clay content

---- Depth = 1000m

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

.... Depth = 2000m

---- Depth = 3000m

Model parameters for Fig. 3 and Fig. 4: sand component - ф03 = 0.45;

CS = 0.00015; SS = 0.10; clay component - ф0с = 0.60;

CC = 0.00037; SC = 10.0

Fig. 4. Permeability vs. clay content trends for anisotropic (a,b) and isotropic (c) models computed for pure sandstone (red curve), pure clay (black curve) and their 50-50 mixture (blue dotted line) in agreement with GPF concept (formulas 8-10) and GHU concept (formula 6)

The plots in fig. 3 represent porosity vs. depth trends and isotropic permeability vs. depth trends computed for pure and mixed clastic rock lithologies. Athy's formula was used here at the computing of porosity vs. depth trend for pure rocks (Magara, 1978):

0(z) = ф0 Exp(-Cz). (11)

Here ф0 and C are lithology controlled compaction constants.

Then the mixture porosity formula (4) was used for mixed lithology. We have used empirically driven compaction constants for pure sands and clay typical for deep water Gulf of Mexico conditions (Reike, Chilingarian, 1974).

Note that the isotropic permeability curve for 50-50 mixed rock lithology is shifted from "pure sand" to "pure clay" lithology side. Thus, a small fraction of clay material in sandy rock significantly decreases the total permeability level.

This effect looks quite similar for lateral permeability components and about one order more significant for vertical permeability component of anisotropic model (8-10) according to the computations displayed in fig. 4. Note that the porosity vs. depth trend (11) was the same as for the results presented above (fig. 3) for isotropic permeability model.

467

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Fig. 5. Stages of generation of the thin-layered reservoir model: a - getting 3D framework; b - getting thin-layered texture; c - getting 3D porosity model initialized; d - getting 3D absolute permeability model specified based on chosen generalized function KA(x;0).

See also fig. 12.

Fig. 12. Comparison of cell (left) and edge (right) transmissibility multipliers approaches implemented for simulation of production history pore pressure maps (b,c) and oil saturation slices (d,e) by the end of synthetic production history.

See Appendix 2 for more details.

The consistency between first (formal edge transmissibility multipliers) and second (cell transmissibility multipliers - 14) approaches has been checked based on reservoir simulation for plausible earth model (see Appendix 1 for more details). Here the industry available reservoir simulator was used for modelling of pore pressure and oil saturation dynamics. The results of this numerical test presented in fig. 11-12 can be considered as confirmation of their consistency.

472

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

6.3. Specification of relative permeability

The presence of several immiscible and non-wetting fluid phases in the pore space of a reservoir matrix diminishes its absolute conductivity for each of these phases.

The essence of this apparent permeability decrease in two-phase case is related to the presence of a new phase of the pore fluid in addition to the initial one within the same pore space volume. This additional fluid portion formally shrinks a pore space volume and hampers passage through it by another phase if the relevant pressure gradient is not zero. However, decreasing of rock conductivity for multiphase fluid flow is not only due to pore space bulk losses but also due to reducing of pore channels' connectivity in their narrow and tortuous parts. This phenomenon is especially strong at low speed flow, when the pore channel connectivity is mainly controlled by capillary forces on contact with non-wetting phases, which could create extremely high entry pressure barriers (Vavra et al., 1992). As it is confirmed by mercury injection tests on clastic rock samples even just a few flow channels locked by residual oil could dramatically reduce the total bulk conductivity (Dake, 1978; Okano et al., 2005).

The process of oil displacement from reservoir by passive or active (injection) water flooding can be roughly treated as a movement of the water invasion front edge toward production well(s). The thickness of the relevant water-oil transition zone will gradually increase with travel time passed from beginning of the process (Bear, Bachmat, 1991; Muskat, 1946). Most of the pore space behind this front edge will be filled with water, but part of it will still hold residual oil. Depending on flow speed and pore space fabric the amount of residual oil could vary in broad saturation ranges [5-50 %]. The nature of such oil pillars consists in losing connectivity between running parts of oil flow, trapping part of oil volume and losing their mobility within the water flow. The opposite is also true. Moreover, there is obviously some virgin water within the pore space of any oil-gas filled reservoir, since HC fluid migration from source rock is considered to be a secondary process which occurs after reservoir matrix dewatering during the primary process of sedimentation-compaction (Jouset, 1986).

So, there are always two consequent processes occurred during any two/multi-phase fluid flow. Namely: first, the given phase permeability decreases sharply as soon as any new immiscible phase appears within the pore space; second, mobility of this phase will be dropped to zero as soon as its saturation approaches a certain limit. In other words, some saturation limit exists, which defines specific volume for residual fluid.

The extension of Darcy Law for volumetric multiphase flow rate Qw/o is obtained by use of a step-down multiplier Krw/o that accounts for the two above-mentioned phenomena. Thus it can be expressed, in the case of two-phase flow ("o" - oil and "w" - water) as follows (Calhoun Jr., 1982):

{

Qo(S,T) = A[KAKro(S)/^o(T)]-[dPJdl; +Po(T)g-sin(0)]; Qw(S,T) = A[KAKrw(S)/Vw(T)\[dPJ% +pw(T)g-sin(0)],

(15)

where A is the cross area normal to the fluid flow vector; KA is the absolute permeability of rock matrix; ^ is the fluid viscosity; dP/^d is the pressure gradient for the given phase; p is the phase density; в is the angle between gravity acceleration and flow rate vectors; T,S are the current subsurface temperature and oil/water saturation respectively.

The relative permeability is treated here as a continuous single-valued function of current saturation. The typical shape of Kr(S) represented in fig. 13 is in good agreement with laboratory data (Muskat, 1949; Kotyahov, 1977; Nikolayevsky, 1996; Zeltov, 1998).

The following features are normally highlighted in connection with the Kr(S) model:

1. sharp decrease of Kr(S) at small departure of relevant phase saturation from 100% level;

2. existence of immobile area (Kr(S)=0) within the saturation interval S < Sw/omin;

3. clear link of the Kr(S) shape with rock matrix lithology;

4. weak link of the Kr(S) shape with fluid viscosity.

There are several analytical Kr(S) models proposed for computing. In particular, the following Chierici's model (Chierici, 1981) was used for data fitting (Christie et al., 2002):

{

Krw(S; S Kro(S; S

< min

min. r< w ; So

min min w ; So

) = exp{-a[(S-Sw ) = exp{-£[(S-Sw’

in)/(1-S

n)/(1-S

min

o

min

o-

S)]C};

S)]D},

(16)

where Swmin, Somin are the residual saturation for immovable water and oil phases respectively; a, b, C, D are the adjustable constants.

Despite good flexibility of this construction, there are four constants in it which require adjustments with no clear physical meaning for each of them.

The robust model in context of data inversion should be tunable by as few model parameters as possible, with clear interpretation for each of them. From this point of view the following functional link proposed by Muskat for two phase fluid flow (Muskat, 1949) looks more attractive:

473

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Г Krw(S; Swmm) = [16S2(S-Swmin)3(1-Swmin)]/[2S2(2-3Swmin)+3S Swmm(3Swmm-2)+ Swmm(4-5 Swmm)]2; l Kro(S; Swmin) = Krw[(1-Sw-S); Swmm].

(17)

Here the same notifications are used as above. Only one new parameter has appeared - a corrector for non-wetting liquid (sw), with a value close to one (sw ~ 0.85-0.95).

The empirical points (markers) correspond to different oil viscosity level. The red-blue and magenta-cyan lines correspond to Chierici's (16) and Muskat's (17) models respectively. Tests made by using Muskat's formula (17) and Chierici's formula (16) proved that the accuracy of approximation of empirical data for two-phase flow system are comparable and limited with 5 % misfit (see fig. 13). At this the applicability of these two formulas for data fitting purposes are evidently different.

Indeed, the only one adjustable model parameter of model (17) is saturation of residual immovable water - Swmin This parameter has clear petrophysical meaning since it is closely related to the peculiarity of specific pore space.

Moreover, the review of published data collected about it on many real HC fields (Kotyahov, 1977; Nikolayevsky, 1996;

Gamatuddinov, 1971) convinces that it has a well-proven lithological interpretation. In connection with the clastic rock materials it is convenient to use summary of lab measurements presented in (Gamatuddinov, 1971) (see fig. 14a).

As it clearly visible from fig. 14a the slope of Swmin(KA) curve becomes steeper as the clay content increases.

From another side the range of variation in absolute permeability for depth interval 1-4 km is about one order less for coarse grained sandstone then for silty sandstone. We believe that the nature of both phenomena is the same. They are both related to the different compaction capacities for pure sandstone and pure clay (Magara, 1978; Mello, Karner,

1996). Thus the clay content (RCC value) seems to be a key factor for definition of residual water saturation in clastic rock.

Fig. 13. Phase permeability measurements for oil and water after (Nikolayevsky, 1996) vs. theoretical curves

18 35 60 100 180 350 600

Absolute permeability [mD]

Fig. 14. The empirical Swmm vs, KA [mD] curves for different lithologies (after Gamatuddinov, 1971) (a) and theoretical trends KA [mD] vs. depth [m] for different quality of sandstones (code and colour of curves)

1-3 - sandstones: (1 - silty; 2 - fine grained; 3 - coarse grained); 4-6 - carbonates; 7 - mudstone. The circle marker indicates Swmin at Ka = 100mD (middle of full range 10-1000mD, in log scale). The relevant depth 1000m was used for fitting

In order to get Swmm linked with reservoir matrix lithology in a more universal way we suggest correcting of empirical relationships Swmin (1-3) displayed in fig. 14a for compaction. This correction can be made by using the following formula to fit with:

s; (v)=(o+vbx )1-

,-z[VCc +(1-V )Cs ]

B

(18)

where z is the depth in [m]; CS,CC are the compaction constants for pure sand and clay respectively, defined according to (Reike, Chilingarian, 1974) as CS = 0.0015m-1 and CC = 0.0037m-1; B0, Bb B2 are the adjustable fitting constants.

Here first linear multiplier in the right part defines the linear term of Swmin(v). The second multiplier is for correction of Swmm(v) for clay compaction with depth.

474

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

We have used the empirical data displayed in fig. 14a and model (18) for identification of constants B0, Bb B2. The following values were obtained based on least square data fitting routine: B0 = 0.13; B1 = 0.52; B2 = 0.5.

So, we suggest formulas (17-18) with the given constants B0, B1, B2 for specification of the relative permeability model for two-phase fluid flow within the clastic rock matrix. It is essentially based on single key model parameter V - clay content (RCC), since the depth is uniquely associated with the given cell. Therefore such specification does not increase the dimension of the relevant inverse problem and simultaneously makes an effective reservoir model more sophisticated. From this point of view this approach to RM specification can be considered to be robust and purposely built.

Practically it is convenient to get relative permeability multiplier values predefined and sorted out in array form KRk:,mwhere indexes k,l,m define depth, RCC and water saturation values respectively. Then specific value of relative permeability multiplier can be found uniquely in accordance with the given cell position (I j, k), attached clay content V iJ>k value and current water saturation (see also fig. 25a). Strictly speaking KRk:lm value should also be updated for the given cell position in agreement with redistribution of saturation field during simulation process, i.e. index "m" is changeable with time.

7. Sensitivity analysis

According to the concept of an effective reservoir model we distinguish in RM specification between fixed structural part (RM framework), fixed petrophysical part (RM initialization) and tuneable part (control model parameter set). The smaller the dimension of the CMP vector space X that can reproduce the production history with the given accuracy, the better the effective RM is. An ideal effective model would be able to reproduce any production history based on adjustment of only one globally universal control model parameter. Practically it is hard to imagine such a situation. Nevertheless it is important to minimize amount of CMP vector components required for adequately accurate simulation. The first step in this direction was made by using a special type of generalized permeability functions (KA, KA, KF and Kr), which all use the clay content related argument (RCC) as a control model parameter among other CMP components and fixed empirical constants.

In the context of data fitting problem it is important to determine how sensitive this parameter is in comparison with other possible parameters, and how stable the relevant sensitivity signature is to the position of the given model parameter vector x within the multidimensional sub-space X of acceptable inverse problem solutions.

The sensitivity vector S(x) on the multi-dimensional model parameter space X is often associated with estimation of discrete Frechet derivative vector:

F(x) = |M[x+Sx] - M[x] | / |Sx|, (19)

where M[] is the operator of the forward problem solution; Sx is the infinitesimal local perturbation of reference model parameter vector xe X.

This estimation is convenient if the analytical solution of the forward problem could replace the numerical solution without losing generality. In the context of the discussed purposes operator M[] implies all stages of reservoir simulation: model gridding, model specification and numerical solution of relevant multiphase fluid dynamic problem. So, F(x) can be estimated only numerically. In connection with the data inversion problems (Menke, 1984) the dimensionless measure of vector F(x) at constant Sx is treated as a sensitivity vector, defined in the same model parameter space X by the set of following i-th components

si(x) = ||M[x+SxJ - M[x] 11/|| M[x] 11, (20)

where || || denotes L2 norm.

If some physical meaning and corresponding unit could be associated with each i-th of model vector components, it is convenient to calculate the local perturbation Sxi of reference vector x as a fixed small fraction of a priori given range of its variation. Then the sensitivity vector (20) can be treated as a combined measure of both model sensitivity and prior uncertainty, which are not the same.

Since all components of CMP vector have clear physical meanings according to the developed specification rules, we will use this last modification and we will understand sensitivity measure as components of the following misfit vector:

Si(x, Axi) = ||M[x+5xi] - M[x] 11/|| M[x] ||, (20*)

where Sxi is the fixed fraction of predefined variation range Axi regardless to the index of component. For example Sx = 0.01Ax.

The multi-run averaged distribution of RM sensitivity vector (20*) can be considered as a signature of the relevant effective reservoir model, which defines the approach to its parametrical upscaling and tuning (see part 2 of this paper for more details). It is important to investigate how stable it is to CMP vector x position

475

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

within the predefined model parameter sub-space Xac X available for adjustments and to the random noise in framework and texture specifications. The more stable the S;(x, Ax,) signature is, the more appropriate the rules of RM specification are and the more universal approaches could be implemented for its tuning (Madatov, Sereda, 2003; Madatov, Sereda, 2005).

7.1. Continuous single-layer reservoir model

Since, the transport properties of continuous reservoir are mainly controlled by its absolute permeability distribution it is logical to use Frechet derivatives of generalised absolute permeability functions in form (19) as a first guess about sensitivity vector (20) for a speculative reservoir model, where relevant CMP specification was implemented. We have used the simplest Athy model to simulate the natural porosity reduction with depth due to compaction of reservoir matrix (11). It was implemented for pure sand and clay lithology components. Then the mixed lithology concept (4) was implemented to simulate isotropic permeability model according to (6). Consequent substitution of (11) into (4) and then into (6) allows getting explicit formula for isotropic permeability as a function of seven parameters: ф0S; CS; SS; ф0с; CC; SC and V. The behaviour of the Frechet derivatives (19) estimated as a function of depth at fixed position of 7-dimensional vector x is displayed in fig. 15.

As it follows from this plot, the first three sensitivity ranks among the seven control parameters specified for single flow element belong to the specific surface constants and clay content components. Their separation from others is several orders at shallow depth and this decreases exponentially with depth. The compaction constants pass ahead in terms of the sensitivity rank below 3000m.

It is important to keep in mind that the presented draft estimation of relative sensitivity ranks gives only very gross estimations of the sensitivity vector for the RM model since the considered F(x) components are not associated with any specific framework, texture and well system data.

In order to come to more specific conclusions we have used multi-run simulation for the speculative single-layer reservoir model (see fig. 16a). The Athy formula (11) and parabolic clay content approximation V(xy) were used for getting synthetic porosity ф (xy,z) be populated according to mixed lithology model (4) for clastic rocks. Then isotropic (6) and anisotropic (8-10) absolute permeability models were used for getting synthetic permeability KA(x,y,z) and KA(x,y,z) be populated within the cellular model (see fig. 16c-d).

The regular and random variations in clay content distributions have been investigated.

Besides the type of the functions V(x,y) which could be potentially used for approximation of lateral variation in clay content the different orientation or just different constants of power function V(x,y) for the same approximation type can be considered as a regular uncertainty in V component of the given CMP vector. The random component of V(x,y) was simulated as an additive white noise 8W with fixed standard variation level o.

Finally, the following model for the lateral clay content distribution according to certain depositional scenario was used:

V(xy) = (A0+SA0) +(A1+SA1) (x+y) + (A2+SA2) (x+y)2+ SW, (21)

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

where SA0, SAl and SA2 are the small uncertain components of the trend constants A0, Aj and A2 respectively; SA0/A0 -SAfAi ~SA2/A2 << 1; SW - the white noise component of the model (19).

The simulations of oil production history were done in agreement with well configuration presented in fig. 16 by using an industry standard reservoir simulator.

Multiple simulations reveal regular features of this model. Namely:

1. The small regular perturbation (SA0/A0 and SA2/A2 < 10 %) generates noticeable distortion of the reference production history expressed in the form of well oil production rate curve (WOPR curve).

2. The same or bigger [o (SW) < 25 %] random perturbation SW does not generate a noticeable distortion of the reference WOPR curve.

Some of the instructive examples are represented in fig. 17 and 18.

.„3

--- Fs 8.687 4

— Fm 5.587 5

... Ss 10.633 3

... Sm 151.902 1

---Cm 2.179 7

---Cs 2.542 8

--- V 15.555 2

500 1000 1500 2000 2500 3000 3500 4000

Depth

Fig. 15. Comparison of range normalized Frechet derivatives of isotropic permeability KGHU as a function of depth at clay content V = 0.3 with the value and rank reported

476

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Fig. 16. Continuous single-layer reservoir model in 3D display views a - framework part; b - porosity model; c,d - vertical and lateral components of absolute permeability respectively; e - oil saturation model by the end of synthetic production history

Fig. 17. The clay content (above) vs. permeability (below) 2D distribution in map views

a, d - reference model:

a - V(xy)= 0.70 - 0.0085(x+y)2; d - KxxydtoVXy)];

b, e - regularly perturbed model:

b - V*(x,y)= 0.60 - 0.0075(x+y)2; e - Kxxyy[ф ;x,V*(xy)]

c, f - irregularly perturbed model:

b - V(xy)= 0.70 - 0.0085(x+y)2 + SW; a(SW) = 0.10

Table 2. Sensitivity rank vs. additive SW noise level for reference reservoir model

CPM component Sensitivity at additive noise level (a): Aver. Rank Comment

0 .01 .05 .1 .25

Fc .003 .005 .026 .053 .108 8 See formulas (4,11)

Cc .006 .010 .049 .101 .213 6 See formulas (4,11)

Sc .175 .170 .150 1 See formulas (8-10)

Fs .007 .004 .027 .059 .123 7 See formulas (4,11)

Cs .017 .024 .061 .124 .268 5 See formulas (4,11)

.169 See formulas (8-10)

.064 .295 Coefficients in V(x,y) approximation by (21)

.022 .023 .085 .132 4

477

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Fig. 19. Averaged sensitivity signatures of the perturbed reference model based on multiple reservoir simulations

a - the variation in the white noise level a(SlV) at V(x,y) approximation is 0.0-0.25; b - the variation in framework error is 0.0-0.05; c - the variation in CMP vector position in normalize metric (xi+Sx,)/xi is 0.0-0.5 (see also formula 18*)

The sensitivity signature of the reference reservoir model was not changed significantly against gradual increasing of random perturbation level (see Table 2). In particular three upper ranks were distributed among components SC, SS and coefficient A0 in V(xy) approximation according to (21). The corresponding sensitivity components are highlighted by red - first rank; orange - second rank and yellow - third rank in Table 2.

The sensitivity signature was proven to be stable also to the tolerant errors in the RM framework and to the random variations in vector x position inside the predefined model parameter sub-space Xa available for model tuning. The averaged sensitivity signatures are displayed in fig. 19.

7.1. Discontinuous single-layer reservoir model

By analogy with the continuous case we first evaluate the Frechet derivative component for control parameter SV (uncertain fraction of clay content) computed in agreement with explicit fault transmissibility model (12-14). The analytical solution reveals highly non-linear behaviour of STm/SV within the range of SV variation [0,1] (see fig. 20).

Fig. 20. The numerical estimation of SV component ol Frechet derivative for transmissibility model (12-14).

Normalized displacement D is curve code

478

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Fig. 22. Synthetic oil saturation at final time step corresponding to locked (left) and open (right) leak zone in Fault 1 in stream line view (a,b) and slice (Z = 2105m) view (c,d) representations

Fig. 21. Combine map of lateral permeability (Log scale) tied to the top of blocked reservoir (see also Appendix 2, fig. A2.1)

Following parameterization of three-fault system is used:

Fault 1:

Segmentl (1.1 - 1.2) D = 1.5; SV1j = 0.4;

Segment2 (1.2 - 1.3) D = 0.5; SV12 = 0.0 - 0.8;

Segment3 (1.3 - 1.4) D = 1.5; 8V13 = 0.4;

Fault 2: D = 0.9; SV2 = 0.25;

Fault 3: D = 0.9; SV3 = 0.15.

The potential leak segment 1.2 is indicated with arrow. The relevant CMP component (SV12) is free for adjustment and the variation range is highlighted with bold font

Fig. 23. Abnormal sensitivity of the reference model to the perturbation in fault leak zone contrast (parameter SV)

a - variations of production history for well P2 against fixed perturbation of CMP components (curve code); b - sensitivity signatures based on synthetic WOPR curves for well P1 and well P2.

479

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

As expected, the transmissibility function associated with some fault barrier becomes dependent on the uncertain fraction of clay content at almost sealing conditions, i.e. when normalised fault displacement D approaches 1 and smear clay potential of the fault zone also approaches the upper limit.

Evidently, an extension of the framework model by introducing of a seismically visible fault system and adding corresponding components to the CMP vector x will change the sensitivity signature of the updated reservoir model. However, such apparent RM complication will not obviously complicate the solution of the associated data inversion problem. Indeed, when the contrast between upper and lower rank sensitivity components increases and the amount of upper level sensitivity components decreases the condition for multidimensional model adjustments could be even more favourable in comparison with a simpler but less sharply sensitive model.

In order to check this insight and to get more specific sensitivity estimations we have complicated the above reservoir model by adding some fault segments into its framework. The well control schedule at simulation of synthetic production history was also complicated in order to make the model more plausible (see Appendix 2 for more details).

The anisotropic permeability model for continuous reservoir interval was preserved (see fig. 16c-d). The fault transmissibility multipliers were computed based on formulas (12-14) in agreement with framework and fault parameters indicated in fig. 21.

Two scenarios were simulated: fully sealed fault 1 (see fig. 22, left) and partially sealed fault 1 with open leaking zone (see fig. 22, right). The locked and open conditions were simulated by putting uncertain clay fraction (SV) for segment 1.2-1.3 on right (SV12 = 0.8) and left (SV12 = 0.0) limits respectively. The differences in pressure and oil saturation 3D distribution at the final time step appeared to be quite dramatic (see fig. 22c-d).

Thus, the leakage zones associated with the sealed fault could have a critical impact on production dynamics especially for the closest production well(s). The relevant CMP component is expected to be highly sensitive to the simulation results and consequently it can be used as one of the key parameters in a history matching. The numerical estimations of the reference model sensitivity to the uncertain clay content (SV12) confirmed this conclusion (see fig. 23).

It follows from analysis of perturbations in leakage zone contrast that the highest impact on the synthetic WOPR curve is connected with the nearest production well P2 (see fig. 23b). The shapes of synthetic oil production rate curves for well P2 are not uniform but reveal the existence of a leak fault zone, which focuses relevant streamlines during second phase, when injector I2 and producer P2 were kept active and partially during third phase when only producer P2 was active (see Appendix 2).

Based on this it appears to be possible to recognise such a particular feature of the blocked reservoir as a leak across the sealed fault through a relatively permeable zone. Indeed, the influences of perturbation in different CMP components on synthetic WOPR curve are not distinguishable from each other as soon as activity of wells is not changeable during production. In contrast, the impact of leakage zone contrast controlled by SV component on synthetic WOPR curve can be easy recognised among other impacts (see fig. 23a). The closer the relevant production well is to the potential leakage zone the bigger difference in WOPR shape due to this phenomenon.

7.3. Discontinuous single-layer reservoir model with lithology depended relative permeability

It follows from definitions (17-18) that the relative permeability for clastic reservoir rock matrix could be presented as a function of its clay content, which in turn is presentable by lateral trend function (21) in agreement with litho-facies correspondence of the reservoir. Thus, there are no additional components for control parameter vector x demanded for an effective RM specification. Still, the meaning of a clay content value specified for the given cell is extended when the lithology-dependent model of relative permeability (17-18) is implemented. This is because it now controls the relative permeability in addition to the absolute permeability. Therefore, the contribution of relevant control parameter components into a sensitivity signature is expected to be increasing (see fig. 24).

The influences of a clay content variation on the shape of Kr(S;Swmm) function and on the result of reservoir simulations are illustrated in fig. 25. The reference reservoir model has been used for simulation (see Appendix 2).

Fig. 24. Sensitivity signatures of the reference model at the constant (a) and lithology depended (b) relative permeability models

480

Вестник МГТУ, том 13, №3, 2010 г.

стр. 457-494

Fig. 25. The influence of constant clay content level specified for the reference reservoir model (code of curves) on relative permeability variations in agreement with model (17-18) (a) and on synthetic water (a) and oil (b) production rates for well P2 (see Appendix 2)

8. Conclusion

Apart from pure formal reservoir descriptions such as Response Surface Model (Salhi et al., 2005; Carrears, Turner, 2006), the complexity of the modern 3D reservoir simulation in terms of model attributing and parameterization remains to be unsuitable as it is for history matching purposes. This situation exists mainly because of the higher priority of the forward problem solution against data inversion associated with history matching, which was historically formed in reservoir engineering. Namely, the various numerical schemes of the forward problem solution formulated for 3D multiphase fluid flow were created first and then the model parameterization routine was adjusted for simulation purposes. As a consequence of this the achieved results were, and remain, helpful for better understanding of simulated phenomena but they are often useless in production practice for prediction of the observed phenomena.

One of the ways to improve this imbalance is reduction of the dimension of the relevant inverse problem by focusing on attributing and parameterization algorithms for such a model parameter space, which presumes data inversion later on. Therefore the main vector of the presented approach to an effective reservoir model generation lies in this direction.

The stochastically driven inter-well population of reservoir transport properties (porosity-permeability) is broadly popular nowadays at the reservoir simulation and history matching (Donselaar, Geel, 2003; Petersen et al., 2003; Raisson, Temple, 2004). This is mainly because of simplicity and soundness of produced synthetic data distributions.

However pure stochastic models ignore underground deposition processes and lead to poor interpretable calibration results.

We suggest a process-based model approach to the reservoir model architecture, in which all available sedimentological knowledge can be used as a guideline (Donselaar, Geel, 2003; Ringrose et al., 2003) for well data interpretation and model parameter specification. In particular, the developed concept implies separation of the fixed structural part (geometry) and the initialized petrophysical part (pore fluid and rock matrix properties) from the unfixed (adjustable) part defined as a texture component of an effective reservoir model.

We also suggest using generalized absolute and relative permeability functions for specification of reservoir transport properties, as they are most sensitive to the simulation results. The construction of the generalized absolute permeability functions (4,6) or (4,8-10) allows reduction of the number of independent arguments demanded for specification of flow units defined in RM texture. Each of these arguments can also be considered as a control model parameter (CMP vector component) in the context of reservoir model identification. In the case where 3D porosity distribution for the reservoir can be initialized based on available wireline log data the set of necessary control parameters includes minimum components per flow unit, which are: Ss, Sc, V. Here Ss, Sc are the specific surfaces for pure sand lithology and pure clay lithology; V is the clay content level (RCC). Normally RCC component allows more comprehensive description then just V-constant level. Particularly in the case of tidal, shore or deep marine facies interpretation the lateral variability of clay

481

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

content within the specific flow unit could be represented by trend approximation. When using the parabolic function (21) it includes two more components (Ab A2) to be specified in addition to A0.

Note, that all of the above-mentioned control model parameters imply clear geophysical interpretation. Namely SS and SC represent integral property of the given pore space fabric; V(x,y) trend function used for the lateral clay content description accommodates geological and litho-facies information about reservoir accumulation and sedimentation.

The generalized function implemented for specification of the lateral flow barrier permeability (fault transmissibility) in agreement with (12-14) implies only one control model parameter - uncertain clay content fraction in a fault membrane - SV.

The relative permeability specification based on extended Musket's model (17-18) does not require any additional components to be used for CMP vector specification, since it uses the same RCC component as it is used for specification of the absolute permeability.

Finally, the proposed RM architecture and rules for its specification appeared to be fitted for purposes of data inversion from one side and is still usable for simulation by industry available simulators from another side. The features of such a model actually coincide with demands to the robust reservoir model, which were summarized by K. Azis in his "Ten Golden Rules for simulation engineers" (Aziz, 1989) as follows:

1. Clearly define the aim of reservoir model setting. Do it against available data

2. The model needs to be facilitated. The simplest model capable to achieve your aim is the best model

3. Subdivide big problem into consequence of solvable pieces

4. Bigger model is not obviously better. Quite often it is opposite

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

5. Stay close to physics. Check the reality of your model parameter values

6. Don't exaggerate the significance of your model. Comparison of possible scenarios is a good result

7. History matching problem has multiple solutions. Most reasonable one should satisfy not only formal criteria but must be based on a priori available geology and petrophysics otherwise a prediction based on it could be wrong

8. Don't average extremes in permeability, since they normally belong to key flow units and/or barriers

9. The local scale measurements are not applicable to field scale. Simple averaging could change origin of field scale phenomenon

10. Model cannot substitute lab and in situ measurements.

The introduced approach to the RM architecture design and population of its transport property allows simulation of various geology settings in terms of tectonics and litho-facies scenarios, which can be used for reproducing observed multi-well production history. The most important results of the sensitivity analysis made for the plausible earth model and multiple reservoir simulations (see Appendix 2) can be summarised as follows:

1. The regular variations in clay-content template (trend type; trend form; trend orientation) appeared to be much more significant than its irregular variations, on condition that an additive white noise level in V(xy) approximation does not exceed 25 %. Therefore, the low-frequency regular component of RM specification such as background clay-content trend V(x,y) is supposed to be a main factor of observed long-term reservoir dynamic and it has to bring a major input for its further prediction.

2. The sensitivity signature of both continuous and blocked reservoir models appeared to be rather stable and invariant to the start position on testing CMP vector within the model parameter space, well configuration and additive high frequency components of lateral permeability distribution (white noise level).

Based on the sensitivity rank analysis it is possible to distinguish between the primary and secondary important CMP components in the context of model parameter adjustments. This can be interpreted as an important prerequisite for further reducing of an initial dimension of the model parameter space.

Acknowledgements

Author is grateful to Schlumberger Moscow Research centre and to Murmansk State Technical University for giving him the opportunity to develop the theory and algorithms applicable at the 3D reservoir geo-fluid system modelling-calibration and prediction. Author also is grateful to Youry Kulikov for scientific cooperation during development of an effective reservoir model concept. Author appreciates the help from Eamonn Doyle during editing this paper and preparing it for publication.

List of references

Advin A., Johnson M. Analysis of faulting in porous sandstones. J. of Structural Geology, v.5, p.19-31, 1985. Al-Busafi B., Fisher Q.J., Harris S.D., Kendall M.U. The impact of faults representation on history match and future generated seismic impedance response in reservoir models: A case study for Pierce Field, North Sea.

Paper presented at the SPE Europec/EAGE Annual Conference held in Madrid, Spain, 13-16 June, 2005.

482

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Al-Deeb M., Bloch H., Salem El-Abd, Sharfeddine M. Fully integrated 3D reservoir characterization and flow simulation study: A field case example. Paper presented at the 10 Abu Dhabi International Petroleum Exhibition and Conference, 12-13 October, 2002.

Alsos T., Eide A., Astratti D., Pickering S., Benabatos M., Dutta N., Mallick S., Schultz G., denBoer L., Livingstone M., Nikel M., Sonneland L., Schlaf J., Schoepfer P., Sigismondi M., Soldo J.C., Sternen

L.C. Seismic application through the life of the reservoir. Oilfield Review, September, p.48-65, 2002.

Amaefule J.O., Altunbay M., Tiab D., Kersey D.G., Keelan D.K. Enhanced reservoir description: Using core and log data to identify hydraulic (flow) units and predict permeability in uncored intervals/wells. SPE 26436, 68th Ann. Tech. Conf and Exhibit., Houston, Tx., 1993.

Aminzadeh F., Barhen J., Glover C.W., Toomarian N.B. Reservoir parameter estimation using hybrid neural network. Computers & Geoscience, v.26, p.869-875, 2000.

Antonellini M., Aydin A. Effect of faulting on fluid flow in porous sandstones: Petrophysical properties. The American Assoc. of Petroleum Geologists Bull., v.78, N 3, p.355-377, 1994.

Ates H., Bahar A., El-Abd S., Charfeddine M., Kelkar M., Datta-Gupta A. Ranking and upscaling of geostatistical reservoir models using streamline simulation: A field case study. Paper SPE 81497 presented at the SPE 13th Middle East Oil Show & Conference held in Bahrain, 9-12 June, 2003.

Averbukh A.G. Investigation of sedimentary rock properties and content by using seismic. M., Nedra, 232 p., 1982.

Aziz K., Settari A. Petroleum reservoir simulation. London, Applied Science Publishers, 489 p., 1979.

Aziz K. Ten golden rules for simulation engineers. J.Petrol. Technol., v.41, N 11, p.1157, 1989.

Baldwin B., Butler C.O. Compaction curves. The American Assoc. of Petroleum Geologists Bull., v.69, N 4, p.622-626, 1985.

Baltrusch S., Bellmann L., Lockner R.K., 0yvann E., Tharahardt R. The Akshalbulak channel model: An interdisciplinary approach to handle complex reservoirs. Paper A011 presented at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Batycky R.P., Blunt M.J., Thiele M.R. A 3D field scale streamline-based reservoir simulator. SPE Reservoir Engineering, November, p.246-254, 1997.

Bear J., Bachmat Y. Introduction to modelling of transport phenomena in porous media. Kluwer Academic Publishers, London, р.553, 1991.

Berryman J.G., Berge P.A., Bonne B.P. Estimating rock porosity and fluid saturation using only seismic velocities. Geophysics, v.67, N 2, p.391-404, 2002.

Bhatt A., Helle H. Committee neural networks for porosity and permeability prediction from well logs.

Geophysical Prospecting, v.50, p.645-660, 2002.

Blair S.C., Berryman J.G. Permeability and relative permeability in rocks. In: Fault Mechanics and Transport Properties of Rocks, p.169-186, 1992.

Boassen T., Palatnik B. Pore- to field-scale multi-phase upscaling for IOR. Paper SPE 94191 prepared at the SPE Europec/EAGE Annual Conference held in Madrid, Spain, 13-16 June, 2005.

Borg H. Fault controlled pressure modelling in sedimentary basins. A thesis presented for the degree of Doctor Ingenior of the Norwegian Univ. of Science and Technology, Department of Math. Science, 156 p., 2000.

Bouvier J.D., Sijpesteijn K., Kleusner D.F., Onyejekwe C.C., Van Der Pal R.C. Three-dimensional seismic interpretation and fault sealing investigations, Nun River field, Nigeria. The American Assoc. of Petroleum Geologists Bull., v.73, p.1397-1414, 1989.

Breton P., Brechet E., Guillon S., Bertini F., Demichelis J. Seismic shale gouge ratio (SGR) evaluation and interpretation for fault characterization. Paper P026 presented at EAGE 67th Conference & Exhibition, Madrid, Spain, 13-16 June, 2005.

Calhoun Jr. J.C. Fundamentals of reservoir engineering. University of Oklahoma Press, Norman, 426 p., 1982.

Calvert R. Insights and methods for 4D reservoir monitoring and characterization. SEG and EAGE issue, Tulsa, Okla, 219 p., 2005.

Carrears P.E., Turner S.E. Tahiti: Development strategy assessment using design of experiments and response surface model. Paper SPE 100656 presented at Western Regional/AAPG Pacific Section/ GSA Cordilleran Section Joint Meeting held in Anchorage, Alaska, USA, 08-10 May, 2006.

Cheong Y.P., Gupta R. Experimental design and analysis methods for assessing volumetric uncertainties. SPE Journal, September, p.324-335, 2005.

Chierici G. Novel relations for drainage and imbibition relative permeabilities. SPE10165, 1981.

Christie M., Subbey S., Sambridge M. Prediction under uncertainty in reservoir modelling. Paper presented on 8th European Conference on the Mathematics of Oil Recovery, Freiberg, Germany, 3-6 September, 2002.

Christie M.A., Mansfield M., King P.R., Barker J.W., Culverwell I.D. A renormalization-based upscaling technique for WAG floods in heterogeneous reservoirs. Expanded Abstracts (Society of Petroleum Engineers), p.353-361, 1995.

483

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Collazos J.P., Nieto J.A., Peet D., Dobbs W. Integrating disciplines in a shared earth model environment: Reservoir characterization of a Nisku Pool in S.Alberta. Paper pres. at CSPG/CWLS Conference, June, 2003.

Corbett P.W.M., Ellabad Y., Mohammed Kh., Pososyaev A. Global hydraulic elements: Elementary petrophysics for reduced reservoir modelling. Paper F-26 presented on EAGE 65th Conference & Exhibition, Stavanger, Norway, 2-5 June, 2003.

Cortez O.D., Corbett P.W.M. Time-lapse production logging and the concept of flowing units. Paper SPE 94436presented at the SPE Europec/EAGE Annual Conference held in Madrid, Spain, 13-16 June, 2005.

Cuddy S.J. Litho-facies and permeability prediction from electrical logs using fuzzy logic. Paper SPE 49470 presented at the 1998 SPE Abu Dhabi International Petroleum Exhibition and Conference, Abu Dhabi, UAE, 11-14 October, 2000.

Dake L.P. Fundamentals of reservoir engineering. University of Oklahoma press, Elsevier, 443 p., 1978.

Davis M., David M. Automatic kriging and contouring in the presence of trends (Universal kriging made simple). The Journal of Canadian Petroleum Technology, 17(1), 1978.

Denesh A. PVT and phase behaviour of petroleum reservoir fluids. Amsterdam-Tokyo, Elsevier, 400 p., 1998.

Donselaar M.E., Geel C.R. Reservoir architecture model for heterolithic tidal deposits. Paper P034 presented at EAGE 65th Conference & Exhibition, Stavanger, Norway, 2-5 June, 2003.

Eide A.L., Omre H., Ursin B. Prediction of reservoir variables based on seismic data and well observations. Jour. Am. Stat. Assoc., v.97, N 457, p.18-28, 2002.

Ellis J.A. Well logging for Earth scientists. Elsevier Scientific Publishing Company, 419 p., 1987.

Evans B. Fault mechanics and transport properties of rock. Edited by Academic Press, 523 p., 1992.

Fisher Q.J., Knipe R.J. The permeability of faults within siliciclastic petroleum reservoirs of the North Sea and Norwegian Continental Shelf. Marine and Petroleum Geology, v.18, p.1063-1081, 2001.

Fisher Q.J. Recent advances in fault seal analysis as an aid to reservoir characterisation and production simulation modelling. Paper SPE 94460 presented at EAGE 67th Conference & Exhibition, Madrid, Spain, 13-16 June, 2005.

Franseen E.K., Watney W.L., Kendall C.G., Ross W. Sedimentary modeling, computer simulations and methods for improved parameter definition. Kansas Geological Survey Bulletin, v.233, 524 p., 1991.

Fristad T., Groth A., Yielding G., Freeman B. Quantitative fault seal prediction: A case study from Oseberg Syd. In NPF Special Publication, p.7, 1997.

Furre A-K., Munkovold F.R., Nordby L.H. Improving reservoir understanding using time-lapse seismic at the Heidrun Field. Paper A-20 pres. at EAGE 65th Conference & Exhibition, Stavanger, Norway, 2-5 June, 2003.

Harding T.P., Tuminas A. A structural interpretation of hydrocarbon traps sealed by basement normal fault blocks at stable flank of foredeep basins and rift basins. The American Assoc. of Petroleum Geologists Bull., v.53, p.812-840, 1989.

Hippler S.J. Microstructures and diagenesis in North Sea fault zones: Implication for fault-seal potential and fault migration rate. Chapter 7 In: "Seals, traps and petroleum systems" / ed. by Surdam, The American Assoc. of Petroleum Geologists Bull., v.67, p. 103-115, 1997.

Hristopulos D.T., Christakos G. Renormalization group analysis of permeability upscaling. Stochastic Environmental Research and Risk Assessment, v.13, p.131-160, 1999.

Hull J. Thickness-displacement relationships for deformation zones. J. of Structural Geol., v.10, p.431-435, 1988.

Idrobo E.A., Choudhary M.K., Datta-Gupta A. Sweep volume calculations and ranking of geostatistical reservoir models using stream line simulations. SPE 6255 7 paper presented at the SPE/American Assoc. of Petroleum Geologists Western Regional Meeting held in Long Beach, California, 19-23 June, 2000.

Johnson H.D., Stewart D.J. Role of clastic sedimentology in the exploration and production of oil and gas in the North Sea. Geological Society Special Publication, N 18, p.219-231, 1985.

Joreskog K.G., Klovan J.E., Reyment R.A. Geological factor analysis. New York, Elsevier SPC, 223 p., 1976.

Joseph P., Fornel E.Du., Eugen T., Guillocheau F., Robin C. Stratigraphic and reservoir architecture of tectonically-confined turbidite ramps: The example of the Gres d'Annot system (Eocene - Oligocene, French Alps). Paper A032 pres. at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Jouset T.P. Migration of hydrocarbons through the sedimentary rocks. M., Nedra, p. 188, 1986 (in Russian).

Kalkomey C.T. Use of seismic attributes as predictors of reservoir properties: Potential risks. 66th Ann. Internat. Mtg., SEG, Expanded Abstracts, 1996.

Kanevskaya R.D. Math modelling of reservoir engineering. Moscow-Igevsk, 140 р., 2002 (in Russian).

King M.J., Datta-Gupta A. Streamline simulation: A current perspective. In Situ, N 22(1), p.91-140, 1998.

Knoth O., Haller N. Stochastic inversion of OBC PZ data at the Eldfisk Field. Paper C026presented at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Knott S.D. Fault seal analysis in the North Sea. The American Assoc. of Petroleum Geologists Bull., v.77, N 5, p.778-792, 1993.

484

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Kooistra J.A., Potter D.K. The links between global hydraulic elements, facies, grain size and sorting. Paper F008presented at EAGE 67th Conference & Exhibition, Madrid, Spain, 13-16 June, 2005.

Krylov A.P., Glogovsky M.M., Mirchink M.F., Nikolayevsky N.M., Charny I.A. Foundations of an oil reservoir engineering. GNTINGTL, Moscow-Leningrad, 416 p., 1948 (in Russian).

Labourdette R., Seguin J., Biver P., Temple F., Hegre J-A. 3D modelling of stacked turbidite channels: Impact on dynamic reservoir simulations. Paper A-13 presented at EAGE 65th Conference & Exhibition, Stavanger, Norway, 2-5 June, 2003.

Lerdahl T.R., Rustad A.B., Theting T.G., Stensen J.A., 0ren P.E., Bakke S. Pore- to field-scale multi-phase

upscaling for IOR. Paper 94191 pres. at the 67th EAGE Ann. Conf., Madrid, Spain, 13-16 June, 2005.

Liu K., Paterson L., Wong P., Qi D.A. Sedimentological approach to upscaling. Transport in Porous Media, v.46, p.285-310, 2002.

Lohne A., Virnovski G., Durlovsky L.J. Two-stage upscaling of two-phase flow: From core to simulation scale.

Paper SPE 89422 pres. at the SPE/DOE symposium on improved oil recovery, Tulsa, 17-21 April, 2004.

Lopez S., Cojan I., Rivoirard J., Fouche O., Galli A. Process-based stochastic modelling on the example of Meandering Channelised Reservoirs. Paper P042 presented at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Madatov A.G., Sereda A.-V.I. The pre-drill and real time while drilling overpressure prediction based on Effective Basin Model concept. Proc. of the Murmansk State Techn. University, v.9. p.361-388, 2006.

Madatov A.G. The overpressure driven seismic velocity response. Review of standard models and methods for extraction in the context of basin modelling approach to overpressure prediction. Proceedings of the Murmansk State Techn. University, v.8, N 1, p.84-119, 2005a.

Madatov A.G., Sereda A.-V.I. A multi-well data inversion purpose-built for calibration of an effective basin model at overpressure prediction. Proc. of the Murmansk State Techn. Univ., v.8, N 1, p.44-83, 2005b.

Madatov A.G., Sereda A.-V.I. The decomposition of 3-D overpressure evolution model in basin scale and its application to the fault seal analysis. Proc. of the Murmansk State Techn. Univ., v.4, N 1, p.79-96, 2001.

Madatov A.G., Sereda A.-V.I. The forward and inverse problems of the fluid dynamics in basin modelling applied to the pore pressure prediction within the sedimentary basins. P.1. Theory aspects. Proceedings of the Murmansk State Techn. Univ., v.3, N 1, p.89-114, 2000.

Madatov A.G., Sereda A.-V.I. Upscaling of a basin model required for pore pressure prediction before and during drilling. Proceedings of the Murmansk State Techn. Univ., v.6, N 1, p.119-144, 2003.

Magara K. Compaction and fluid migration. Elsevier Scientific Publishing Company, p.319, 1978.

Manzocchi T., Heath A.E., Walsh J.J., Childs C. The representation of two-phase fault-rock properties in flow simulation models. Petroleum Geo-science, v.8, p.119-132, 2002.

Mello U.T., Karner G.D. Development of sediment overpressure and its effect on thermal maturation: Application to the Gulf of Mexico basin. The American Assoc. of Petroleum Geologists Bull., v.80, N 9, p.1367-1395, 1996.

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Menke W. Geophysical data analysis. Discrete inverse theory. Academic press, New York, 312 p., 1984.

Merriam D.F., Davis J.C. Geologic modelling and simulation, sedimentary systems. Kluwer Academic/Plenum Publishers, New York, 352 p., 2001.

Miall A.D. Principles of sedimentary basin analysis. Springer-Verlag, Berlin, p.616, 2000.

Moushin I.A., Brodov L.U., Kozlov E.A., Hatjanov F.I. Structure-formation interpretation of seismic data. M., Nedra, 299 p., 1990 (in Russian).

Myers R.H., Montgomery D.C. Response surface methodology: Process and product optimization using designed experiments. New York City, Wiley, 381 p., 1995.

Nikolayevsky V.N. Geomechanics and fluid dynamics. M., Nedra, 447 p., 1996 (in Russian).

Nordahl K., Ringrose Ph., Wen R. Petrophysical variation as a function of scale and implications for core-log integration. Paper F-25 pres. at EAGE 65th Conference & Exhibition, Stavanger, Norway, 2-5 June, 2003.

Norris R.J., Lewis J.J.M. The geological modelling of effective permeability in complex heterolithic facies. Paper SPE 22692 pres. at the 66th Annual SPE Techn. Conf. and Exhib., Dallas, TX, 6-9 October, 1991.

Oil and Gas Geology handbook. M., Nedra, 480 p., 1984 (in Russian).

Okano H., Pickup G.E., Christie M.A., Subbey S., Sambridge M., Monfared H. Quantification of uncertainty in relative permeability for coarse-scale reservoir simulation. Paper SPE 94168 presented at the SPE Europec/EAGE Annual Conference held in Madrid, Spain, 13-16 June, 2005.

0ren P.E., Bakke S. Reconstruction of Berea sandstone and pore-scale modelling of wettability effects. J. Petroleum Science & Engineering, v.39, p.177-199, 2003.

Petersen S.A., Farrelly B.A., Braathen B.I. Modelling of fluid effects in compartmentalised Areaso: An Oseberg south case study. Paper A-25 pres. at EAGE 65th Conf & Exhib., Stavanger, Norway, 2-5 June, 2003.

485

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Raisson F.H., Temple F.M. Impact of sedimentary heterogeneity on reservoir monitoring in a Turbiditic Channel Complex (Angolan deep offshore Girassol field case). Paper D017 presented at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Reike H.H., Chilingarian G.V. Compaction of argillaceous sediments. London, New York, 424 p., 1974.

Reynolds A.C., He N., Chu L., Oliver D.S. Reparameterization techniques for generating reservoir descriptions conditioned to variograms and well-test pressure data. SPE Journal, December, p.413-426, 1996.

Ringrose P.S., Skjetne E., Elfenbein C. Permeability estimation functions based on forward modeling of sedimentary heterogeneity. SPE 84275 presented at the SPE Annual Techn. Conference and Exhibition, Denver, CO, 5-8 October, 2003.

Robertson E.C. Relationship of fault displacement to gouge ratio and breccia thickness. Institute of Mining Engineers, v.35, p.1426-1432, 1983.

Salhi M.A., Van-Rijen M., Wei L., Dijk H., Alias Z., Upadhyaya A., Lee H. Structured uncertainty assessment for a mature field through the application of experimental design and response surface methods. Paper SPE 93529 presented at 14th SPE Middle East Oil-Gas Show Conference, held in Bahrain International Exhibition Center, 12-15 March, 2005.

Sandjivy L., Leron A., Boudon S. The reservoir qualification of seismic data using geostatistics. Paper P273 presented at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Sarma P., Durlofsky L.J., Aziz K. Efficient closed-loop production optimization under uncertainty. Paper SPE 94241 presented at the SPE Europec/EAGE Annual Conference held in Madrid, Spain, 13-16 June, 2005.

Schelkachev V.N. The reservoir engineering worldwide and in Russia. Institute of computer research, Moscow-Igevsk, 132 p., 2002 (in Russian).

Shang B.Z., Yang H.D., Bourgeois S.A. Stochastic fault seal analysis and reservoir simulation for field with thin layers and many faults. Paper F010 presented at EAGE 68th Conference & Exhibition, Vienna, Austria, 12-15 June, 2006.

Talukdar S., Brusdal L.H. Reservoir management of the Njord Field. Paper SPE 92419 presented at the SPE Europec/EAGE Annual Conference held in Madrid, Spain, 13-16 June, 2005.

Tarantola A. Inverse problem theory: Methods for data fitting and model parameter estimation. Elsevier, Netherlands, p.386, 1987.

Ter-Sarkisov R.M. The development of natural gas fields. M., Nedra, 659 p., 1999 (in Russian).

Tetzlaff D.M., Harbaugh J.W. Simulating clastic sedimentation, van Nostrand Rheinold Series on Computer.

Methods in Geosciences, New York, p.202, 1989.

Tipping D.E., Deschenya M., Deimbacher F., Kovyazin D. Increasing of production efficiency for a HC field in Siberia based on quantitative production risk estimations, SPE 101808. Paper presented on SPE Conference held in Moscow, Russia, 3-6 October, 2006 (in Russian).

Van Wagoner J.C., Mitchum R.M., Campion K.M., Rahmanian V.D. Siliciclastic sequence stratigraphy in well logs, cores, and outcrops: Concepts for high-resolution correlation of time and facies. Am. Assoc. Pet. Geol. Meth. Explor., Ser. 7, p.55, 1990.

Virnovsky G.A., Friis H.A., Lohne A.A. Steady-state upscaling approach for immiscible two-phase flow.

Transport in Porous Media, v.54, p.167-192, 2004.

Walsh J.J., Watterson J., Heath A., Gillespie P.A., Childs C. Assessment of the effects of sub-seismic faults on bulk permeabilities of reservoir sequences. Geological Society, London, Special Publications, v.127, p.99-114, 1998.

Walsh J.J., Nell P., Yielding G. Fault transmissibility multipliers for flow simulation models. Petroleum Geoscience, v.5, p.53-63, 1999.

Wang P., Yotov I., Wheeler M., Arbogast T., Dawson C., Parashar M., Sephernoori K. A new generation EOS compositional reservoir simulator: Part I - formulation and discretization. Paper SPE 37979 presented in 1997th SPE Reservoir Simulation Symposium, Houston, TX, USA, 1997.

Wen R-J., Martinius A.W., Naess A., Ringrose P. Three-dimensional simulation of small scale heterogeneity in tidal deposits: A process based stochastic method. In Buccianti et al. (Eds) Proceedings of the 4th Annual Conference of the International Association of Math Geology (IAMG), Ischeia, p.129-134, 1998.

Williams K.E., Madatov A.G. An effective basin model for pore pressure evaluation. Paper presented at American Association of Petroleum Geologists Hedberg Conference, The Netherlands, May 6-9, 2007.

Wong P.M., Boerner S.T. Quantifying uncertainty in statistical regression techniques for reservoir mapping. Paper C042 presented at EAGE 66th Conference & Exhibition, Paris, France, 7-10 June, 2004.

Yielding G., Freeman B., Needham D.T. Quantitative fault seal prediction. The American Assoc. of Petroleum Geologists Bull., v.81, p.897-917, 1997.

Zeltov Y.P. Reservoir Engineering. M., Nedra, 365 p., 1998 (in Russian).

486

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Appendix 1. The multi well data analysis for 3D porosity models initialization and evaluation of lateral

trends in the day content distribution. Real case study

Some results of the case study performed for the huge gas-condensate field in Western Siberia with many years production history is used here to check the possibility of initialization of the 3D porosity distribution and 2D (lateral) clay content trend evaluation based on available standard set of wireline logs.

The 31x12 km sub-area belonging to the upper Cretaceous (Cenomanian) age interval (depth interval 1000-1200m) includes 255 wells with 80 of them available for data analysis (see fig. A1.1). For reasons of confidentiality the real coordinates and names of analyzed wells are hidden on the map and section views here and below.

The reservoir model framework was represented on a fine regular grid by top and bottom depth maps. Its internal heterogeneity (textural model) was recovered based on wireline log correlation between wells only, since there were no seismic data available. Following well logs data were considered for correlation:

1. Gamma Ray (RGK-1 bore);

2. Sonic (BCS bore);

3. Density (DC bore);

4. Resistivity;

5. Gas saturation.

A1.1. Porosity model

The sonic and density logs serve as a main log data input for porosity estimations. The core porosity measurements available on coarse well grid were used for data transformation control.

The sonic data transformation into porosity scale was based on Wyllie’s formula (Ellis, 1987)

ф (вв) = (в - вг)Щ-fy, (A1.1)

where ф is the reservoir porosity estimation; в/ and вг are the transit time (acoustic slowness) in fluid and rock matrix respectively: 190 ps/ft for fresh water; вг ~ 51-55 ps/ft for sandstones; 43.5-46 ps/ft for siltstones (Oil

and Gas..., 1984); в is the measured acoustic transit time.

The validity of Wyllie’s formula is widely varied from lithology type, confining pressure regime and current level of porosity. The sonic-driven porosity for shallow intervals (unconsolidated rock with porosity greater than 30-35 %) tends to overestimate the "true" formation porosity in the Western Siberia conditions. Thus, the relevant empirical corrections were required. The final check was done against available porosity measurements.

The density log (GRR) measures a secondary emitted gamma ray radiation sensitive to the density of electrons in the porous media. This kind of density estimation essentially coincides with a measure of bulk density for most of the rocks. Although some of the rock basis minerals such as coal, evaporite, sylvite etc. cause significant perturbation and relevant bulk density estimations require corrections too (Oil and Gas., 1984). Still, the Wyllie type of formula (A1.1) remains to be robust for GRR log based estimation of the porosity in Western Siberia conditions:

ф (Pf, Pr) = (pGRR - Pr)/(Pf -pr), (A1.2)

where ф is the porosity estimation; pf, pr are the fluid and rock matrix densities respectively; Pggr is the measured GRR density related with bulk porosity (pB = фр/ + (1-фр) as follows: Pgrr = Pb + Vc(pc-pr). Here Vc is the volumetric fraction of rock occupied with clay, evaporite or other distorting mineral which have density pc.

The value of GRR driven porosity requires more corrections at shallow non-consolidated part of the clay reach sections, where the difference (pc - pr) cannot be ignored. The dissolved or second (third) phase oil (gas) may also lower the GRR density values and hence lower relevant porosity estimation. Again, the fluid phase density differences have more contrast at shallow intervals. By analogy with sonic - porosity transformation (A1.1) the final adjustment of pf, pr was done against available porosity measurements.

The results of combined sonic and density transformation into porosity scale for the wells along section line (see fig. A1.1) are displayed in fig. A1.2a and fig. A1.3b (red curves).

i ЩЛ P11025 l о 0 \ \ ЩШ 12 о A °\ ° - 1050 * \ \ ' 1075 “Trl ° \ I3V |\° *b \° \ 1100 a l\o °y °llo °фо 0 1125 <P I \ 2 j \J о 9 о P __ 1150 i / Ш 1175 of ° \ / j|H [m] J 0 <гП° 0 ЯР 7 ° i° ЛШ Ру! 68 0 2.5 5.0 7.5 10 (km)

Fig. A1.1. Field area used for reservoir model initialization with the top depth contours, well positions and section trace indicated. The speculative production well system indicated with underlined name fonts (red -for production wells; blue - for injectors)

487

Depth [km] Ч о .,p Depth [km] Depth [km]

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Fig. A1.2. Log data integration (a) and inter-well correlation of possible flow units (yellow) in section (b) and 3D cube representation (c) Well-log data legend:

Sonic-Density log driven porosity (red); Resistivity (black); Gas Saturation (light blue)

ig. A1.3. Porosity data cube based on mixed 3D Kriging interpolation andXYZ linear trend approximation (a)

and cube data check (b) at control well positions

Control data legend: Sonic-Density log driven porosity - green; data extracted from cube (a) - red

1.05

— 1.10-

Cl

Ф

О

1.15

1.20

0,2 0.3 0,4

Porosity

02 0.3 0,4 C!2_0;3_Q.4

Porosity Porosity

V— '

S b l

0 2 ' 4 Distance [km] 6 8

Fig. A1.4. Fitted to the data Athy's compaction model (c) and relevant synthetic porosity data vs. real porosity data at the control well positions in single well (a) and section (b) views

Data legend:

Porosity data extracted from cube: 3D Kriging interpolation - dark red, mixed 3D Kriging interpolation - dark blue; Synthetic porosity profile based on Athy's model: grin - for poor flow units PFU1-5, dark yellow - for fine flow units FFU1-5 according to the fitted compaction constants FS, CS, FC, CC (c)

488

Oil Production Rate [STBfday] Oil Production Rate fSTB/day1)

Вестник МГТУ, том 13, №3, 2010 г.

стр.457-494

Fig. A1.5. Porosity data cube based on mixed Kriging interpolation (a) and relevant isotropic permeability 3D model (b) against corresponding synthetic porosity data cube (c) and relevant isotropic permeability 3D model (d)

> f WaII PI

'T ■ Computed poros Porosity data cub Mis match

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

/ by m odai e _

./ t * E= a, 062 ,

; ‘

г * : j: * -*

■, - - ... - .(я. РМ-- -

4000

Time [days]

II P2 -

Г 1 - г We

~r

Г“Г - -

1 л

..... -■:.

E = 0.027

i

—*— Porosity data cube Mismatch

1 к—,

0 2000 4000 6DOO 8000

Time [days]

Fig. A1.6. Final mismatch between Well Oil Production Rate curves simulated based on synthetic porosity cube (black) and based on initialized porosity cube (red)

The total misfits (E) are indicated in boxes (see also fig. A1.1 to get positions of the production wells P1-3). See also fig. A1.5

489

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Fig. A1.7. The WOPR simulation errors (value in parentheses) due to uncertainty in initialized porosity model (a) against the WOPR simulation errors (value in parenthesis) due to CMP component variations (b).

The abbreviation for curve code on the upper plot (a):

1. Xy(z) - simulation result based on relevant linear trend approximation implemented at the 3D porosity model initialization;

2. Kr - simulation result based 3D Kriging interpolation at the 3D porosity model initialization;

3. Kr xy(z) - simulation result based 50-50 % mixture between 1 and 2 at the 3D porosity model initialization;

4. Calc - simulation result based on synthetic porosity at the 3D porosity model initialization;

5. MEAN - multi-curve averaging result considered as a reference curve.

The abbreviation for curve code on the lower plot (b):

1. A0+0.1A0 - simulation result based on -10% perturbation of reference model in constant clay content level (component A0 in CMP specification);

2. Ss+0.1Ss - simulation result based on -10% perturbation of reference model in upper rank sensitive specific area constant (component Ss in CMP specification);

3. Calc - simulation result based on synthetic porosity at the 3D porosity model initialization;

4. MEAN - multi-curve averaging result considered as a reference curve.

A total of 53 wells with the best data sets were selected for 3D porosity data interpolation on cube 100x100x100 (see fig. A1.2c; fig. A1.3a). The porosity transformation value related to each specific cube cell was obtained by some smoothing due to coarser z-grid for porosity data cube than for porosity well data. Then the 3D mixed Kriging algorithm was implemented for the whole cube to get the draft 3D porosity model. This algorithm includes W1-W2 mixture between classical 3D Kriging interpolation (Davis, David, 1978) and power xy,z approximation over all included porosity data in the area, where W1, W2 are the weight coefficients adjustable against real data.

The final texture interpretation was based on the draft 3D porosity model and log correlation. It allows allocating of five sub-seismic elements of reservoir texture, each of which reveals relatively higher porosity level and allows merging of relevant depth interval into a simply connected 3D objects (see yellow layers in fig. A1.2b-c; and A1.3b). These objects were recognized as fine flow units (FFU1-5) separated from each other by five poor flow units (PFU1-5) indicated with green colour in the lithology columns, sections and cubes.

490

Вестник МГТУ, том 13, №3, 2010 г. стр.457-494

The feed back control was based on blind tests. Some of the wells with the porosity data available were excluded from the mixed Kriging interpolation and then included back into the project for comparison between "data cube porosity" and well-log driven porosity (see fig. A1.3b).

The blind test reveal that the best result in initialization of 3D porosity cube is based on mixed 3D Kriging and linear porosity trend approximation, where both fractions are taken in equal proportion, i.e. W1=W2=0.5.

We used also Athy's formula (11) to get porosity models at each of the 53 well positions. The compaction coefficients FS, CS, FC, CC were treated as compaction constants specified for each components of RM texture. Then the log-derived porosity data were used as real data at their adjustments during tuning of the relevant compaction models. The results are represented in fig. A1.4.

Therefore, two porosity cubes were initialized on reference fine grid: one - based on real well log driven porosity data; second - based on synthetic porosity model (compare porosity cubes in fig. A1.5a,c). Then two isotropic permeability cubes were computed (see fig. A1.5b,d) based on relevant specification routine (see section 6, formula 6).

Finally two synthetic production histories were simulated: first - based on the real data driven porosity model (fig. A1.5a); second - based on the synthetic porosity model fitted to the real porosity data on control points (fig. A1.5c). The aim of this numerical experiment was to check how big the difference between these two reservoir simulations is.

During simulations the gas-condensate was replaced with the black oil and initial OWC position was tied to the common depth level of 1080m for each flow unit. The scheme for three production wells plus four injection wells was implemented together with three-stage production schedule by analogy with the speculative example (see Appendix 2):

The synthetic well oil production rates (WOPR curves) computed for two reservoir simulations were compared and analysed (see fig. A1.6-7).

The misfit between WOPR curves simulated based on the above porosity initialization modes appears to be reasonably low. In contrast, the model parameter errors converted into the relevant synthetic WOPR misfit turns out much bigger (see fig. A1.7).

The general conclusion based on the given data analysis and simulations could be formulated as follows:

There are two possibilities for 3D porosity model initialization:

1. It can be done based on well log data transformation and later on 3D interpolation (see fig. A1.5a).

2. It can be done based on synthetic porosity model fitted to the real porosity data at control points (fig. A1.5c).

The misfit between relevant synthetic production histories simulated by industry available reservoir simulator is comparable with the tolerant production data error (e ~ 5 %). At the same time the plausible uncertainty in the specific surface (SS, SC) or in the clay content (V) components of the model parameter vector leads to the much bigger errors in reservoir simulation (see fig. A1.7b). Hence, the compaction constants can be fixed and excluded from the data fitting process on condition that there is enough static well log data to initialize consistent 3D porosity model.

A1.2. Clay content model

The gamma ray log measures a primary emitted gamma ray radiation sensitive to the content of clay minerals such as coal, illite, montmorilonite in clastic rocks. Thus the intervals with higher averaged level of gamma ray log are normally associated with greater clay content (Oil and Gas..., 1984). The Wyllie type of formula (A1.1-2) remains to be robust for rough gamma-ray based estimation of clay content in Western Siberia conditions:

V = (Y - Yc) / (Yc - Ys), (A1.3)

where V is the clay content fraction; ys, Ym are the gamma ray log values for pure sand and clay or mud-rock intervals of the host rock. Formula (A1.3) has been calibrated against core description in 24 wells based on analysis of gamma ray data y (GVK1 bore) as follows: ys = 0.95; yc = 16.0.

The results of gamma-ray log transformation into shaliness scale are displayed in fig. A1.9-10. The good correlation between high permeable, gas saturated intervals and relatively low V levels is clear visible on detailed plot (see fig. A1.9a). The relevant section intervals (yellow in fig. A2.10) coincided with "fine flow units" detected by porosity data analysis (see section A1.1). Thus the 3D clay content population inside the reference fine grid reservoir model was built in agreement with above-introduced texture (see fig. A1.2c) by using 3D Kriging procedure.

Fig. A1.8. Explanation of raw clay content фс term after (Ellis, 1987)

491

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

Log data legend:

- Log derived permeability;

- Gas Saturation (scale unit);

- Log derived porosity (scale unit);

- Gamma Ray (Bore GK1);

- GK1 derived clay content after smoothing (scale unit)

- Clay content pseudo log extracted from relevant data cube

Fig. A1.9. Log derived porosity vs. gamma-ray derived clay content data in single well (a) and section (b) representations

The relevant V-cube (see fig. A1.10) was considered as a rough real data driven estimation of reservoir shaliness to be used for investigation of inherent lateral trends.

The following four-step procedure was independently implemented for every fine flow unit for recovering of associated V(xy) trend and analysis of relevant approximation reliability and unit-by-unit consistency:

1. Getting 3D Kriging interpolation of gamma ray log driven clay content values within the given unit.

2. Getting clay content map for the given unit from (1).

3. Getting V(x,y) trend estimation from (2).

4. Estimation and analysis of error from (2-3).

Fig. A1.10. Clay-content cube based on 3D Kriging interpolation of gamma-ray derived 0c-logs

Five units (FFU1-5) included in RM texture model (see fig. A1.4c) were analyzed in agreement with this procedure. The results are represented in fig. A1.11

The kick errors on each differential maps were analyzed in connection with the nearest well locations (see pink coloured indication in right maps in fig. A1.11). Normally, such kick errors were caused by shortage of gamma-ray data in the well.

Finally, the general conclusion based on gamma ray data analysis could be formulated as following.

The clay content area distributions for all textural components of the given reservoir model reveal clear recognizable lateral trends, that could be approximated by linear function V(x,y) with layer averaged accuracy from 0.5 up to 2.5 %. The close to west-east orientation of relevant trends V(x,y) is consistent for all layers except the deepest one (FFU5), which is most likely due to shortage of data available.

492

Вестник МГТУ, том 13, №3, 2010 г. стр. 457-494

•0.48

,0.46

;0.44

:Q 42

0.4

0 38

О 36

0.34

0.32

0.3

0.28

0.26

0.24

0.22

0 2

0.18

|FFU1|

|0.44

0.42

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

0.4

0.38

0.36

0.34

0.32

0.3

0.28

0 26

0.24

0.22

0.2

0.18

FFU2

•0.48

0.46

•0.44

0 42

0.4

0.38

0.36

0.34

0.32

0.3

О 28

0.26

0.24

0.22

0.2

0.18

0.4

0.3©

О 36

0.34

О 32

■0.3

:0_28

■0.26

О 24

0.22

0.2

PL

,0.38

0.36

0.34

0.32

0.3

0.28

0.26

0.24

0.22

0.2

0.18

0.16

FFU3

,0.36 0.34 0.32 0.3 0.28 0.26 0.24 0.22 0.2 0.18

FFU4

,0.34

0.32

0.3

0.23

0.26

0.24

0.22

0.3

0.29

0.28

0.27

0.26

0.25

0.24

0.23

О 14 0 13 0 12 0.11 0 1

9.6-002 sue З.е-002

7.3- 002 б.е-002

Fig. A1.11. The lateral V(x,y) trends extracted from clay-content integration and analysis for reservoir flow units FFU1-4

1. Map of V(x,y) distribution for indicated layer based on Kriging interpolation for layer-related data set (left);

2. Map of V(x,y) background x,y trend based on relevant approximation for the same data set (middle);

3. Map of standard error |Map1 - Map2|/|Map11 (right)

493

Madatov A.G. An effective field scale reservoir model for history matching... (Part 1)

APPENDIX 2. A

single-layer blocked reservoir model architecture and specification

Fig. A2.1. The contour map tied to the reservoir top with oil-water contact level marked (a) and relevant absolute permeability cube (Kxx) with well positions indicated (b)

Cube parameters: Grid (X,Y,Z): 40 40 40

Xmin Xmax [m]: 166247 172648 Ymin Ymax:[m]: 200 6659

Zmin Zmax:[m]: 1000 3600

Well control: Phase 1 (15 time steps):

Production well P1 - open BHP = Production well P2 - shut; Production well P3 - open BHP = Injection wells I1 - shut;

Injection wells I2 - open with rate Injection wells I3 - open with rate Injection wells I4 - open with rate

1020 psia;

= 50000 stb/day; = 50000 stb/day; = 50000 stb/day;

Phase 2 (7 time steps):

Production well P1 - open BHP = Production well P2 - open BHP = Production well P3 - shut; Injection wells I1 - open with rate Injection wells I2 - open with rate Injection wells I3 - open with rate Injection wells I4 - shut;

1020 psia;

= 70000 stb/day; = 70000 stb/day; = 70000 stb/day;

Phase 3 (7 time steps):

Production well P1 - open BHP = Production well P2 - open BHP = Production well P3 - open BHP = Injection wells I1 - shut;

iНе можете найти то, что вам нужно? Попробуйте сервис подбора литературы.

Injection wells I2 - open with rate Injection wells I3 - open with rate Injection wells I4 - open with rate

1020 psia;

= 90000 stb/day; = 90000 stb/day; = 90000 stb/day.

Fig. A2.2. Consequent phases of the synthetic oil production history simulated in agreement with presented well

Fig. A2.3. Synthetic well oil production rates for reference reservoir model simulated in agreement with a three-phase well control schedule (see also fig. A2.1-2). The phase terminations are marked with red arrows

494

An effective field scale reservoir model for history matching and reservoir prediction. Part 1. Model architecture and specification Текст научной статьи по специальности «Науки о Земле и смежные экологические науки»

Аннотация научной статьи по наукам о Земле и смежным экологическим наукам, автор научной работы — Madatov A. G.

Похожие темы научных работ по наукам о Земле и смежным экологическим наукам , автор научной работы — Madatov A. G.

Текст научной работы на тему «An effective field scale reservoir model for history matching and reservoir prediction. Part 1. Model architecture and specification»