Lessons Learned in a Large Full-Field Simulation Study

The Safaniya field is the world's largest offshore field. The Safaniya reservoir is a Cretaceous-age fluvial deltaic deposit consisting of stacked sandstones and interbedded shales. The reservoir is a large shallow-relief anticline comprised of three domes: South, Central, and North (Fig. 1).

The Safaniya reservoir is typically divided into three zones: the Upper Stringer sands, Main Sand, and Lower Stringer sands (Fig. 2). Net sand thickness is 300 feet in the south, thinning to 150 feet in the north due to increasing shale content. The structural dips range from ½° to 3°. Minor faulting has been observed in the North Dome. Average porosity is 28% and core permeabilities range from 500 to 5,000 md. The reservoir contains undersaturated 28° API oil with an average solution GOR of 126 scf/STB.

The Safaniya reservoir was placed on production in 1957 when only partially delineated. Simulation modeling began in the late-1950s and primarily addressed facility requirements. The first models were coarsely-gridded, 2D areal full-field models (Fig. 3). The first 3D full-field model was composed of two layers, the upper layer modeling the oil zone above the oil/water contact and the lower layer modeling the aquifer. This energy model was used to determine long-term facility requirements.

A companion seven-layer window model was created to address completion strategy and evaluate the multi-zonal water fingering observed in the field. The dilemma was that the full-field model could not simulate the water fingering mechanism, while the window model could not be scaled up to account for field-wide reservoir heterogeneity and structural completion locations. Therefore, a combination energy and front-tracking full-field model was needed.

A 3D simulation model was built to predict reservoir performance in Safaniya field. This was done in a team environment by geologists and reservoir engineers over a ten-year period. The engineers and the geologists were supervised by their own functional supervisor. As shown by the timeline in Fig. 4, work was done in a sequential manner in that the completion of the geologic model preceded the construction of the simulation model followed by the history match.

The initial engineering work involved interpreting water encroachment from a large regional aquifer. The primary problem in the interpretation was developing a zonation to explain the water movement and reconciling it with a geologic model that was based on a stratigraphic interval concept.

Simulation work began with the building of cross-sectional models to develop interblock relative permeabilities. The initial plan called for deriving shifting factors and assuming a straight line relationship between end points. However, the cross-sectional modeling demonstrated the need for a non-straight line approach in a multi-zone model.

A detailed review of all well data (e.g., completions, work-over history, and flowmeter logs) followed. A project database and reporting programs were developed for use with Saudi Aramco's database to speed construction of well data. During this period an eight-layer geologic model was completed incorporating maps of structure, porosity, net thickness, and permeability.

A full year was taken to construct the model. The eight-layer 500-meter by 500-meter gridded model was made to fit within memory requirements of a Cray 1M computer. Geologic maps were digitized and gridded including basal shale maps for assignment of vertical permeability between the layers. Well data incorporated interpreted PI and flow splits for assignments of initial layer productivities.

Two years were taken to history match the model. Problems occurred in matching water influx in the North Dome. Additional geologic evaluations never achieved a satisfactory match probably due to more tortuous flow paths as a result of undescribed stratigraphy or faulting. It pointed to the need of either a finer layered model or a new geologic description. Problems also occurred in matching water encroachment in assumed homogeneous sections of the massive sands. Cross-sectional modeling revealed a need to increase residual oil. This was accomplished with a more heterogeneous description using geostatistics. Adjustments to the end-points in the interblock relative permeability curves were made via a sweep factor.

The final six months of the project was used to calibrate the model with a surface facility simulator and document the study.

As seen from the recounting of work, building and history matching of the simulation model took four years. Cycle time could have been reduced by using pseudomodels in a parallel planning approach as suggested by Saleri¹  since a number of problems would have been encountered and resolved sooner.

The layering scheme developed for the model resulted from an integration of geologic and engineering data. This was the first zonation effort applied to this reservoir beyond the categories of Upper Stringers, Main Sand, and Lower Stringers.

The Safaniya member was subdivided into 24 geologic intervals using chronostratigraphic markers or time lines to represent a similar time of deposition. Eight of these markers represented major transgressions of sea level that separated eight deltaic sequences. These major markers strongly corresponded to shales which were identified as significant locations of water fingering. These major markers are shown in the type log in Fig. 2.

Three types of shales have been identified in the Safaniya reservoir:

1. Correlatable shales, which are detected by well logs and can be readily correlated between two or three wells.

