Introduction

The X field located in northern of Iraq, a primary evaluation on field potential had been performed at exploration phase using 3D seismic data and information from two wells, accordingly a light field development plan proposed for a potential gas development of the X field reservoir. The target gas production rate provided by operator Company is 175 Mscf/day to cover the demand of power plant, the plan proposed to drill another 18 wells with potential similar to existing two drilled wells, the result of the gas production forecast shows that 20 wells are required to achieve an 7 years plateau period of 175 MMscf/day. At first phase, the target rate can be achieved with 7 wells in total (2 existing and 5 new wells), but reservoir re-evaluation at early production phase by using dynamic data shows lower reserve potential because of Asphaltene (Bitumen) occupying part of pore space which effectively lead to revise field development plan; especially on reducing target rate, number of producer wells, and managing Bitumen production problem..

Correctly and continuously evaluation of Hydrocarbon initial in place by using different techniques is worthy effect field development plan and strongly influences project economic feasibility. Generally, in carbonate gas reservoirs or any other reservoir type, volumetric calculation is performed at exploration phase, depending on structure geometry and both rock & fluid property distribution without taking any consideration for reservoir dynamic energy, which put high possibility to have an overestimate Hydrocarbon initial in place. While at development field phase, both simulation history match analysis and material balance calculation techniques can be managed to re-assessment reservoir potential and narrow down actual fit value with valid evidences completely explain any disparity reflect reservoir dynamic flow energy. Material balance is one of the fundamental tools of reservoir engineering, the plot of z/p versus cumulative gas production, Gp, is a widely accepted method for solving the gas material balance (Dake 1978; Pletcher 2000). In this study case, carbonate gas reservoir potential evaluated by three different techniques; volumetric calculation, simulation initialization with & without history match (HM), and material balance calculation. Volumetric calculation estimates GIIP value almost two times bigger than other techniques estimation. Simulation history match result confirms material balance estimation for reservoir potential of four years production and pressure history data, the overestimated value from volumetric calculation explain by Asphatena “Bitumen” existing. In gas reservoir, continuously updating of simulation analysis and material balance calculation is required using every last production and pressure data to optimize field managing. Finally, both volumetric calculation and simulation analysis should not be viewed as a replacement to material balance because the latter can yield valuable insights that can be unobserved from other two techniques. Performing a separate material balance study will usually improve overall reservoir understanding and enhance any subsequent reservoir evaluation and simulation study as well as material balance should be viewed as a complement to simulation, not as a competing approach, and using both to improve analysis of hydrocarbon reservoirs (Dake 1978).

3D geological static model

The field has a clear surface expression, manifested by a sinuous anticline. The general stratigraphy framework of the area is simply characteristic by Folded Zagros, the sedimentary sequence has more than 10 km thick; start with a ductile Upper Precambrian series, generally shows an alternance of shallow-water and often dolomitized carbonates with more ductile marly intervals and evaporitic layers, that can act as decollement surfaces (Jassim and Goff 2006). The study zone is sited in the north-eastern margin of the Arabian Plate. The Nabitah Orogeny affected this area from 640 to 680 Ma (Stern and Johnson 2010), in the Late Precambrian, the NW–SE trending Najd rifting event (520–610 Ma) started. This formed the NW–SE trending Najd fault system in Iraq, which is prominent in the basement and influences the boundaries of tectonic zones in northern Iraq (Jassim and Goff 2006). N–S trending grabens and tilted fault blocks were formed in the Arabian Platform during the Cambrian–Carboniferous.

Petrel geo-modelling software Schlumberger, 2015 version is used to build a static reservoir model over target reservoir interval of study field. The static model is result of a major incorporates of interpretation outcomes of 3D seismic survey acquired tied in with acquired well log data and actual well formation tops from six wells that penetrated target reservoir.

Framework of modeling 3D Structural Model using Petrel software consists of fault modeling, gridding, top and bottom of reservoir surfaces/horizons/zones, and rock petro-physical properties modeling. The uncertainly evaluation of geological model is concerned on representing the complexity of existing fault.

A fault network plane was used as the base of the fracture architecture. Major faults were defined as active to set up horizons structure and displacement. The minor faults were defined as non-active and were used to control the dynamic simulation. Totally fifteen faults were picked, most of which were parallel to the fold axis of orientation, NW–SE, as shown in Fig. 1. The figure clearly indicates that the SW plunge of the structure was deformed by a main garben sealed fault.

Fig. 1
figure 1

Top structural depth 3D map and fault network

The geo-cellular pillar gridding was generated using a fault model with a zig-zag trend attribution. The grid cell was defined by following the boundary of the area of interest and clipping the non-interest area, particularly in the southern flank.

