The decomposition of 3-D overpressure evolution model in basin scale and its application to the fault seal analysis
A.G. Madatov1, V.-A. I. Sereda2
1 Institute of Marine Geophysics; Higher Mathematics and Computer Science Chair of the MSTU
Abstract. A new 3-D fluid dynamics model of multilayered media has been developed, with complex burial history and capable for calibration. This new method extends the 1.5-D approach to the solution of forward and inverse problems for excess formation pressure evolution at the basin scale, earlier proposed by the authors. The concept of a laterally conductive reservoir carrier among massive permeable rocks has been introduced. The system of sealing faults has been considered as pressure barriers in the relevant lateral drainage system. The vertical and lateral directions of pore fluid migration have been accounted recursively in two relevant mass conservation models. It allows building of an effective numerical scheme capable of modelling formation pressure compartmentalization and calibrating the sealing properties of surrounding faults. We introduce as a goal of calibration the set of independent fault attributes sensitive to the seal properties and sufficient to describe the lateral barrier. A software program package based on this approach has been developed, and some results of its practical application for fault analysis in the North Sea basin have been discussed.
Аннотация. Разработана эволюционная трехмерная модель динамики геофлюидов, отличающаяся от классических бассейновых моделей возможностью калибровки. Метод развивает предложенный ранее авторами подход к решению прямой и обратной задачи флюидодинамики для ряда одномерных, формационно увязываемых моделей генерации и разгрузки сверхгидростатических давлений в шкале времени седиментации. Вводится концепция латерально проводящих каналов - пластов и латеральных экранов - геологических разломов. Для латеральных экранов определяется набор основных эффективных параметров, влияющих на латеральную миграцию флюида в соответствующем латерально проводящем пласте. Уровни эффективных латеральных оттоков увязываются с уровнем эффективной вертикальной миграции порового флюида в процессе уплотнения пород, температурного расширения и генерации геофлюидов. Калибровка единой трехмерной модели генерации и разгрузки АВПД позволяет определять экранирующие свойства разломов, включенных в развивающуюся геофлюидальную систему в качестве латеральных флюидодинамических барьеров. Метод реализован в виде компьютерной технологии. На примере анализа экранирующих свойств разломов в бассейне Северного моря детально обсуждаются элементы и этапы соответствующей обработки. Детали численной реализации алгоритмов моделирования и результаты практической калибровки приведены в приложениях.
1. Introduction
The faults cutting the clastic section interfere with lateral pore fluid movement. Numerous investigations indicate that they can act both as a pore pressure/fluid flow barriers and also as drainage zones depending on relative permeabilities of fault and host rock. In both cases the fault permeability parameter is more uncertain and less predictable due to the more complex origin of this petrophysical characteristic (Manzocchi et al, 1998a). Indeed, rock permeability is mainly a function of compaction state (porosity level) and porous channel connectivity at the micro-scale, whereas fault permeability is dependent on its morphology (length, displacement, deformation zone thickness, juxtaposition, fault trace density, tortuosity, fault segments' linking and orientation), plus host rock properties -porosity, clay content, net-sand connectivity, etc. In addition, the time scale of fluid flow investigation is also important.
For example, production time scale implies reservoir locality and steady state assumptions are taken for the fluid dynamic description. The history of sedimentation and associated faulting is not important here, but fault internal structure should be described in more details. In particular, the fractal nature of the fault-fracture system geometry existing in the present-day reservoir should not be ignored (Childs et al., 1997). The stochastic model for the fluid flow simulation across the fault system can be appropriate and the relevant model parameters will be calibrated. If there are some exploration aims directed for an object, then other fault model parameters appear to be more important. In particular, burial and associated rifting history influences the result of the lateral migration through the carrier beds crossed by faults (Antonellini, Aydin, 1995). Thus, the time scale of fluid
transport model to a certain degree determines the space scale and sensitive parameter set for the fault system to be assigned.
The effect of faulting on fluid flow dynamics in porous carrier beds for both production and exploration time scales is actively discussed in literature (Knott, 1993; Walsh et al., 1998a,b; Manzocchi et al., 1998a,b, etc.). Numerous models have been proposed for different applications: to simulate lateral fluid flow patterns, account for pore pressure compartmentalization and predict fault sealing capacity or HC field integrity. All of them could be treated as fluid dynamic models of the fault zone or, more generally, fault system models. Depending on aims and time scales, the fault zone/system fluid dynamic models could be stochastic - based on homogenization principles (Childs, 1997; Manzocchi et al., 1998a), or empirical - based on "wise old man" approach (Watterson, Walsh, 1994). In both cases a certain model upscaling and calibration are required, i.e. the attributes of the fault system most sensitive to fluid dynamics should be established in order to replace the unpredictable nature of this natural phenomenon with an effective numerical device suitable for reproduction of observed field data. The range of variation and most appropriate values of the fault models in each particular case should be properly calibrated before use. So, collecting new fault models applicable for fluid dynamics in both production and exploration time scale is an endless process unless some common calibration tools can be implemented.
A gap exists between the tools developed for upscaling of the fluid dynamic fault model and its calibration. If the set of fault parameters most sensitive to lateral fluid dynamics are established and even ranged quantitatively, the calibrations are poorly understood and the corresponding approaches are often reduced to simple statistical analysis.
In this paper we focus on developing a simple fluid dynamic model of the fault system appropriate for simulation of pore fluid lateral migration and its associated problems in the basin evolution time scale. In contrast to the existing fluid dynamic models of the faulting, we propose a model aimed for convenient and unique further calibration to be performed. The way of calibration is based on the inverse problem solution theory developed for highly nonlinear forward model operators (Madatov, Sereda, 2000). The input for the calibration should include pore pressure, porosity and formation lithology-stratigraphy data collected from the series of representative wells within the calibration area.
Evidently, there is a certain range of possible fluid dynamic fault models, which could be proposed and calibrated for the same area. Which type is best? What is the optimal parameter set capable of data reproduction and prediction? The answers to these questions lie far beyond the scope of this article. Nevertheless, we maintain that any model type, i.e. a set of fault system parameters and rules for computing of synthetic data, is suitable as soon as it can simulate all of the real data prototypes observed in the offset wells with an accuracy not less than field data quality. So, if the given fluid dynamic model of the fault system is capable of reproducing and predicting the data that have been used for its calibration then such a model is viable and can be used for the chosen application.
From this point of view, the fluid dynamic model of the fault system presented herein can be suggested for the application to pore pressure prediction, compartmentalization, seal analysis and HC field integrity.
2. The 1.5-D pore pressure evolution model with effective lateral outflow carriers
The mass conservation principle for a simple porous rock unit in sedimentation time scale could be reduced to the following differential equation (Madatov et al., 1996):
AcP/dt = divq + cB/ct. (1)
Here the current excess pore pressure (fluid dynamic potential) are time derivatives and current source term - time derivatives as well (¿P/ct and cB/3 respectively); q is the vector of pore fluid flow speed; A is the function of pore fluid and rock matrix compressibility. The terms A and B represent continuously changeable during burial time interval properties of rock and pore fluid vs. current depth and/or temperature. Normally, they are given empirically in the form of "normal trends" or by PVT relationships. Conversion from depth or temperature to time scale is being performed via sedimentation and/or heating rates.
The numerical solution of (1) with respect to excess pore pressure is accompanied with porosity and litho-stratigraphy units' distribution in subsurface.
The fluid flow divergence term could be represented via the Darcy law for the nearly whole section depth except may be upper metres where the active dewatering from nonconsolidated sediments could have quite high rate. On condition of mostly applicable slow drainage rate the base equation (1) can be rewritten as follows:
A cP/ctt = ViCVP) + dB/dt, (1*)
where the pore fluid volumetric flow rate q is represented according to the Darcy law in the form:
q = CVP
and C - is the fluid conductivity function of the simple volume. For single phase flux this function is given by the ratio k/u - rock matrix permeability/pore fluid viscosity.
The instant vertical component of the formation pressure gradient could be reconstructed at any specific time slice of burial history by using back stripping technique. Thus, this information is recoverable from available well data and normal compaction trends established for clastic and carbonate sediments in many regions. In contrast, the lateral component of the formation pressure gradient is less certain at present time and is hardly possible to be reconstructed back to burial history. So, the following simplification of the general model (1*) could be taken in order to upscale the relevant forward model in the way more suitable for the single well data inversion (Madatov et al, 1997).
A cP/dt = d(C3P/âz) /& - SqL + cB/ct.
(3)
«B iriiil Tim e
C
FirP
X ■: ¡i omps ielea ction
£
a QJ
solid phase space 2
1 1 1
The fluid flow divergence term is subdivided here on the explicit vertical drainage term - d(CdP/dz)/dz and the implicit lateral drainage term - SqL. This model has been introduced (Madatov et al., 1996) for purposes of the calibration of the pore pressure evolution model against the single offset well data.
The lateral Sink term could be interpreted so far only effectively. Indeed, it represents the total rate of pore fluid outflow from specific volume averaged for single time integration step and integrated over all lateral azimuths.
Thus, a Steady State assumption is implemented here implicitly for lateral outflow. The main advantage of this model consists in reducing of the corresponding inverse problem down from 3-D to 1-D scale. This approach allows simple and stable calibration of the sensitive formation parameters mostly contributing to the relevant synthetic pore pressure and porosity. As it was shown above (Madatov et al, 1998) from the practical point of view the so called "1.5-D" model (3) is capable to reproduce the excess pore pressure profiles observed in most offset wells. Still it is simple enough for calibrating by using inversion technique (Madatov, Sereda, 1997). A set of calibrated 1.5-D models within the area forms necessary framework model yielded for the quantitative predicting of abnormal pore pressure profile at any target well location.
A numerical scheme developed for solving forward problem (3) implies subdivision of the 1.5-D formation model (well model) onto several formation units having each own time gate and rate of sedimentation as well as lithology type with associated "normal" compaction and fluid conduction relationships (Dickey, 1976). The burial history for this model is subdividing on the corresponding number of sedimentation events - macro steps lasting with the constant sedimentation and lateral drainage rates (Fig. 1). The solution is achieved numerically on the common for all formation units time - a depth grid (Fig. 1a). The depth grid is chosen with the density sufficient to get at least three nodes for the thinnest formation. The terms A,B and C of model (3) are currently re-computing at each time slice for a refreshed depth grid in agreement with the sedimentation rate and "normal compaction" trends predefined for the given formation lithology, PVT properties of the pore fluid and phase conversion potential (HC-generation rate).
Corresponding inverse problem is being
Fig. 1. Burial history sketch for multiformation 1.5-D pore pressure evolution model given in global coordinates (a) and flow rate vector sketch given in local space coordinates coupled with single volume of burying rock (b)
Burial time - T
Depth fn]
rK- Sandstorm, m
"-<- o*, m
J_I
[1.001 0.01 0.1 10 10 101 Fluid conduction [illinn 000 year] at single vertical DC re pressure drop [1 Pa/in]
Depth [rr]
Fig. 2. Relative lateral to vertical pore fluid outflow estimation based on the simple 1.5-D model (2). a -burial history sketch with fluid lateral and vertical outflows through single sand volume marked; b - pore fluid conduction capacities for normally compacted rocks; c - lateral / vertical pore fluid outflow ratio R(D) vs. depth for thin sand layer burying among massive clays
solved with respect to five adjustable parameters most sensitive to the synthetic porosity and pore pressure assigned for every formation unit (Madatov et al, 1998). These parameters include two compaction constants and one fluid conduction constant defining relevant normal trends (Fig. 2b); HC-generation potential included into Source term B and effective lateral Sink term - SqL.
At the calibration of model (3) in the first approximation there aren't any constructive background mechanisms to be implied regarding to the lateral Sink term. So, it's treated as a constant for the given formation unit during every macro step occurred from its deposition till present time. Evidently, the lateral Sink term assigned and calibrated for every formation has an effective (non-physical) meaning, unless the lateral pressure gradients and lateral fluid conductivities are "hidden" in it (see formulae (1) and (3)). Indeed, according to the Darcy law, the term qL combines some lateral draining properties of the given formation with the evolution of the lateral gradient of the fluid dynamic potential.
The simple 1.5-D model represented in Fig. 2a can give some gross estimate of the lateral to vertical outflow rates relationships during normal compaction in the clastic section. The thin sand layer buried inside the thick clay mass with constant sedimentation rate 50m/MY (million of years) is assumed to be laterally continuous and connected to the surface at some remote point (Remote Discharge Zone - RDZ, see also Fig. 3). In this particular model, the lateral offset distance to RDZ is about 200 km. During the compaction history the sand bed acts as a lateral drainage channel for both underlying and overlying clays. Escaping pore fluid under constant overburden increment can use both vertical and lateral directions to provide the normal compaction of the host rock matrix (Fig. 2b). Due to much lower permeability of clay (Fig. 2b) the contribution of the vertical fluid outflow to the total fluid volume escaped through the top of sand will gradually decrease with time and depth (Fig. 2c). Depending on sedimentation rate and relative sand/clay permeability, the lateral to vertical outflow rates ratio will reach some maximum and then will gradually decrease due to reduction of total pore fluid volume to be released.
The contribution of lateral pathways to full pore fluid divergence during compaction in nature cannot be overestimated. This is especially true for the oldest part of the sections overlapped with the thick mudstone and shale formations that are poorly permeable in the vertical direction. According to (England et al., 1987; Bredehoeft et al., 1998), only about 10 % of the total section thickness provides about 90 % of the total volumetric lateral outflow from the whole sediment column during its compaction. Our modelling and calibration results in the North Sea basin also support this conclusion (Madatov, Sereda, 2000).
3. The 1.5-D pore pressure evolution model for single lateral outflow carrier. General approach
The lateral Sink term in (3) can have some "physical" meaning as soon as the lateral pressure gradient has been recovered through the whole burial history. Being linked to the top of some lateral channel, mapped for some area and considered at each time slice of the burial history this effective formation property could deliver valuable information about lateral conduction capacity of the whole section. In particular, the set of calibrated Sink terms associated with the regionally good reservoir could be used after ranging and sorting out for quantitative analysis of the main lateral pressure seals normally connected with the main faults in the area. This becomes possible within the range of some framework lateral draining model established for the calibration area to describe charging and discharging of the fluid dynamic potential versus burial time and space coordinates.
It is important to stress that the full-scale 3-D charging-discharging (CD) forward model inevitably leads to the classical non-resolvable calibration problems associated with inversion of offset data (Madatov, Sereda, 1997). In addition to the natural non-unique and unstable solution of corresponding 3-D inverse problem the complexity of this task will create huge computing problems. Thus, the development of the calibratable effective framework lateral draining model is vitally important for the establishment of a proper link with the available area data. This model can be based on the calibration results obtained against single offset wells. In this case it should imply further incorporation of the calibrated "Sink" parameters SqL with tectonic information from the seismic data. In order to establish such a model some basic assumptions should be made:
1. The predominant direction of the lateral component of pore fluid outflow vector would be possible to evaluate. Generally, this direction coincides with the averaged gradient of total thickness increasing and is preserved through the sedimentation history time gate if the global stretching of the crust didn't change the azimuth. Thus, for the rift basins, like the North Sea, the lateral outflow vector could be taken on average coinciding with the centre to the flanks direction.
2. Any of the lateral drainage carriers selected within the area for calibration have a stationary connection to the surface at a certain remote discharge zone associated with the outer basin edge (Fig. 3).
Fig. 3. Basin wide scheme of the compaction induced pore fluid drainage discharging excess pore pressure. The higher overpressure level the darker is colour indicating the corresponding pressure compartment. Sub-vertical white strips mark sealing faults - lateral pressure barriers
The distance Л between this zone and the centre of the calibration area is much bigger than the area size S. So, the lateral component of the gradient of the current fluid dynamic potential can be easily defined from the current excess pore pressure at any point inside the calibration area. Of course, the local distortions of the lateral patterns must be small enough to be ignored.
3. The space model of lateral fluid outflow is associated with the main drainage carriers and with the main faults crossing the corresponding pathways. This is true for most post-depositional fault systems observed in the sedimentary basins. For example, the North Sea rift basin contains a membrane type of fault seals with permeability 3-5 orders lower than the permeability of the Jurassic sandstone which is the main lateral flow connector during sedimentation (Knott, 1993). At the same time the accommodating mud rocks have even lower average permeability. The typical porosity-permeability ranges according to (Matthai, Roberts, 1996) are as follows:
Sandy carriers: permeability: 10-14 - 10-12 m2; porosity: 0.1 - 0.35
Holding shales: permeability: 10-19 - 10-17 m2; porosity: 0.1 - 0.35
Intersecting fault: permeability: 10-21 - 10-12 m2; porosity: 0.1 - 0.35.
Thus, the lateral drainage from the calibrated area affected by the fault system could be described by some network space model, where lateral carriers and flow connectors are isolated from each other vertically by the poor permeable mud-rocks and laterally by the intersecting faults.
4. Despite that present day fault throw represents the cumulative result of the complex strain-stretch history of the basin evolution, it is hardly possible to recover all the displacement increments related to the corresponding reactivation phases. Faulting is a more catastrophic geological event than compaction-conduction, which can be treated as continuous coupled processes for the given sedimentation period. So, the main faults normally included in seal analysis and calibration should be interpreted as occurring synchronous to the sealing properties, which are not time-dependent and do not change rapidly during the burial history. The start point of the fault system history associated with development of lateral barriers should be a matter of investigation. In particular for the North Sea rift basin, this point could be taken to coincide with a reverse slip that started in Late Cretaceous time (Knott, 1993).
Thus, the connector properties of the network drainage model are controlled mainly by the rapid changing of the lateral carrier conductivity related to the fault activity. So, some trigger mechanism can be applied to describe the evolution of the network model through burial time.
Let us now consider the approach to the fault transmissibility definition and its calibration based on the pore pressure evolution forward problem formulated for a single lateral flow carrier intersected by a single fault.
According to assumption 2 the lateral channel has a stationary connection to the surface in some remote zone outside the calibration area (Fig. 4). Let the geometry and section properties be constant in Y direction through all of the burial history and the set of the single 1.5-D forward models of pore pressure evolution models described by (3), have been calibrated against the relevant offset well data. The position of these wells are highlighted in the X, Z plane by the formation columns.
Evidently, the basic equation (1*) written for sedimentation time scale is still applicable for modelling of pore pressure evolution inside the lateral channel. Moreover, it is convenient also to split the fluid divergence term into lateral and vertical components. In contrast to the 1.5-D vertical model, we will now express the lateral component explicitly and the vertical one via the cumulative vertical Sink term:
SqZAcP/dt = c(ClcP/cL) /cL - SqZ + cB/ct, (4)
where: ¿P/dL is the currently acting lateral component of excess pore pressure gradient; CL - the fluid conductivity of carrier in the lateral direction; A,B have the same meaning as in the previous expressions (1), (3).
Fig. 4. Combination of the 1.5-D vertical and lateral approaches for calibration of the lateral draining CD model for single fluid flow carrier crossed with the fault. a - flow units' sketch; b - three cells' model used for lateral conductivity approximation
It is clear that the vertical component of the fluid flow divergence in this lateral model has no direct relation to the considered formation unit, but it has to be introduced for current account of mass balance. In contrast to 1.5-D vertical case (3), the coefficients of equation (4) belong to the same formation unit except the local fault zone, where CL can be dramatically changed according to the fault zone thickness and permeability (Antonellini, Aydin, 1994). Thus, the space variation of the model parameters is expected to be less abrupt than in the vertical case and the integral conductivity of the lateral channel during the whole post-activation phase (assumption 4) will be mostly dependent on the sealing property of the fault system developed in the area. On the other hand, the pressure data available for calibration of pressure evolution model (4) are normally distributed much less densely. Indeed, the calibration of the 1.5-D vertical CD model (3) is based on formation pressure evaluation in the offset well, whereas the calibration of the space model (4) can only be based on the set of pressure data (usually RFT measurements) coupled with X, Y well coordinates.
The accuracy of the numerical solution of (4) should be in reasonable agreement with the natural variability of the lateral carrier conduction and with the possibility to calibrate it. From this point of view the space grid for the introduced draining 1.5-D lateral model can be much less dense than the one used for the 1.5-D vertical CD model. In particular, the length of the X cell can cover the distance between neighbouring faults along the predominant lateral fluid flow direction (assumption 1) associated with the cross size of the possible pressure compartment. In the case of two adjacent pressure compartments separated by a single fault (see Fig. 4b) the following approximation for the carrier conductivity is valid (Manzocchi et al., 1998b):
Cl1,2 ={m/(Li+L2)[(L1 - e/l)/Kl + 6/Kf + (¿2- 0/2)/Kt\Y
(5)
where: n is the pore fluid viscosity; Ki,2 and KF are the intrinsic reservoir and the apparent fault permeabilities respectively; L1+Ll is the interval between the adjacent pressure compartments including thickness of the fault zone 6. Taking into account, that Ki « Ki,2» KF and L1,L2 » 0, we can reduce (5) into a more homogeneous
form:
C\2 = {^ /(L1+Ll)[1/K + 1/TRF ]} -1,
(5*)
where TRF = 0/KF is the fault transmissibility.
The constructions (5) or (5*) present the obligatory element of the 1.5-D lateral CD modelling. Since it contains the fault seal parameter to be calibrated, we consider it in greater detail in the next section, and here just emphasise that the fault transmissibility explicitly entered into model (4) is the only new model parameter that must be calibrated as a most sensitive to the synthetic pore pressure to establish the goal framework CD model for the calibration area.
The single combination of the consequent solutions of the 1.5-D vertical forward problem (3) and 1.5-D lateral forward problem (4) can be considered as a simplified 3-D solution of the relevant fluid dynamics problem formulated for the same lateral carrier. Let us call the corresponding CD model Combine CD model further in the text. The Combine CD model must be in agreement with the mass balance conditions (1). This could be done in a few loops of the consequent vertical and lateral solutions based on the same set of the adjustable section model parameters (see the scheme in Fig. 6).
A - compartments' sketch;
B - excess pore pressure gradient (vertical component) for deepest, middle and shallowest pressure compartments. The depth positions of the sand are highlighted;
C - lateral-to-vertical components ratio (absolute value) for the fluid flow divergence in the sealed lateral channel (sand) vs. single fault
permeability.
Fig. 5. Results of combine CD pressure evolution and fluid migration modelling for the single lateral carrier burying between massive clays
The simple 2-D model displayed on Fig. 5 illustrates the different capacities of overpressure retention and lateral-to-vertical release fluid flow ratio at different sealing and burial conditions for the single lateral channel (sand) buried inside massive clays. Three hypothetical tectonic blocks in this model (Fig. 5a) are separated from each other by the equally permeable faults. The deepest, the middle and the shallowest blocks are 50, 40 and 30 km respectively distance from the RDZ. The relatively good permeable sand interval is sealed vertically by clay with 3-4 orders lower permeability (see also Fig. 2b) and laterally by sealed fault(s). The sedimentation rate is taken as a constant during the final sedimentation episode for each block. So, the deepest block 1 buried with the highest sedimentation rate (35 m/MY) is connected to the RDZ via three consequent lateral barriers associated with fault 1, 2 and 3. The lateral inflow into this block from the left side is supposed to be equal to zero.
As it was expected the relevant pressure compartment number 1 retains the highest level of excess pore pressure since it has the worst lateral communication capacity and reaches the maximum temperature. Conversely, the shallowest and coolest block 3 has the best lateral communication to the RDZ and close to zero excess pore pressure.
Let us consider now the general scheme suggested for the calibration of the Combine CD model of pressure evolution and fluid migration for the lateral carrier.
As it has been shown earlier, the solution of the 1.5-D lateral forward problem (4) implies a mass balance for pore fluid inflow and outflow at the current pore pressure potential. This balance is computed for each cell of the lateral flow unit at each node of the time grid coinciding with the one used for calibration of the relevant vertical models (3). The inflow charging the system incorporates some portions of the pore fluid released from the pore space during the matrix compaction, thermal extension of pore fluid and solid-to-fluid phase degradation (if there are any). The outflow discharging the system is governed by the Darcy law (2) and depends on the vertical and lateral fluid transport capacity of the rocks-fluid system. Since equations (3) and (4) are applied to modelling the same pressure evolution process for the same area they could be coupled for the intersected cells. Namely, the Source and Sink terms related to the lateral channel could be estimated for each location of the offset well after calibration of the corresponding 1.5-D vertical CD models (3) by the inverse problem solution. Consequently, they can be used as the first approximation
mínffiíllí?*
Fig. 6. Scheme of combining of the 1.5-D CD vertical (above) and lateral (below) inverse problem solutions during calibration of the relevant Combine CD model. p" - Pz/ll5-d - pressure data and 1.5-D models
(Z - vertical, L - lateral) respectively; <5qZ, SqL - vertical and lateral Sink terms; c, f - adjustable formation and fault
parameters; n, m - number of elements for the lateral (fault) and vertical (formations) models; min{ } - operator of non-linear minimization
of the misfit functional
input to define the lateral outflow part of the Sink term more exactly, based on its explicit form (4). In addition, the sealing properties of the lateral barriers related to the faults could be calibrated according to the model assumption assigned above to describe the connector model (5) of the target lateral carrier. The result allows explicit computation of the Sink term for each calibration well. This output in turn can be used as a feedback input for updating the calibrated "Sink" parameters for the relevant interval of the 1.5-D vertical CD model. This looping can proceed until some pre-defined stop criteria are achieved. Finally, the calibration routine for the Combine CD model can be illustrated as shown in Fig. 6.
Note, that both vertical and lateral 1.5-D inverse problem solutions are reduced to fit the synthetic and field pore pressure data collected for the same space coordinates. In general, this fitting problem can be reduced to non-linear minimization of the misfit function given, for example, in the least square norm (Madatov, Sereda, 2000). In general, this problem requires implementation of the optimization technique to be resolved. However, a more simple way could be suggested for the first approximation estimates of fault transmissibility parameters.
4. The fault attributes required for CD modelling and calibration
When talking about seal properties of the fault system it is important to remember, that within the seismically visible space scale (tens of metres of displacement and tens of kilometres of trace length) the only fault attributes that are dominantly sensitive to the observed pressure drop can be calibrated and then used for further prediction. This implies a gross space property upscaling to be done in advance.
It is convenient and commonly accepted to consider fault transmissibility as a combination of the fault zone permeability (a seal membrane characteristic) and disturbed connectivity of the flow unit (carrier juxtaposition). In general, both of these factors develop over geological time and are quite changeable in space.
In addition, 6 and KF attributes of the fault transmissibility (Fig. 5) cannot be treated as completely independent. According to (Manzocchi et al., 1998b) the carrier bed displacement across of the fault is assumed to be the main control of the fault zone thickness, d, and the secondary control on the fault zone permeability, KF. The shale content of the fault zone mainly determines its permeability (Blair, Berryman, 1992).
Using Shale Gouge Ratio (SGR) as a representative measure of the shale content (Manzocchi et al., 1998b) suggested the following empirically based relationship:
Log (KF) = - 4 SGR - 0.25 Log(D )(1-SGR)5, (6)
where: KF is the fault permeability in mD; D is the cross fault displacement in metres; SGR is given as a fraction.
The studies reported by (Walsh et al., 1998b) confirmed that SGR value greater than 0.15-0.2 results in a membrane seal. The fault permeability value obtained from (6) has to be interpreted as a median value of the log-normal permeability distribution covering about two orders. In particular, for D = 20 m, SGR = 0.2, the predicted KF lies between 0.012 and 1.2 mD (Yielding et al., 1997).
Evidently, calibration of the non-linear model (6) applied to the definition TRF from (5*) will require the introduction of one more additional attribute for the fault zone - SGR. The advantage of this approach is that the KF value has here no functional link with the carrier permeability and consequently no systematic errors can follow from the above-calibrated flow unit parameter K.
A more common and simple way is to recognise the permeability of the fault zone as a fraction of the host rock permeability (Antonellini, Aydin, 1994). For example, the following model has been suggested by (Knott, 1993):
KF = KNGR (1.0 - ND), (7)
The NGR (Net to Gross Ratio) here is the lithology mixture characteristic assigned for the lateral carrier. Values can vary from 1 for pure sandstone down to close to 0 for nearly 100 % clay content. ND is the normalised displacement, i.e. ND = (fault throws/lateral carrier thickness). The validity range of (7) is restricted within about 0.99 > ND > 0.01. For faults with throws greater than the reservoir thickness (ND > 1) an appropriate full lateral seal should be assumed. According to (Knott, 1993) 100 % of all faults analysed in the Northern North Sea, 90 % of those analysed in the Southern North Sea (a total of over 400 faults) appear to be sealing.
In general, the fault zone thickness d as a geometry attribute is normally more changeable in space than the membrane seal attribute KF. On the other hand, d has good correlation with the fault displacement, which could be estimated from the seismic data. Sub-seismic fault population modelling (Walsh et al., 1998a) reveals that some power laws or fractal scaling laws on logarithmic scales exist to connect d with the fault displacement and/or with fault segment length. One of the popular empirical models suggested by (Hull, 1988) looks like:
0 = D/x, (8)
where ^is the adjustable coefficient defined in the interval [0.5 - 1.0].
It is important to stress that total lateral outflow from any cell of the carrier fault system is the result of 3-D integration over all possible pathways. In particular, the lateral spill points of the fault zone, acting as the fault terminals, are especially important and can not be ignored in the situation close to the full lateral sealing. Fault thickness could vary markedly either over quite a short distance. So, simply averaging KF and d along the whole fault trace could lead to a severe flaw in the framework lateral draining model. Some segmentation of the slalom fault trace is a necessary element of the fault upscaling. The arithmetic averaging of KF and harmonic average of # within the single fault segment is strongly recommended (Manzocchi et al., 1998b).
Both fault permeability and throws can be treated as a function of geological time. The easiest way to do this is to bind them to the fault displacement history. According to assumption 4 above, the displacement of an existing fault was most likely developed as a series of positive finite increments coinciding with each new reactivation time. The time when each new reactivation occurred is negligible compared to any of the macro time steps involved in the burial history. Substituting the displacement with the corresponding piece-wise time function D(T) into formulae for permeability and throws (see equations 6-8), we can get the time-space model of the lateral seal controlling factors. The uniqueness and accuracy of the fault system calibration will be fully controlled by a number of pressure compartments included in the common CD model, the pressure contrasts between them and the quality of the data required for calibration.
5. The effective CD model and calibration of fault properties for the 3-D multi compartments case
It is clear that the formulation of the fault control problem for the lateral drainage modeling and calibration of the relevant seal properties demands some pore pressure difference observed across the faults in adjacent pore
pressure compartments. Otherwise the pore pressure evolution modeling and inversion problems would not require addressing 3-D fluid dynamics. The lateral variations in fluid conductivity related to lateral facies changes are normally not so abrupt (Blair, Berryman, 1992) and are poorly visible from seismic data. Hence, the framework CD pressure evolution modeling for the calibration area could be based on interpolation of sensitive formation parameters between the wells (Madatov et al., 1998). This allows a reasonable compromise between the data uncertainty and complexity of the model.
Pore pressure compartments in sedimentary basins are associated with blocks of good permeable reservoirs isolated hydraulically from each other. The surrounding faults, which alter the connectivity of the reservoir treated as a flow unit, are usually considered as lateral boundaries of the pressure compartment. A complete lateral isolation during the whole interval of the faulting history is hardly possible because the hydraulic behavior of any fault systems is extremely inhomogeneous (Helth et al., 1994). So, the pressure drop presently existing between adjacent compartments indicates not just "lack of hydraulic connection" (Dickey, 1976), but more precisely, the different rate of the fluid potential discharge via lateral pathways. From this point of view the pore pressure compartment is an excellent natural cell required for numerical simulations of lateral draining within the range of the Combine CD model (3-4). The number of cells, their size and details accounting for describing of their fault boundary are essentially a matter of field data quality and interpretation.
Further we will consider two possibilities suitable for gross and full scale descriptions of pressure compartments required for Combine CD modeling and fault seal analysis.
The first one is based on a steady state assumption about lateral fluid flow during the whole fault activation phase. As shown in item 3 above it can be suggested as a first approximation for rank evaluation of the transmissibility of the main sealing fault within the area. The second one includes the more sophisticated Combine CD model and is reduced to coupled numerical solutions of the relevant forward problem (equations 3 and 4) according to the recursive scheme given in Fig. 6.
In general, there are no restrictions on the details in the description of the fault geometry for the second approach. One should remember, however, that the more blocks and boundary segments that are introduced in the Combine CD model definition, the less unique the inverse problem solution should be expected. Some regularization can only be possible by constraining the variation range for goal transmissibility parameters. Conversely, the steady state approximation provides a unique but quite coarse inverse problem solution with respect to unknown transmissibility parameters in conditions of strict limitations on the number of pressure cells and number of fault segments between the adjacent cells. The steady state approach to fault seal analysis can be recommended as a reasonable start point in further recursive calibration scheme.
In both cases the boundary and initial conditions in terms of pore pressure distributions could be formulated for CD forward modelling of lateral drainage.
Let the calibration area Q include the set of good permeable formations (reservoirs) interpreted by lateral drainage channels and the system of intersecting faults potentially acting as lateral drainage barriers be established in it. Let also the joint interpretation of the reservoir pore pressure data available within the fault system result in a set of pressure compartments (Fig. 7a). Each compartment co is bounded with a finite number of fault segments S and the reservoir thickness can be continuously interpolated to any point inside the given compartment.
Let also the membrane seal quality KF be assigned as the given fault property and the fault zone thickness d vary from segment to segment according to displacement (juxtaposition) variation.
Formally speaking, the frames of the calibration area could also be considered as apparent lateral drainage barriers which include an unknown number of faults located outside the area. Indeed, according to the first model assumption there is a remote discharge zone outside the calibration area. Hence, some fictive compartments have to be attached across the outside boundary of the area along the direction of the expected outflow vector (Fig. 7a). Conversely, some "bottom" compartment has to be introduced to contribute to inflow rate into the CD system Q.. We must also make the following assumptions about the outside Q area continuation for the rest of the compartments:
1. By analogy, to the 1.5-D vertical CD model, we introduce an impermeable lateral barrier between the outside ("bottom") compartment and the adjacent highest pore pressure cell a>i of the system Q.. That is, a zero flow input is expected to be attached through the relevant area edges to the compartment cot. (Fig. 7a).
Pi >P2 >P3>P4>Pi>Ps
— Fault trace Imp emu sib le fictive input fiul ^
- - -: Generalised fictile output fiul
Fig. 7. Six compartment model of the calibration area. a - the pressure compartmentalisation and lateral
outflow ray diagram; b - the spatial model of single pressure compartment; c - the network scheme of the lateral flow connectors
2. By the same analogy, the lowest pore pressure compartment (at on Fig. 7a) has to be connected with some dummy compartment containing remote discharging zone (RDZ). The corresponding part of the external frame (dotted line) could be interpreted as a generalised dummy output fault acting as the last apparent lateral drainage barrier between the area and RDZ.
3. The rest of the pressure compartments in the area are considered as intermediate pressure cells. Since the lateral pressure communications are assumed to be possible only via the lowest pore pressure compartment, the rest of the area frame is considered as transparent for fluid flow, with zero pore pressure gradient across to the given area edge directions. It means that neither inflow nor outflow can be expected via the relevant frame parts.
Let us consider a single activation phase coinciding in time with the relevant macro step used for 1.5-D vertical CD modelling within the area. According to the model assumption 4 it means that the fault transmissibility will be represented by the time-independent functions of KF and 6. It is clear also, that any of the particular solutions achieved for the conjugated 1.5-D vertical CD model calibrated inside the given compartment are applicable on the same time-space grid. So the initial conditions for CD forward modelling of lateral drainage come from the incorporation of individual results of calibration of the vertical 1.5-D models.
Now two possibilities for the next steps in fault transmissibility estimation are considered.
5.1. Calibration in steady state approximation
As shown earlier, the lateral outflow elements of the pore fluid divergence 8qL - lateral Sink term are one of the results of calibration of the 1.5-D vertical CD model. Since it is assumed as constant for a single macro step of fault activation the steady state assumption about lateral flow should be satisfied. Indeed, according to the lateral draining model assumption 1, the location and zero fluid potential of the remote discharge zone are stable in time. On the other hand, the sealing properties of lateral barriers are kept constant during the activation time macro step (assumption 4). So, the differences in pore pressure evolution rates between compartments in first approximation could be compensated for by the differences between Sink terms determined for each relevant compartment. Experience in calibration of the 1.5-D framework models indicates (Madatov, Sereda, 2000) that the effective Sink terms for good permeable formations normally vary by orders from one compartment to another. Normally the corresponding values increase from deep to shallow and stay close to a constant level inside each single compartment (see also Fig. 9). So, the map of lateral Sink terms calibrated in the area and attached to the specific lateral channel can be used directly for gross estimation of the transmissibility of the lateral barriers. Some other simplifications should be added to the above described multi-compartment model in order to do it:
1. The displacement D of the fault zone between i and j compartments should be recovered along each n single segment SriiJ . This can be controlled for a given carrier formation top/bottom from the relevant contour maps.
2. The length of each segment should be constrained by the range of displacement variation along it in order to use the average value of displacement < Dnij > as a representative value to estimate the fault zone thickness < (fij > for a given fault segment. The empirical relationships used in such an estimate should be unique for the whole calibration area and be one of aims of the calibration.
3. Each single fault Fl contains a finite number of the segments N and can cross more than two pressure compartments. The membrane seal parameter KlF is a common attribute of the single fault that must be calibrated in addition to the area-wide empirical relationships between the fault zone thickness and displacement given as an example in (8).
4. The number of unknown transmissibility elements of the fault system within the area includes the number of single faults involved - L plus one area-wide unknown coefficient x in the empirical formula (8). So, the number of pressure compartments must not be less than L+1.
Let us assume the faulting activation phase has begun simultaneously in the sedimentation history of the area for all the faults involved in the calibration and let the corresponding macro time step be equal to AT. The initial (at the start of activation time interval) excess pore pressure P0 and the vertical components of the fluid flow divergence <SqZ> averaged along this macro time step can be determined and attached to each compartment after calibration of the relevant 1.5-D vertical model for area. The present day excess pore pressure PT is assumed to be known in the form of measured values at the reservoir level for each compartment.
We can get the following simplified formula by using steady state assumptions to express the Combine CD model as it follows from equations (3) and (4):
[<A> (Po - Pt) / AT] = [<^> + <q>]„ (9)
where i - is the number of compartment, and <*> - is the time averaging operator.
Note that the real pressure data PT are included in the equation describing balance between the vertical and lateral components of fluid divergence averaged in time (during fault activation interval) and space (for the given compartment). On the conditions introduced above it obtains a unique solution for the inverse problem with respect to the goal properties of lateral sealing faults by using explicit formula for the < SqL> elements. We will show it on the six compartments' example displayed in Fig. 7.
First of all the number of faults to be involved in the membrane seal parameter calibration should be restricted to five. The compartments have to be sorted according to decreasing initial fluid potential in order to form the corresponding connector network model (Fig. 7c).
This network model includes 10 flow connectors (Fig. 7c). The configuration of the flow ray diagram for the given system of compartments and attached initial fluid potentials should not be changed during discharge phase due to steady state assumptions about lateral flow elements.
The following approximation is suggested to describe the average lateral outflow element attached to
any
i-th element of the connector network model introduced in (9):
< S(t > i = Zm=1,M (<qLi,m >/Aim), (10)
where the index m is attached to all M adjacent compartments having fluid connection with i-th one; Xim - the distance between the centres of the i and m cells.
In particular, for the given six-pressure cell model in 10 connector steady state representation the following system of algebraic equations can be written:
<SqL>i = <qL >1,2+ <qL>i,3 <SqL>2 = <qL>2,3+ <qL>2,4+ <qL>2,5 - <qL>1,2
<SqL>3 = <q L >3,4+ <qL>3,6 - <qL>1,3+ <qL>2,3 (10)
<SqL >4 = <q L >4,5- <q L >2, 4" <q L >3,4 <SqL >5 = <qL >5,6 <q L >2,5- <q L >4,5 <SqL >6 = <q L >6,rdz <qL >3,6~ <qL >5,6.
The value of the average lateral flow element <qLi,m > can be expressed based on the Darcy law and definition of the fluid conductivity in the lateral direction (5, 5*):
<qL i,m >=&t,m<Pt,i— Pt,m>/^i,m, (12)
where: <Ptii- Pt,m> is the pressure drop between compartments i and m, averaged over time; Áim is described above; a,m is the fluid conductivity corrected for the fault segmentation as follows:
&t,m ^n=1,N (^i,mC i,m), (13)
where: C"im= (^ [(Ai,m - d"im) / Kt¡i + d nim /KFi,m] / Ai,m}'1 is the fluid conductivity of the n-th fault segment; ^im = S"im / S i - the normalized length of the n-th fault segment; S"ik - the length of the n-th fault segment; S, -the perimeter of the i-th compartment; dnhm and KFm are the thickness and membrane seal of the fault between m and i cells; Ku and ^ - the lateral carrier permeability and fluid viscosity respectively.
Note, that according to assumption (3) the membrane permeability of the fault zone KFim is a common unknown parameter for the given fault. The thickness d nim according to assumptions (2-4) is defined from segment to segment by using measurable average displacement < D"im > and the unknown parameter x of the relevant empirical equation (8) which is unique for the whole area.
5.2. Calibration by non-linear minimization of the misfit functional
It is clear that the draining model introduced above for the calibration area Í2 in steady state approximation is quite gross in terms of pressure cells and fault segment limitations. From this point of view it is only recommended as a first approximation for the full-scale lateral drainage modelling and fault seal properties' calibration. A more general approach implies more detailed segmentation of the chosen faults. As soon as the number of segments with unknown seal properties exceeds the number of the pressure cells where the direct pressure measurement data are available, even steady-state simplification won't provide a unique solution for the relevant inverse problem. In such a situation it is more appropriate to apply 3-D CD modelling of the excess pore pressure field evolution and tuning of membrane permeabilities to provide the minimum misfit between the synthetic and real data for all available calibration well locations.
The numerical scheme for 1.5-D lateral forward modelling with the effective implicit scheme is given in detail in Appendix 1, where this problem is considered in a two-dimensional statement, i.e. represented as
more possible from the computing point of view. Here we introduce the goal function of the relevant fitting problem (Fig. 6) in the operator form and discuss possibilities for its minimization. Essentially, the non-linear optimization philosophy developed for the inverse problem solution in the 1.5-D vertical CD model definition remains the same, and more recent works (Madatov, Sereda, 2000) can be recommended for more details.
Let the detailed segmentation of the fault system be done based on tectonics and pressure data and steady state approximation of the compartment model for the target carrier. Let also the synthetic pore pressure field be calculated at any i,j node of a regular 2-D grid that is tied to the top surface of the target carrier and fully covers the calibration area Q (see Appendix 1). The relevant synthetic pore pressure value can be represented in the operator form as follows:
p% = Pl15D [Sqz; fi,f2,...,fn)], (14)
where: PLL5~D[ ] is the operator of pore pressure forward modelling in the 1.5-D lateral definition (4); SqZ is the vertical Sink term available at each control node after the single well model calibration; f1,f2,...,fn are the adjustable seal parameters (permeability) of the fault segments; n is the total number of the fault segments involved in the Combine CD model of the target carrier.
Let the pressure data available at M offset wells in the calibration area for the target carrier interval be estimated on average at each of the offset well positions and tied to the relevant top surface. The offset well coordinates can be upscaled within the 2-D grid system introduced for 1.5-D lateral modelling. It allows us to form a set of control nodes - 1,2,...M, where comparison with the synthetic pressure is possible, and to define on it the real data vector p*={ p*h p*2...p*M}.
The subset of synthetic pore pressure values collected and sorted out for control nodes forms the model vector pA={ pA1: pA2, ...pAM} suitable for estimation of misfit with the real data. The aim of calibration at this stage is minimizing if this misfit by tuning of adjustable seal parameters f1,f2,...,fn of the fault segments.
In order to formulate the relevant fitting problem generally we introduce the dimensionless misfit vector function R(f) in the form:
R(f) = W-[p* -pA(f)] / ||p*||, (15)
where: W is the diagonal (MxM) matrix of weight coefficients, regulating confidence for the pressure data quality; ||p*|| - the Euclid vector norm for the data vector.
Then the fitting problem can be formulated as follows:
f = argmin (||R(f)||) and f e FacF. (16)
Here Fa denotes the a priori set being searched within the model parameter hyperspace F, which represents the full range of possible variations for goal fault seal parameters (see above).
Since the forward modelling operator PLL5~D[ ] is highly non-linear with respect to the vector f, the analytical solution of the fitting problem (17) cannot be designed. By analogy to the solution of the inverse problem at the calibration of the 1.5-D vertical model we suggest the application of the combined gradient (steepest descent) and the Hauss-Newton search strategies (Madatov et al, 1998) for solving the given fitting problem.
6. The real data test
The real data test was based on the high-pressure region located in the northern part of the Viking Graben (the North Sea), comprising parts of both east and west margins with some extension to the Horda platform. The region has been actively explored over the last 20 years. So, the geology and tectonics of this basin are widely described in literature (see, for example, Introduction to the petrolium..., 1984; Thorne, Watts, 1989; etc.). The abrupt changes of pore pressure across the faults at the Jurassic-Triassic reservoir levels are commonly experienced in exploration drilling. That phenomenon is commonly interpreted within the range of the set of the hydraulically isolated pressure compartments with different drainage regimes during the Tertiary phase of the sedimentation history (Chiarelli, Duffaud, 1984). From this point of view the area can be treated as an ideal natural polygon for testing Combine CD modelling and inversion approach to the fault seal properties estimation.
6.1. Geology, tectonics and overpressure review
The sedimentary rocks in the Viking Graben form two megasequences separated by the Base Cretaceous unconformity.
The Jurassic-Triassic sediments (below Base Cretaceous unconformity) have a complicated tecto-sedimentologic history, resulting in a highly variable lithology where the major unconformities are presented at several stratigraphic levels. Their structural complexity and sealing conditions are controlled by faults with downward rotation, listric geometry and salt tectonics. The fault tectonic dominates in the present day 3-D structure image of the
Jurassic-Triassic part of the section. According to (Johnson, Stewart, 1985) the major fault movements in the area are post-depositional and can be dated from the middle of the late Tertiary time interval (up to 64 m. years bpt.)
The present day excess hydrostatic fluid pressure distribution at the Base Cretaceous level and below reveals several compartments each with its own hydrodynamic condition. The main oil and gas targets are linked to the Middle Jurassic and Triassic sediment groups. The closer these targets are to the pre-axial part of the Graben the higher the risk of overpressure occurrence at their levels.
The Tertiary part of the section is a relatively homogeneous sequence with high content of clay and silt. There are no significant accumulations of hydrocarbons in the Tertiary because of the absence of the high permeable channels for hydrocarbon migration and the lack of proper reservoir beds. So, the main lithological and petrophysical features of the formations are more constant and correlated with the formation level.
The Cretaceous part of the section is extremely changeable in thickness from the margin to the central part, and the thickness of the near-axial Cretaceous sediments exceeds 1000 m. Here they are represented mainly by very compacted shales, which can act as regional sealing barriers for vertical fluid outflow from the deeper part of the section.
Finally, the main geological features having impact on the present day pore pressure distribution could be listed as follows:
1. Fault-block tectonics on the rift megasequence level, which can create significant lateral seals between different pressure compartments.
2. Two different types of subsidence during rifting and post-rifting time which caused quite isolated pressure regimes for upper and deeper Cretaceous unconformity megasequences.
3. Several periods of non-deposition or uplift and erosion during Jurassic-Early Cretaceous time which could cause simultaneous local release of the overpressured pore fluid.
4. Higher temperature gradient in the pre-axial part of the Graben due to different sedimentation and cooling rates. Consequently different start times for the maturation and HC-generation from source rocks and different phase content of the pore fluid observed in different pressure compartments at the present day.
As a result the mechanisms responsible for the pressure development vary from one pressure compartment to another and from the upper part to the deeper part in the section.
More detailed review of the overpressuring mechanisms acting in the region is available in the published literature (see, for example, Chiarelli, Duffaud, 1984; Buhrig, 1989). The following is a summary of the most important points:
The lateral conduction of sealing faults and gas impulsion into highly permeable sandstone channels are dominant among overpressuring factors in the Jurassic-Triassic parts of the megasequence. The effects of temperature and gas generation are crucial in closed compartments located in the axis part of the Graben, i.e. to the south-west. The overpressuring here is mainly restricted by the fracture limits of the rocks or the pressure release capacities of the adjacent overburden faults acting as a pressure valve. In contrast to this, the pressure compartments open to the East (belonging to the Horda platform) show almost hydrostatic pore pressures at the reservoir level.
The origin of potential overpressuring in the Tertiary part of the section is proven to be related to the high sedimentation rate at the late Miocene-Pleistocene time, including glaciation. It could create undercompaction and overpressuring for the poorly permeable Green Eocene clay.
Fig. 8. Geology and pressure data input.
Structural map (a) and 3-D surface view (b) for the target carrier top with the fault traces and calibration well marked. The filled squares (a) and triangles (b) highlight offset wells involved in both vertical and lateral loops of the CD model calibration. (c) - section through the offset wells with pore pressure data attached. The RFT data available at the target carrier intervals were involved in both vertical and lateral loops of the CD model calibration traces are highlighted with bold lines in Fig. 8.
According to the general scheme (Fig. 6) the 1.5-D vertical CD models were introduced for the 19 offset wells in the area at the first phase of the pressure study. This model was based on a formation column which was common for the area. The 14 formation units were considered to be sufficient to describe all the pressure variations observed in the wells; for the Tertiary, Cretaceous and Jurassic--Triassic sections comprising 4, 3 and 7 formations respectively. Each of the individual 1.5-D vertical models was calibrated based on the 1.5-D inverse problem solutions. Finally, the appropriate framework model for the area was established after the model parameter regularization (Madatov, Sereda, 2000).
The calibrated properties of the Cretaceous-Tertiary and Triassic-Lower Jurassic intervals of the section allowed the calculation of the vertical component of the pore fluid divergence (vertical Sink term) for the target lateral channel. The lateral barriers in pressure communication become visible within the lateral Sink term map (Fig. 9a). The six pressure compartments were confirmed for the fluid dynamic unit of the target carrier. The lateral and vertical Sink term map reveals agreement with the depth and pore pressure discontinuity across the relevant fault boundary of the compartments (Fig. 9). In particular, the low average levels of the lateral fluid divergence component inherent to highly overpressured compartments 1 through 3 were in good correlation with the
From about 2800-3000 m depth the quality of the vertical conduction of Cretaceous sediment could decay dramatically due to diagenetic porosity collapse (Byrd, 1990). In some places this process led to high overpressure values through the whole interval of Cretaceous formations. In the more easternmost areas some strings of fractured carbonates and siltstones possibly provided a lateral drainage for the releasing pore fluids, which can explain the near-hydrostatic level for the same formation in the neighbouring.
6.2. Establishing and calibration of the Combine 3-D CD model in the area
The commonly good permeable Middle Jurassic oil and gas reservoir was used as a main lateral draining carrier to establish the Combine CD model and to calibrate the sealing properties of the relevant fault system.
The structure map of the reservoir top reveals several terraces separated by steeply dipping faults (Fig. 8). Some of these faults could be recognized as sealing faults by joint interpretation of subsurface topography and actual formation pressure data available for the target carrier. The relevant fault
Fig. 9. Lateral (a) and vertical (b) "Sink" parameter maps for target carrier after calibration of 1.5-D framework model. The relevant variation range is given in Log scale
1.0 1 2i 1.5 175 2/J Excess pore pressure scale [iv'co] | | | | | ]-^-1
Fig. 10. Combine 3-D CD modelling of pore pressure evolution for 6-compartment model of the lateral draining system of target carrier.
Sealing faults are marked with striped lines; squares show offset wells; current time slice of the numerical solution is shown at the low left corner in million years bpt.
onent (sedimentation rate) and vice versa.
H.
fir
The 3-D distributions of the top depth, vertical Sink term, matrix compressibility and fluid conductivity were available for the target carrier interval at any paleo time slice during the whole fault activation interval. Since their variations vs. time were less than the limits of numerical errors for modelling the relevant average values were taken as an input for the second phase of calibration. According to the general scheme, this phase included the fitting of the fault seal properties to the lateral draining model.
The pressure data involved in calibration were collected from the six offset wells representing each pressure compartment (Fig. 8). The fault system included 5 faults subdivided into 20 segments with known carrier juxtaposition. So, the number of adjustable independent fault seal properties was six - five unknown membrane seal parameters for five independent faults (Fig. 11a) plus one unknown empirical constant in displacement vs. fault thickness relationship (8). Based on formula (6) that allowed gross estimation of the fault permeabilities segment by segment in steady state approximation (9-11). More rigorous fitting was based on the 3-D CD numerical modelling performed on the regular grid by using the implicit scheme described above. At this stage the segment permeabilities were adjusted locally to their steady state approximation. The fitting was stopped and relevant permeabilities were fixed as the stop criteria were reached simultaneously at each offset well.
1[J
it
mW.
E F
I
-1 -2 -3 -4 -5 -6 -7 -8 -3
Fault seal permeability [mp] (Log
Fig. 11. Pressure compartment scheme
(a) and steady state approximation of the lateral fluid drainage
(b) during time interval of the CD modelling for target carrier. Numbers of compartments are shown inside the white circle. Sealing faults are marked with capital letters
The fault system established for the area by geometry of the segments and their calibrated permeabilities can be used to account exactly for the lateral Sink terms at any 1.5-D model calibration attached to the offset wells. Moreover, this formation parameter for the given target carrier was included in the next loop of the Combine CD model calibration (Fig. 6) as a constant for each offset well. It allowed more exact tuning of the rest of the adjustable formation parameters against single well data. Updated vertical Sink term, matrix compressibility and fluid conductivity for the target carrier interval enter again into the next iteration of lateral drainage modelling.
The looping process (Fig. 6) was continued until the data stop criteria were satisfied for both vertical and lateral CD models. Some fragments of the Combine 3-D forward problem solution for excess pressure evolution in the six-compartment drainage model of the target carrier are represented in Fig. 10. Note that this modelling is based on calibrated segment-by-segment permeabilities specified for all of the faults involved in the 3-D framework model. The details of the calibration are available in App. 2.
7. Discussion
The basin modelling approach to fault seal analysis implies a clear understanding of the burial and fault activation history in addition to general knowledge of local geology, lithology and tectonics. It is likely that different tectonic and historic frameworks established for the calibration area can lead to different fault segmentation and permeability estimations, in which case some effective Combine CD model parameters will be obtained after converging loops of calibrations. From a pragmatic point of view any model capable of calibration within a physically reasonable range of model parameter variations is viable and can be used for prediction. Of course, the larger the systematic errors in the description of local tectonics, lithology, stratigraphy, etc., the more remote the target model parameter values are from the real ones. In particular, certain historical and structural framework provide for realistic estimation of fault seal properties and allows for interpretation of 3-D Combine CD modelling for the given lateral carrier in terms of the fluid migration system analysis.
From this point of view the calibration of the Combine CD model of excess pressure evolution inside the target carrier allows the following reconstruction to be done for this flow unit evolution within the basin model scale:
1. Six existing pore pressure compartments had full pressure communication at the start of fault activation history. At those times (early Paleocene) the fluid dynamic potential of the target flow unit was equal to zero, i.e. pore pressure was close to hydrostatic level everywhere within the area.
2. Different rates of sedimentation in the centre (compartments 1, 2 and 3) and at the flanks (compartments 5 and 6) resulted in the generation of a lateral pore pressure gradient. The biggest fault in terms of average displacement, D-D', acted as a main lateral seal on the drainage pathway. Other faults serve as local lateral barriers and cause redistribution of fluid dynamic potentials in the central and flank parts of the area.
3. The highest lateral pressure gradient is kept across the main sealing fault D-D' through the whole Tertiary time interval. So, the main lateral flow pattern is focused on the weakest section of the main seal, located in the southern part of the area (see lateral flow, Fig. 11b). This relatively narrow pattern stayed approximately constant and provided the main volume of pore fluid migration during the whole time of the fault activation and cell-by-cell pore pressure distribution.
Authors have pleasure to thank Eamonn F. Doyle for many years' fruitful collaboration in the pore pressure prediction field. In particularly, authors are grateful for his close attention to this article and some wise notices and practical help that Eamonn Doyle kept provided us. Authors also like to thank Hans B. Helle for his fantastic efforts in coordination of joint research with Fault Analysis Group headed by John Walsh (Liverpool University, UK).
References
Antonellini M., Aydin A. Effect of faulting on fluid flow in porous sandstones: Petrophysical properties. AAPG
Bulletin, v.78, N 3, p.355-377, 1994. Antonellini M., Aydin A. Effect of faulting on fluid flow in porous sandstones: Geometry and spatial
distribution. AAPG Bull, v.79, N 5, p.642-671, 1995. Blair S.C., Berryman J.G. Permeability and relative permeability in rocks. Fault Mechanics and Transport
Properties of Rocks. ISBN 012-243780-2, p.169-186, 1992. Bredehoeft J.D., Djevanshir R.D., Belitz K.R. Lateral fluid flow in compacting sand-shale sequence: South
Caspian Basin. The American Association of Petroleum Geologists Bulletin, v.72, N 4, p.416-424, 1998. Buhrig C. Geopressured Jurassic reservoirs in the Viking Graben: Modelling and geological significance.
Marine and Petroleum Geology, v.6, p.35-48, 1989. Byrd W.D. Geology of the Ekofisk Field, Offshore Norway. Geological Society Spec. Publ., N 55, p.446-451, 1990.
Chiarelli A., Duffaud F. Pressure origin and distribution in Jurassic of Viking basin (United Kingdom -
Norway). AAPG Bulletin, v.64, N 8, p.1245-1266, 1984. Childs C., Walsh J.J., Watterson J. Complexity in fault zone structure and implications for fault seal
prediction. NPF Special Publication, v.7, p.61-72, 1997. Dennis J.E., Scnabel R.B. Numerical methods for unconstrained optimization and nonlinear equations.
Prentice-Hall, Inc., Englewood Cliffs, New Jersey, 440 p., 1983. Dickey P.A. Abnormal formation pressure. Discussions. AAPG Bulletin, v.60, N 7, p.355-377, 1976. England W.A., Mackenzie A.S., Mann D.M., Quigley T.M. The movement and entrapment of petroleum
fluids in subsurface. Journal of the Geological Society, London, v.144, p.327-347, 1987. Helth A.E., Walsh J.J., Watterson J. Estimation of effects of sub-seismic sealing faults on effective permeabilities in sandstone reservoirs. North Sea and Gas Reservoirs. Norwegian Institute of Technology (NTH), p.327-347, 1994.
Hull J. Thickness-displacement relationships for deformation zones. Journal of Structural Geology, v.10, p.431-435, 1988.
Hunt J.M. Generation and migration of petroleum from abnormally pressured field compartments. AAPG Bul., v.74, p.1-12, 1990.
Introduction to the petroleum geology of the North Sea. Edited by K.W. Glenkie, Oxford, London, 236 p., 1984. 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. Knott S.D. Fault seal analysis in the North Sea. AAPG Bulletin, v.77, N 5, p.778-792, 1993. Madatov A.G., Sereda A.-V.I. Computing aspects of the solution of some inverse problems in geophysical investigations. International Conference "Northern Universities" (Oct. 16-18, 1997). Abstracts, Murmansk, p.57-58, 1997.
Madatov A.G., Sereda V.-A.I. The forward and inverse problems of the fluid dynamics in basin modeling applied to the pore pressure prediction within the sedimentary basins. P.1 Theory aspects. Proceed. of the MSTU, v.3, N 1, p.89-114. 2000.
Madatov A.G., Sereda A.-V.I., Doyle E.F., Helle H.B. The "1.5-D" inversion approach to pore pressure evaluation. Concept and application. Paper presented on workshop "Compaction and overpressure current research", 9-10 December 1996, Institute Francais du Petrole, Paris, France, 1996. Madatov A.G., Sereda V.-A.I., Doyle E.F. Integration and inversion of well data into the basin evolution model: Way to the new generation of pressure prediction technologies. Paper presented on American Association of Drilling
Engineers Forum "Pressure regimes in sedimentary basins and their prediction", 2-4 September 1998, Houston, USA, 1998.
Madatov A.G., Sereda V.-A.I., Doyle E.F. Pore pressure prediction by using inversion before and during drilling. Paper presented on workshop "New methods and technologies in petroleum geology, drilling and reservoir engineering", 19-20 June 1997, Krakow, Poland, 1997. Manzocchi T., Ringrose P.S., Underhill J.R. Flow through fault systems in high-porosity sandstones.
Geological Society, London, Special Publication, v.127, p.65-82, 1998a. Manzocchi T., Walsh J.J., Nell P., Yielding G. Fault transmissibility multipliers for flow simulation models.
Fault Analysis Group. The University of Liverpool Department of Earth Sciences, p.1-37, 1998b. Matthai S.K., Roberts S.G. The influence of fault permeability on single-phase fluid flow near fault-sand intersections: Results from steady-state high-resolution models of pressure-driven fluid flow. AAPG Bull., v.80, N 11, p.1763-1779,
1996.
Richards Ph.C. The early to mid-Jurassic evolution of the northern North Sea. Oil and Gas Reserves,
Geological Society Special Publication, N 55, p.191-205, 1993. Spencer A.M., Larsen V.B. Fault traps in the Northern North Sea. Oil and Gas Reserves, Geological Society
Special Publication, N 55, p.281-298, 1993. Thorne J.A., Watts A.B. Quantitative analysis of North Sea subsidence. AAPG Bul., v.73, N 1, p.88-116, 1989. Walsh J.J., Watterson J., Heath A., Childs C. Representation and scaling of faults in fluid flow models.
Petroleum Geoscience, v.10, p.1-11, 1998a. 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, Spec. Publ, N 127, p.99-114, 1998b. Watterson J., Walsh J. Fault seal prediction: Impossible, or just very difficult? Geological Society Special
Publication, N 18, p.219-231, 1994. Yielding G., Freeman B., Needham D.T. Quantitative fault seal prediction. AAPG Bulletin, v.81, p.897-917,
1997.
Самарский A.A., Гулин A.B. Численные методы. М., Наука, 430 е., 1989.
APPENDIX 1. The numerical scheme of the 2-D pore pressure forward problem solution on regular grid
The forward problem for the excess pore pressure evolution inside the lateral drainage system Q was formulated in (4). For numerical solution reasons this problem can be represented by the following differential equation written for unknown P(x,y,t) distribution.
A(x,y)dP/dt = d(kx(x,y)(cPlcX))lcX + d(ky{x,y){dP/dy))/dy + f(x,y,t). (A.1)
Here the time argument t scans through the burial history interval [0,7] and the spatial arguments x and y are snapped to the finite rectangular grid: xe{0,Lx], ye[0,Ly] stretched over the calibration area. The time averaged functions A(x,y), kx(x,y) and kH(x,y) represent combined rock matrix and pore fluid compressibility and x,y component of the fluid conductivity terms respectively. The total source functionfx,y,t) determines the rest of the pore fluid injection term due to vertical migration and inside cell generation into the drainage system. All these terms should be available as a result of calibration of the 1.5-D vertical CD models against offset wells and creation of the framework model within the area limits - {Lx,Ly}. The initial distribution of excess pore pressure can be available by the same means. Thus, the initial conditions for the problem (A.1) come from calibration and further interpolation for the rectangular grid:
P (x,y,0) = Po(x,y). (A.2)
The boundary conditions were introduced based on three assumptions about continuation outside of grid frame Q (see item 5). In a more formal way they could be formulated as follows:
cP(x,y,t)/cx = 0, cP(x,y,t)/cy = 0, V(x,y) e Gb
cP(x,y,t)/cx = P(x,y,t)/A, cP(x,y,t)/cy = 0, V(x,y)e Go (A.2)
dqx/dx + dqy/dy = 0 V(x,y)eG,-
Here Gb, Go, Gt - represent the "bottom", "outside" and "intermediate" parts of the Q edge respectively (Fig. 7).
The "bottom" condition concerns the zero flow regime across the deepest compartment edge. The "outside" condition fixes the lateral pressure gradient in the x direction by the small fraction of the current pore pressure value achieved at the highest compartment edge. Finally, the last condition is applicable for the rest of compartments - the "intermediate" parts of the Q area edge. It implies that the lateral components of the fluid flow divergence are
constantly equal to zero for the elements of this edge, i.e. there is neither increase nor decrease of the pore fluid mass at the single cells adjacent to G1.
In general, the boundary conditions can be represented continuously along of Q area edge by the function regulating the specific value of x, y components of the pore fluid divergence.
In order to construct the final difference numerical scheme let us introduce the time grid ®t with the constant simple step ht as follows:
a, = (tn = nht, n=0,1,2,...,N-1,N), Nht = T. (A.4)
Let us also introduce the space grid along x and y with the constant steps hx and hy respectively:
A,y=((x!,yj)/x=ihx, i=0,1,2,...,Nx-1;yj=/hy, /=0,1,2,.,.,Ny-1}, hxNx=Lx, hyNy=Ly. (A.5)
Let <sx,y be the set of internal nodes of the spatial Qj grid and ^xy= tfxy u ^°x,y u the set of edge nodes related to the corresponding edge elements Gb, Go, Gi.
We suggest applying the implicit scheme to get the numerical solution of the forward problem (A.1). The final difference approximation of this scheme can be written as follows:
Ani,i (Pn+1M - PnM)/ht = (l/hx){o"x,i+i (Pn+11+i,1 - Pn+11,1)/hx -ox,i (Pn+1M - Pn+11.i,1)/hx}+
+(l/hy){ oj (Pn+11,1+1 - Pn+11,1)/hx -oyj (Pn+11,1 - Pn+11,1-1)/hy} + fn1j ;
(x^j) e ax,y ; n=0,1,2,.. .,N-1;
, n=0,1,2,N-1 ;
P01,1 = Po(x1 ,yj), (X1,yj) E^jy ; Pn+11,1 = P1(x1 ,yj,tn+1), (x1,yj) e^x,y
(Pn+11,1 - Pn+11-1,1)/hx= P(x,y,t)/A; (Pn+11,1 - Pn+11,1-1)/hy= 0, (X1,yj) e v °x,y °n+1 ■ Pn+11-1,1)/hx = (Р°1,1 - Pn1-1,1)/hx; (Pn+11,1 - Pn+11,1-1)/hy = (Р°1,1 - Pn1,1-1)/h.
n=0,1,2,N-1.
(Pn
(X1,yj) e 77 "x,y , n=0,1,2,N-1.
The terms ox,1 and cryj 1n (A.6) are def1ned by:
o-x,1= (kx(x1-1 ,yj)+ kx(x1 ,yj))/2;
cry- (ky(x1 ,уя)+ ky(x1 y))/2.
(A.6)
(A.7) (A.8) (A.9) (A.10)
(A. 11)
Fig. A.1. The band type matrix structure of the target system (A.8) for the Nx= Ny = 5. The nonzero elements are filled
It is well known that the implicit scheme (A.6-11) is stable at the arbitrary ratios for the single steps ht, hx, hy (Самарский, Гулин, 1989). Still, the numerical solution of the system (A.8) requires a special technique to be developed in order to handle the large dimension of the practically important space - ax,y and time - сц grids.
This technique is reduced to using the Pieceman-Rachford final difference scheme of the variable directions method (Самарский, Гулин, 1989) with some adoption to the specific structure of the relevant matrix (Fig. A.1). This five-diagonal structure of the target system matrix can be easily created by the following reorganization of indexes within the space grid ax,y:
k = (z-1)(Ny-1)+/, /=1,2,...,Nx-1, /=1,2,..., Ny-1,
1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 3 7 4 4 414 0 0 0 0 0 0 0 0
1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 3 7 4 4 414 6 6 6 6 6 6 6 0
1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 3 7 4 4 414 6 6 6 6 6 6 6 0
1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 3 7 4 4 414 6 6 6 6 6 6 6 0
1 1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 7 4 4 414 6 6 6 6 6 6 6 0
1 1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 7 4 4 414 6 6 6 6 6 6 6 0
1 1 1 1 1 1 1 3 3 3 3 3 3 3 3 3 3 7 4 4 8 515 6 6 6 6 6 6 0
1 1 1 1 1 5 2 2 3 3 3 3 3 3 3 3 3 7 4 4 8 5 515 6 6 6 6 6 0
1 1 1 1 1 5 2 2 3 3 3 3 3 3 3 3 3 7 4 4 8 5 5 515 6 6 6 6 0
1 1 1 1 1 5 2 2 2 3 3 3 3 3 3 3 3 7 4 4 8 5 5 516 6 6 6 6 0
1 1 1 1 1 5 2 2 2 3 3 3 3 3 3 3 3 7 4 4 8 5 5 516 6 6 6 6 0
1 1 1 1 5 2 2 2 2 3 3 3 3 3 3 3 3 7 4 4 8 5 5 516 6 6 6 6 0
1 1 1 1 5 2 2 2 2 3 3 3 3 3 3 3 3 3 7 4 8 5 5 516 6 6 6 6 0
1 1 1 1 5 2 2 2 2 3 3 3 3 3 3 3 3 3 3 9 5 5 5 516 6 6 6 6 0
1 1 1 1 5 2 2 2 2 3 3 3 3 3 3 3 3 3 3 9 5 5 5 516 6 6 6 6 0
1 1 1 6 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 9 5 5 5 516 6 6 6 6 0
1 1 1 6 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 9 5 5 5 516 6 6 6 6 0
1 1 1 6 2 2 2 2 2 3 3 3 3 3 3 3 3 3 9 5 5 5 5 5 517 6 6 6 0
1 1 1 6 2 2 2 2 2 3 3 3 3 3 3 3 3 3 9 5 5 5 5 5 517 6 6 6 0
1 1 6 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 9 5 5 5 5 5 517 6 6 6 0
1 1 6 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 9 5 5 5 5 5 5 518 6 6 0
1 1 6 2 2 2 2 2 2 3 3 3 3 3 3 3 310 5 5 5 5 5 5 5 518 6 6 0
919 2 2 2 2 2 2 2 3 3 3 3 3 3 3 310 5 5 5 5 5 5 5 518 6 6 0
2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 310 5 5 5 5 5 5 5 5 5 518 60
2 2 2 2 2 2 2 2 2 2 4 3 3 3 311 5 5 5 5 5 5 5 5 5 5 518 60
2 2 2 2 2 2 2 2 2 2 4 3 3 312 5 5 5 5 5 5 5 5 5 5 5 518 60
2 2 2 2 2 2 2 2 2 2 4 3 312 5 5 5 5 5 5 5 5 5 5 5 5 518 60
2 2 2 2 2 2 2 2 2 2 4 3 312 5 5 5 5 5 5 5 5 5 5 5 5 518 60
2 2 2 2 2 2 2 2 2 2 4 313 5 5 5 5 5 5 5 5 5 5 5 5 5 518 60
2 2 2 2 2 2 2 2 2 2 4 313 5 5 5 5 5 5 5 5 5 5 5 5 5 518 60
APPENDIX 2. Test results from real data
Table A.1. Fault seal calibration ^ summary.
Fig. A.2. Segmentation of the calibration area on regular 30*30 grid. The fault segments are marked with bold numbers
Fault Fault Adjacent Permeability
unit1 segment compartments [mD]
F-F' 0 1 2 .450e-02
B-B' 1 1 3 .450e-05
B-B' 2 2 3 .450e-01
B-B' 3 2 3 .450e-04
B-B' 4 2 3 .225e-04
A-A' 5 1 2 .450e-03
A-A' 6 1 2 .675e-02
C-C' 7 3 4 .225e-04
D-D' 8 3 5 .450e-06
D-D' 9 4 5 .450e-08
D-D' 10 3 5 .450e-08
D-D' 11 3 5 .450e-06
D-D' 12 3 5 .450e-08
D-D' 13 3 5 .450e-09
D-D' 14 4 6 .150e-07
E-E' 15 5 6 .150e-02
E-E' 16 5 6 .450e-02
E-E' 17 5 6 .150e-03
E-E' 18 5 6 .675e-02
A-A' 19 6 0 .150e-02
1 See also Fig. 11 in the main text