2. Non-correlatable shales, which are detected by well logs but are not readily correlatable between wells.

3. Fine-scale shales, which are undetectable by well logs. These thin beds are detectable in cores only, but may occupy as much as 6% to 10% of what is normally considered net sand. They act as baffles to vertical flow (i.e., limited areal extent).

The major shales appeared in wells throughout the field, yet their occasional absence suggested that there were probably holes in the shales between the wells. These major shales were strongly corroborated by pressure differentials observed in the reservoir. Downhole crossflow detected in shut-in wells also indicated that vertical drainage in the reservoir was not evenly distributed, and that these shales could create areas where pressure support from the aquifer was restricted. This engineering data strengthened the geologic correlations and provided an indication of the shale continuity between wells, although interwell shale continuity ultimately became a critical water saturation match variable.

The Safaniya model was designed to model fluid flow through the net effective pore space, not within shales or poor quality sands. The idea was to break the reservoir into correlative sand packages, or flow units. The development of the model layers, called the engineering layers, were designed to:

1. Fit within the depositional model and stratigraphic framework developed for the zonation.

2. Provide improved modeling of reservoir effects such as vertical variation of rock properties, pressure discontinuities, water encroachment (particularly water fingering), and bypassed oil sections.

3. Provide improved well modeling by describing vertical distributions of flow in the wellbore, crossflow observed in many wells, and water coning.

Eight layers were deemed necessary to reasonably capture these effects, meet the study objectives, and accommodate then-current model size constraints.

These layer tops were usually the eight major geologic markers, except in wells where there was local development of a non-correlative shale. Therefore the eight geologic layers were sometimes shifted to accommodate local shale development in a well. These adjustments are illustrated in Fig. 5.

The criteria for net sand determination was based upon historical cut-offs of 15% porosity and, in the oil zone, oil saturation greater than 50% (2.5 ohm-meters resistivity). In the field no sands below a resistivity of 10 ohm-meters were ever determined to be flowing by flowmeters so the actual contributing pore volume was probably less than represented by the net sand maps.

Maps were prepared for the eight layers which included structure tops, net sand thickness, porosity, permeability, and basal shale thickness. These maps were created by a combination of computer mapping and hand adjustments by the geologists. This preserved depositional trends and sedimentological interpretations that could not be captured by the mapping software.

Porosity maps closely followed the net sand trends although it was recognized that the contouring of average layer porosities generally washed over the fine-scale geologic variations. Porosity variations ultimately had insignificant effects on reservoir performance.

Permeability was calculated from the porosity values using a porosity-permeability transform developed from core. The permeabilities were mapped as scalar values and were likewise used as a starting point for the model permeability matrices. This approach obviously made no provision for directional variation at each point and would eventually require extensive history match modifications.

Maps of shale thickness at the base of the layers were created in a fashion similar to the net sand isopachs (Fig. 6). These maps provided the basis for initial vertical communication assignments between the layers. The level of history match adjustments reflected a significant departure from the original maps in the interwell areas.

A water encroachment study²  was conducted at the outset of the study concurrent with the geologic modeling. The objectives of the water encroachment study were to:

1. Establish the stratigraphic location of the various water intrusions throughout the field.

2. Co-develop with geology an appropriate zonation for the reservoir to capture the observed water encroachment on a stratigraphic basis.

3. Determine the effective remaining oil saturation behind the encroachment front through volumetrics and material balance.

4. Establish the basis for the model water saturation match.

Three mechanisms were identified by which water enters the producing area:

1. Vertical rise, in which the oil/water contact rises directly below a producing well from its original depth.

2. Water fingering, in which the water travels laterally over a shaly interval towards the producing wells, usually bypassing oil in underlying sands in the process.

3. Water coning, a local rise of the oil/water contact into the completion intervals.

Three sources of data were used to map where water had encroached: (1) openhole logs of wells recently drilled or deepened into swept areas, (2) cased-hole pulsed neutron logs (TDTs) and, (3) water breakthrough dates from the well histories.

The well logs spanned nearly 30 years of production history providing snapshots in time of where water had encroached. To perform meaningful mapping and volumetrics these data were projected to a common date (January 1989). This was done by estimating vertical rises of the OWCs since the most recent log date.

The breakthrough of water into a well was also used to establish when water had reached the deepest completion layer in the well. Lacking TDTs, water occurrences were assigned by reviewing completion intervals and local shale geometry from well logs. This effort required an extensive review of the well histories for over 450 wells.