Cell size was dimensioned as medium to coarse. A resolution grid of 150 × 150 m was fixed to represent the horizontal heterogeneities, preserve the complexity of the geology, and facilitate practical simulation times.

Reservoirs are generally grouped according to their common features. The zonation in this research primarily used lithological and petro-physical criteria from data of six wells. In particular, the most important differentiating criterion was the vertical distribution of the porosity (Adnan et al. 2010; Nooruddin and Enamul Hossain 2011). Three zones were identified vertically along the structure: The upper zone was characterized by a rich porous unit, with porosity higher than 20% and dominated by vuggy dolomite; the middle zone had a porosity range between 15–20% and was dominated by porous dolomite; and the lower zone had a poor porous unit of tight limestone, as shown in Table 1 and Fig. 2.

Table 1 Reservoir zonation and rock petro-physical properties
Fig. 2
figure 2

Reservoir zonation base on well correlation data

The Geo-cellular pillar gridding model had following dimension;

  • I × J × K = 27 × 121 × 3,

  • Dx × Dy = 150 × 150 m,

  • Dz = 33–105 m.

  • Total grid blocks = 9801, see Fig. 3.

Fig. 3
figure 3

3D Geo-cellular gridding model

Fractures were introduced by discrete fracture network model (DFN) (Dershowitz et al. 2004), the mean value of the fracture porosity was 0.06%; low fracture storage can be considered. Additionally, the mean value of the sigma factor was 6 m−2; however, fracture properties populated in 3D gridding model to support reservoir intake permeability. This means that every grid-block is associated with 1 fracture element, which effectively doubles the total number of grid-blocks to 19,602.

Geological daily drilling report, wells log data and core data are used to evaluate target formation and establish matrix property model. Wire-line jobs performed for all the six wells that penetrated target formation by Schlumberger Company, integrated raw logs data had been collected, processed, and calibrated with available core analysis data from four wells.

Rock property modelling, porosity, saturation, and permeability models were scaled up into the grid cells adjacent to each well and then distributed throughout the whole model using an appropriate algorithm; Gaussian random function simulation method. See Figs. 4, 5, and 6.

Fig. 4
figure 4

Matrix porosity model

Fig. 5
figure 5

Water saturation model

Fig. 6
figure 6

Matrix permeability model

Volumetric calculation run for 3D geological model using following input parameters: matrix and fracture property models, set up gas down to (GDT) at depth 1495 m below sea level, condensate (vaporized oil) Gas ratio (CGR) = 30 stb/MMscf, and 0.005 ft3/stander ft3 as gas formation factor value. The calculation result introduces the value of original gas in place (GIIP) equal to 682 Bscf, and original oil (condensate) in place (STOIP) equal to 21 MMSTB, as shown in Table 2.

Table 2 Volumetric calculation based on static model

Dynamic flow simulation model initialization

The simulation model for this case study was a compositional model with a 3D Cartesian grid system. The dynamic flow grid had the following dimensions which is the same dimensions of Geo-cellular gridding model:

  • Nx × Ny × Nz = 27 × 121 × 3

  • Dx × Dy = 150 × 150 m

  • Dz = 33–105 m

  • Total grid blocks = 9801

Fractures were introduced by preparing the model in dual permeability mode. This means that every grid block was associated with one fracture element, which effectively doubled the total number of grid blocks to 19,602.

The faults in the dynamic model were copied from the static model. Major groups of fault were modelled from partially to fully transmissible, while minor faults were modelled to be fully transmissible owing to their relatively low throws.

Simulations were conducted for the dual-porosity/dual-permeability model. The matrix permeability was not enough to explain the sufficiently high productivity observed in the wells throughout the field (Nelson 1992). The current understanding of the reservoir is that of a typically fractured carbonate in which the matrix acts as storage and the fractures as conduits to the flow (Ben 2017).

The equation of state (EOS) compositional fluid model created using PVTi software was imported to the dynamic model and tied in with initial reservoir and separator conditions.

Three matrix saturation function curves (relative permeability and capillary pressure for the dolomite, dolomitic limestone, and limestone units) were imported to the dynamic model. Additionally, one saturation function curve was used for the fracture system.

Fracture properties were determined as part of the fracture network modelling exercise. Prior to the history match, the fracture properties are defined as follows:

The mean value of the fracture porosity was 0.06%, as no fracture storage was considered; hence, dual permeability simulation was conducted. Additionally, the mean value of the sigma factor was 6.5 m−2.

  • The mean value of fracture permeability was approximately 5420, 3300, and 4894 mD for the I, J, and K directions, respectively.

  • Eclipse (E300) simulator used to initialize the dynamic model to adopt gas condensate compositional fluid model.