The earliest producers were clustered along the crest of the reservoir. Water-free (dry) oil production was a necessity due to lack of water handling prior to 1982. When wells produced water they were plugged back above the nearest competent shale.

The combination of plugback workovers and discontinuous shale geometry created a labyrinth of water encroachment into the producing area. Water movement was dominated by these factors and vertical conformance suffered in the process. Water movement followed a stair-stepping path into the crest, moving from one stratigraphie interval to another. Large areas of bypassed oil were identified in the middle massive sand portion (Main Sand) and the less continuous portion of the reservoir sands (Lower Stringer).

One of the earliest observations from the water encroachment study was that local shale geometry and completion interval locations were the primary factors affecting water movement. Some wells contained several shales and two or more water fingers. These shales were reasonably correlative in the upper, less continuous portion of the reservoir (Upper Stringer) but less so in the Main Sand; therefore, marker depths could not be based solely on water finger locations, i.e., they were not always stratigraphie equivalents.

The presence of shale also retarded vertical water rise in the Main Sand and Lower Stringer. This permitted the dry producing strategy for many years. Some thick shales appeared competent in many wells which did not retard water rise, presumably due to limited extent. In nearly every case, however, the shales would trap some amount of oil beneath them, reducing vertical conformance. Although the thickness of the trapped oil appeared related to the shale extent there was poor correlation between shale thickness and shale extent. Shale extent was often inferred between wells which exhibited similar water encroachment.

Swept thickness isopachs were prepared for the eight major reservoir layers. These maps illustrated the thickness of the swept sand in the wells and portrayed the encroachment pattern interpreted in the interwell areas. The swept sand maps also provided the basis for the model water saturation match.

By comparing the bulk volume of swept reservoir from these maps to the cumulative oil production in the South, Central, and North Safaniya areas, the effective remaining oil saturation was calculated behind the front. An average value of 35% was calculated for each area, compared to the corefiood measurements of 14% to 22%.³  The increase in ROS was attributed to trapped oil and reduced conformance in areas which were mapped as swept. These results assumed greater significance during the history match.

There were about 400 wells that had been completed in the reservoir at the time of the study. All of these wells were worked over several times during their producing lives, requiring frequent model updates.

The completion histories for these wells were recorded in a database in the form of perforation depths (in MD and VSS) versus time. A computer routine was used to create the model input cards for completion layers versus time. Wells were often hand-corrected to account for data entry errors, leaky bridge plugs, and other extraneous data which affected the distribution of flow and were not accounted for in the database, underscoring the need for detailed engineering control in reservoir simulation.

The productivity indices for each well were assigned using PI tests from the field when available. If PI tests were not available for a particular configuration, the PIs were calculated using the completion cell properties. Distribution of the PIs to the individual layers were assigned by flowmeter data when available. If flowmeter data was not available for a particular configuration, then the allocation of PI was based upon permeability-thickness weighting of the processed well log data. Due to the large number of PI and flowmeter tests, and the frequent well configuration changes, a FORTRAN routine was used to compile these data and create the initial PI assignments to the completion layers (Fig. 7).

To capture the common occurrence of water coning, an analytical coning correlation⁴  was implemented in the model. This correlation calculated WOR as a function of rate and several geometric factors, such as distance to oil/water contact, vertical permeability, and fluid properties. This required additional CPU time for each time step but provided the rate sensitivity that common well psuedos could not provide. Coning parameters were eventually adjusted to match breakthroughs and water cut versus time performance in the wells as needed.

The eight model layers corresponded to the eight engineering layers developed in the water encroachment study. The regional aquifer was described by a single-layer description, while the Safaniya oil reservoir was described with eight layers. The adjacent fields were described by a three-layer description.

At the edges of both the three-layer and the eight-layer description (Fig. 8) are rows of cells used to provide a connection between the two different layering schemes. The intent of these chimney cells is to abruptly increase vertical resolution without introducing discontinuities to reservoir properties or behavior (i.e., allow pressure equilibrium).

The areal grid for the model (Fig. 9) is 93 x 186. The Safaniya reservoir was defined by a fine grid of 500 m by 500 m square cells. The cells outside the Safaniya reservoir were coarsened outward to cover the entire aquifer and fit within model size constraints. The total number of cells is 138,384 of which 61,257 are active.