Saturation curves from rock model used to define the water distribution in the dynamic model have been adopted from the core analysis. Three saturation zones have been defined in the dynamic model following the geological zonation, i.e. Upper dolomite, middle dolomitic limestone, and limestone. The relative permeability and Pc curves for the different unit are presented for the reference case of residual water saturation in the matrix is 0.16, 0.35 and 0.7 per each zone. The relative permeability curves defined for the fractures assume segregated flow and have been set up accordingly with a zero capillary pressure.

The reservoir pressure is initialized by setting it to 2688 psia at the GDT corresponding to the dew point pressure. A ternary plot of the reservoir at production start is provided in Fig. 7.

Fig. 7
figure 7

Ternary plot

The dynamic model initialization presents GIIP value equal to 714 Bscf and original oil (condensate) in place (STOIP) equal to 21 MMSTB, and this value is more than GIIP value from static model by about 26 Bscf which is considered an accepted margin range. This disparity is based on several factors, and some of the factors are enumerated in Table 3 (Schlumberger Manual 2014).

Table 3 Volume calculation versus simulation case initialization (Schlumberger Manual 2014)

Material balance (MBAL) model description and calculation

Using gas material balance traditional p/Z (Dake 1996; Ahmed and McKinney 2005) versus cumulative production plot in attempting to history match the performance of gas reservoirs, which is supposed to be one of the simplest subjects in the whole of reservoir engineering and, indeed, the physics and mathematics are quite trivial, if unappreciated, it can lead to a complete misinterpretation and serious error in assessing the reservoir drive mechanism (water-drive or volumetric depletion) of drive mechanism and a serious overestimation of the GIIP. In fact, it is the intention to mention this technique in conjunction with other methods. So the main purpose of the MBAL model is to determine reservoir type and compare estimated GIIP with GIIP from geological and reservoir models.

Petroleum exports IBM tools used to build up material balance model: consist of reservoir tank and six wells, as shown in Fig. 8.

Fig. 8
figure 8

Material balance model: reservoir tank and six wells

The PVT fluid model description used in the model originates from EOS composition model generated by PVTI tool/ of the tuned PVT sample from well X-02, giving a gas condensate fluid type. The required input parameters were define for initial reservoir condition; temperature, pressure, rock and fluid compressibility, and average petro-physical properties were used as defined previously in Table 1, see Fig. 9.

Fig. 9
figure 9

Material balance model: reservoir parameters

The model used individual well by well production and pressure history data of four years, the base case in place volume has been incorporated in the model and no aquifer support has been specified, and Fig. 10 shows wells production and pressure history data.

Fig. 10
figure 10

Material balance model: reservoir parameters

The traditional gas material balance concept of p/Z versus cumulative production plot is widely utilized as a powerful tool for evaluating the performance of gas reservoirs, as this method produces linear behavior;

$$ \begin{gathered} {\text{P}}/{\text{Z }} = {\text{ Pi}}/{\text{Zi }}{-} \, \left[ {\left( {{\text{Pi}}/{\text{Zi}}} \right){ 1}/{\text{G}}} \right]{\text{ Gp}} \hfill \\ {\text{OR}} \hfill \\ {\text{P}}/{\text{Z }} = {\text{ Pi}}/{\text{Zi }}\left[ {{1 }{-} \, \left( {{\text{Gp}}/{\text{G}}} \right)} \right] \hfill \\ \end{gathered} $$

where Pi; initial reservoir pressure psi, Zi; initial gas deviation factor, P; average reservoir pressure psi, Z; gas deviation factor at P, G; gas initial in place mmscf, and Gp; cumulative gas production mmscf.

The using of Cole (1969) plot provides a history matching approach to determine if the reservoir is volumetric depletion reservoir or it is supported with water drive, also to get a good estimate of gas initial in place. The plot is derived from the generalized MBE as given in an expanded form by;

$$ \begin{gathered} \left( {{\text{Gp Bg }} + {\text{ Wp Bw}}} \right) \, / \, \left( {{\text{Bg}} - {\text{ Bgi}}} \right) \, = {\text{ G }} + \, \left[ {{\text{We}}/\left( {{\text{Bg }}{-}{\text{ Bgi}}} \right)} \right] \hfill \\ {\text{OR}} \hfill \\ {\text{F }}/{\text{ Eg }} = {\text{ G }} + \, \left[ {{\text{We }}/{\text{ Eg}}} \right] \hfill \\ \end{gathered} $$