The model grid includes all fields with oil and water withdrawals from the Safaniya member.

The initial model permeabilities were derived from core-based transforms scaled-up with permeabilities from pressure buildup data. The scale-up factor derived from this analysis was 1.156 with a high degree of variance due to the combination of sand types in the reservoir.

Vertical permeabilities were assigned to model cells based on a K_v/K_h ratio from whole core and plugs and the fraction of area covered by basal shale. Since the largest lithology modeled as pay was clean sandstone, a K_v/K_h ratio of 0.5 representing this lithology was chosen as the initial assignment. This ratio was used for sand on sand connections; however if a fractional area of the cell face was covered by shale as indicated from basal shale maps, then that same fraction was used to reduce the vertical permeability.

Resistivity logs and capillary pressure data indicated that the oil-water transition zone is negligible.

The rock relative permeability curves needed to generate interblock relative permeability were created from a set of laboratory-measured curves with a 9% residual oil saturation. During the process of generating the interblock curves, results of recent native-state core indicated a residual oil saturation of 14%, resulting in the work having to be redone.

Originally interblock oil/water relative permeability was to be generated by assuming laboratory measured residual oil and a straight-line relationship for both. The holdup factor for the water saturation was determined by running finely gridded cross-sectional models. Using the Kyte and Berry method,⁵  the shape of the curves proved important in matching the performance between coarsely-gridded and finely-gridded cross-sectional models. The difference between the straight line approximation and a curve generated by the Kyte and Berry method is shown in Fig. 10. Six variables were shown to have significant effect on the shape of the curve:⁶  net thickness, reservoir dip angle, oil viscosity, irreducible water saturation, basal shale, and the location of the cell from the original water/oil contact.

The major assumption in generating the interblock curves was that of homogeneity, in that properties were kept constant in each cross-sectional model. This was in keeping with the geologic interpretation; however, this proved to be a critical flaw during the history match. The 3 man-years effort and sequential approach in generating the relative permeability curves ultimately proved to be time-consuming, costly, and was probably over-engineering of a probabilistic parameter.

The history match was separated into a pressure match and a saturation match. The pressure match used shut-in bottomhole pressures and RFT measurements. Flow splits were matched for each well using flowmeter data. The overall pressure match was achieved first by adjusting regional aquifer properties and then locally adjusting properties within the Safaniya reservoir.

The saturation match was achieved in three steps. First the water encroachment pattern from the water encroachment study was matched by layer, and then water arrival times from TDT and wet production observations were matched on an individual well basis. Finally water cuts were matched for each well.

Approximately 5,630 shut-in bottomhole pressures were used in the pressure match. Bottomhole pressure data was not adjusted with a Peaceman correction⁷  since the 24-hour shut-ins are beyond the range of noninterference effects.

The initial estimate of the aquifer size was about 2,000 times the oil volume reflective of an infinite aquifer. Using RFT data from surrounding aquifer wells, the aquifer pore volume was adjusted. To achieve the volumetric average pressure match in the oil reservoir, permeabilities were increased by a factor of three from the core-based transform.

A more precise pressure match was obtained by altering the lateral and vertical permeabilities through the use of multiplying factors that tied to geologic interpretation. One type was a sand quality factor related to the initial water saturation of an area. Layers tended to group themselves by area into distinct and definable initial saturation ranges. These saturation ranges were thought to be closely related to sand quality. The sand quality factors ranged from 3.0 for the highest quality sand to 1.0 for the lowest quality sand.

An additional concept, that of a connectivity factor, was developed as a measure of sand continuity, since the ability for fluid to flow from one cell to another still depends on the length or tortuosity. Reservoir sand in south Safaniya field was considered to have the best connectivity (a factor of 1.0) and the sand in the North Dome the poorest (a factor of 0.5). The total multiplying factor for a given region was then the product of the two factors.

To match individual well pressures both horizontal and vertical permeability were changed on a localized basis. The resulting pressure match statistics are shown in Fig: 11. Since permeability was the most altered variable the resulting permeability map did not have the same character as the original mapped permeability from the geologic model. Given the extent of these adjustments, the permeability maps were another example of costly over-engineering of a probabilistic parameter.

The flow split match was achieved to within a difference of 5% with very little changes. Generally vertical permeability was changed locally between layers. The ease of matching the flow split is attributed to the time spent in assigning PI during the model building step and the lack of much differential depletion in the reservoir.