where: Bgi; gas formation volume factor at initial reservoir pressure bbl/scf, Bg; gas formation volume factor at certain reservoir pressure bbl/scf, Wp; cumulative water production STB, Bw; water formation volume factor bbl/STB, We; cumulative water influx bbl.

Both traditional P/Z versus cumulative production and Cole plots of four years production and pressure history data show history match of gas reservoir with no evidence of aquifer existing or limit aquifer supporting, see Figs. 11, 12 and 13, the methods estimate the GIIP equal to 478 Bscf. Material balance GIIP (478 Bscf) is significantly lower than of volumetric calculation (682 Bscf) and simulation initialization (714 Bscf); this significant uncertainty on this disparity needs more investigation and explanation. Asphaltene issue can be managed by introduce Asphaltene amount as cut off of pore volume, to incorporate model real potential.

Fig. 11
figure 11

Pressure history match plot

Fig. 12
figure 12

Traditional P/Z plot

Fig. 13
figure 13

Cole methods

History match and simulation result

The purpose of a history matching is to build a simulation model consistent with actual data to capture the level of detail necessary for production forecasts with a high predictive power.

In this study, the fluid flow in the naturally fractured reservoir was described and presented in dual‐porosity/dual‐permeability model (Kazemi et al. 1976; Warren and Root 1963). The simulation was conducted in gas rate control mode, which means that the gas rate was exactly matched. The other parameters to be matched were condensate rate, reservoir pressure, productivity index, basic sediment and water, and flowing tubing head pressure.

The initial simulation run (Base case) result does not show matching for reservoir pressure, but shows a higher pressure trend comparing with actual measured reservoir pressure trend, see Fig. 14, that could be explained as a reflection of smaller reservoir potential comparing with the potential that calculated from volumetric calculation and simulation initialization.

Fig. 14
figure 14

Initial simulation result

Re-run performs for the initial simulation (Base case) after applying the usual reservoir engineering approach, modifying the porosity value (pore volume) to match the depletion levels for each well. A match is obtained by multiplying the overall porosity by 0.725, after the existence of bitumen within matrix pores and fractures was confirmed. The model initialisation indicated an overall gas initially in place value of 476 Bscf. Note that this value corresponded with the preliminary material balance calculations.

As shown in Fig. 15, a good history match achieved for dual‐porosity/dual‐permeability model for all parameters; both matrix blocks and fracture network storage and conductivity are significantly contributed to equilibrium deliver reservoir potential.

Fig. 15
figure 15

Initial simulation result after applying potential reduction

Reservoir potential affected by Asphaltene “Bitumen”

From the core samples, Asphaltine “Bitumen” existing confirmed with significant amount occupied a series of matrix pores and filling some fractures, as well as Asphaltine accumulation observed at production facility while production from wells X-05 and X-06. Figures 16, 17, and 18 show asphaltine existing in recovered core samples, and received samples at production facility.

Fig. 16
figure 16

Core samples shows Bitumen occupied matrix pores

Fig. 17
figure 17

Lab core analysis and description shows fracture filled by Bitumen

Fig. 18
figure 18

Collecting Bitumen samples at production facility

So bitumen occupation for matrix pores and fractures within reservoir interval is used to explain the disparity of reservoir potential that getting from volumetric calculation, simulation initialization without history match (HM), simulation initialization with history match (HM), and MBE calculation as shown in Table 4.

Table 4 Reservoir potential

Asphaltene “Bitumen” issue is managed to be introduced as cut off of pore volume to incorporate reservoir actual potential and simulation history match, and this reduction in pore volume (reservoir potential) is considered as valid for reservoir real recording intake. This reduction also can be a combination of higher water saturation and/or lower NTG.

Conclusions

  1. 1.

    Both material balance and reservoir simulation with history match estimated same reservoir potential (GIIP = 478 Bscf), while volumetric calculation and reservoir initialization assessment give higher figure of GIIP (682–714 Bscf), almost twofold of the actual reservoir potential.

  2. 2.

    Reservoir simulation and material balance are complementary rather than competing tools; reservoir simulation does not eliminate the need for classical material balance analysis.

  3. 3.

    Core analysis confirms Asphaltene “Bitumen” occupation part of matrix pores and filling some fractures within Upper study formation interval, and it is used to explain and confirm that the lower actual figure of reservoir potential for gas reserve (GIIP = 478 Bscf).

  4. 4.

    Material balance in line with simulation history match can provide valuable understandings into reservoir potential, mechanisms and processes that unnoticed by volumetric calculation or simulation analysis.

  5. 5.

    Material balance is used to assist for setting up correct reservoir potential for run reservoir simulation initialization and for getting simulation history match.