The water encroachment patterns were matched at a specific point in time (January 1, 1989) whereas both the arrival times and water cuts were matched throughout history.

Three basic variables had the greatest impact on matching the encroachment patterns and water arrival times. The first was the cell vertical transmissibility used to direct the encroaching water into the proper layers (i.e., matching vertical or lateral movement). The second was the volumetric sweep efficiency used to control the velocities of the encroaching water once it was in the proper layers. The third was the cell vertical permeability which was used to fine-tune water cuts.

The concept of conformance was introduced for the purpose of history matching water encroachment patterns. During the pressure match, the model was not able to duplicate the extent of advancement in the encroaching edge water observed in the field even with variation in pore volume (e.g., porosity). Fig. 12 compares the water encroachment pattern for layer 6 at January 1, 1986, between the observed values (from the water encroachment study) and the one from HM064. The model underpredicted the profile. Any amount of buckling introduced in the curve (to provide higher mobility to water) did not achieve the desired results. Lowering the holdup likewise did not show any improvements. The only factor which moved the water farther into the reservoir was an increase in residual oil.

When the interblock curves were generated the layers were assumed to be homogeneous. Shale was assumed to be only basal shale. A study⁸  using cross-sectional models containing dispersed shale indicated that the presence of fine-scale shales gave rise to higher residual oil saturations, with oil trapped mainly below and behind shale streaks.

A solution to the problem would be to run new cross-sectional models incorporating different distributions of fine-scale heterogeneities to generate new interblock relative permeabilities. However this would have delayed the history match and the resultant curves would have been similarly adjusted to match specific behavior. The conformance or un-swept-factor approach was a more expedient approach. It not only moved the encroaching edge water significantly into the reservoir but also accelerated the progress of the history match.

The unswept-factor is a measure of how much of the reservoir will be volumetrically isolated. It was incorporated into the simulator as another attribute of a grid cell to reflect reservoir heterogeneities. The simulator internally assumes a portion that is unswept within the grid cell remains at an irreducible saturation while the remaining portion can be swept.

The unswept-factor affects the residual oil saturation which in turn depends on fine-scale heterogeneities. The unswept-factor can be used to represent the degree of heterogeneity not initially incorporated in the geologic model. The residual oil saturation is assumed to remain dispersed in the reservoir, being trapped below and behind the fine scale heterogeneities.

Reductions in vertical permeability also had a significant impact on matching the water encroachment patterns. The cross-sectional modeling suggested that a reduction in vertical permeability was necessary to account for the effect of fine-scale heterogeneities. Accordingly, the vertical permeabilities within the sands were reduced by a factor of 25 giving a final assignment of K_v/K_h of 0.02.

The water-cut match was achieved by adjusting vertical permeabilities in the well cell and coning parameters used in the coning correlation. Fig. 13 shows the overall quality of the match for the Safaniya field.

1. A thorough review of water encroachment is vital to understanding reservoir behavior and modeling historical water movement.

2. The scalar mapping of permeability was not beneficial to the history match because it did not capture the tortuosity of flow and the reservoir anisotropy.

3. The effect of shale geometry between and within layers on water movement is more significant than porosity variability.

4. Changes to endpoint saturations of relative permeability and pore volume accounted for the effects of reduced net pay and fine-scale shales in the history match of water movement.

5. The distribution of rock properties beyond control and in interwell areas are more efficiently captured in the history match process than in the geological modeling. This is evidenced by the need for connectivity and sand quality factors and the level of adjustments to the basal shale and permeability maps.

6. The independent a priori generation of detailed pseudo-relative permeability functions did not facilitate the history match process, hence they were time and cost ineffective.

7. The 500-meter x 500-meter cells were not suitable for area! front-tracking and near wellbore coning despite their relative size to the reservoir.

8. Sequential project management allowed setbacks which extended each phase of the study.

9. The results of history matching should be applied to the geologic and water encroachment interpretations in a cyclic manner to achieve a fully-integrated reservoir description.

We thank all the members of the Safaniya team for their contributions to the study, in particular: Messrs. H. H. Al-Zahrani and A. I. Al-Towailib from Saudi Aramco; and Messrs. U. K. Acharya, W. R. Coffelt, D. H. Jones, J. L. Stockey, J. H. Oldfield, R. B. Haulenbeek, R. H. Kirby, T. G. Brickner, R. P. Lockwood, M. F. Mendeck, and Ms. D. L. Adamczak from Chevron.