Lithology and fluid discrimination using rock physics-based modified upper Hashin–Shtrikman bound: an example from onshore Niger Delta Basin

The Ozifa reservoir is proven reservoir that cuts across the Northern and Greater Ughelli depo-belts of the Niger Delta Basin. This reservoir possesses heterogenous character southward of the field, making elastic properties, lithologies and fluid types difficult to describe accurately. In this study, rock physics template was applied to porosity and acoustic impedance (AI) crossplot clusters to illustrate rock–fluid relationships using modified Hashin–Shtrikman upper bound, Voigt upper bound and Reuss lower bound, as an input in the template. Values of acoustic impedance and porosity were used as lithofacies classification parameters for discrimination of lithofacies and fluid types. Our result showed that modified Hashin–Shtrikman upper bound line when applied in acoustic impedance (AI) and porosities (φ) crossplot domain discriminated gas-filled reservoirs from brine filled reservoirs and shale effectively. Similarly, results from crossplot showed clear separation of shale, heteroliths filled with brine and gas bearing sand, which was not plausible using conventional petrophysical analysis. This approach was successfully applied in analysing lithofacies and fluid relationship in different well locations and serves as a model for successful prediction of different lithology and fluid types, a major requirement for determining effects of geological variables such as sorting, clay distributions on the reservoir connectivity and optimum production using time-lapse (4D) seismic interpretation.


Introduction
The understanding of the changes in reservoir properties has helped in understanding the nature of the geological variables such as sorting, grain contact and facies distributions controlling reservoir quality (Bjørlykke 2010). In recent times, rock physics crossplots and modelling techniques involving acoustic impedance (AI) and porosity (φ) have been used to predict fluid type distributions, changes in reservoir properties and perform feasibility study for time-lapse (4D) seismic interpretation (Toqeer and Ali 2017;Avseth et al. 2005;Mavko et al. 2009). It is also known that rock physics-based models, when combined with acoustic impedance inversion, will improve the accuracy of reservoir quality predictions from seismic data (Avseth et al. 1999, Veeken andDa Silva 2004;Sen 2006). In this study, we cross plotted acoustic impedance with porosity and constrained sediment distributions with modified Hashin-Shtrikman upper bound, Voigt upper bound and Reuss lower bound to establish a trend that is consistent with local geology (Han 1986;Avseth et al. 2005;Mavko et al. 2009;Bjørlykke 2010). The objective of this study is to use modified Hashin-Shtrikman upper bound to discriminate lithofacies and fluid types distribution in the southern part of the field using well log data. The methodology used to accomplish the objective is presented in the workflow shown in Fig. 1.
According to Reijers (2011), the study area ( Fig. 2) has been characterized to be shalier in the southern part than the northern part of the field, and the main productive sand interval known as the Ozifa sand cuts across the Northern and Greater Ughelli depo-belts.
There are reports that three wells are producing from the sand while three were marginal, one abandoned and two dry wells in the southern part of the field (Shell 2007). The stratigraphy of the study area (Fig. 3) is detailed in Short and Stauble (1967); Doust and Omatsola (1989); Lawrence et al. (2002) and Reijers (2011). The extensional features favourable for hydrocarbon accumulation according to Evamy et al. (1978) and Stacher (1995) are simple rollover faults, multiple growth faults, antithetic faults and collapsed crest faults. The local stratigraphy of the Ozifa sand began with the deposition of proximal deltaic deposits and channel units that are separated by laterally extensive shale packages that represent flooding episodes, as shown in Fig. 4.   Fig. 1 Simple flowchart of petrophysical and rock physics models used in this study. The end product is lithology/fluid discrimination and final localized rock physics template for the study area

Study methodology
Suites of well log from three selected wells were quality checked (QC) before loading into geological and geophysical software such as Interactive Petrophysics (IP), Rokdoc 6.1.4 and Petrel 2014 for various interpretation (see Table 1). The study adopted petrophysical analysis and well log rock physics template (RTP) analysis of three producing wells (D, A and C). The limitation of this study mainly lies in the fact that there is no core petrophysical data to establish relationships between well log-derived petrophysical parameters and core-derived petrophysical parameters. Secondly, there is no laboratory data on cored sandstone interval to infer the nature of the sandstone (unconsolidated or cemented).

Petrophysical analysis
The petrophysical analysis was primarily achieved by using Interactive Petrophysics (IP) tool to quantitatively determine the reservoir properties across wells D and A. The Well C reservoir properties were determined manually following equations stated below. The reservoir properties studied include total porosity ( ), water saturation (Sw), average porosity av , average effective porosity Ave eff , average water saturation S av and the average volume of shale ( Vsh av ) . The total porosity, was calculated using the density-porosity relationship as shown in Eq. (1).
where φ is the density derived porosity, ma is the grain or matrix density, b is the bulk density and fl is the density of fluids residing in the pore spaces. The effective porosity was calculated by removing the effects of clay on the grain matrix and which was determined using Eqs. 2 (Asquith and Krygowski 2004).
(1) where sh is the bulk density of adjacent shale, V sh is the volume of shale, IGR is the gamma ray index, GR log is gamma ray log value, GR max is the maximum value of gamma ray in the data and GR min is the minimum value of gamma ray in the data. Also, water saturation, S w was obtained from Eq. (3) as expressed by Asquith and Krygowski (2004): where R w is the resistivity of the formation water, R t is the true resistivity of the formation as indicated on the ILD log, "a" is the tortuosity factor with an average value of 1, m is cementation factor that determines the pore geometry and connectivity of pore spaces, and n is saturation exponent determines the spatial distribution of water in the pore space; m and n are assumed to equal 2. The averages of porosity av , water saturation Sw av , and the volume of shale ( Vsh av ) was determined using Eq. (4-6) as expressed in Darling (2005): where n is the number of sample interval, i is the i th input value of porosity, h i is i th input of thickness interval, S w is water saturation while Vsh i = ith input of volume of shale.

Well log rock physics template analysis
Rock physics template (RPT) analysis provides an important framework for predicting the effects of rock fabrics and pore geometry in a mixture of rocks and fluids using well log data. The rock fabric and pore geometries are very important in understanding the stiffness (and velocity) of rocks (Mavko et al. 2009). Rock stiffness has a key influence on the rock frame, which is critical in determining the velocity and rate at which fluid affect compressional velocity (Vp)  and Reuss lower modulus describe the upper and lower bounds for quartz and water mixture, and it is expressed as: where K qtz = modulus of quartz, K w = modulus of water, Vol qtz = volume fraction of quartz and Vol w = volume fraction of water.

Hashin-Shtrikman and modified Hashin-Shtrikman bounds
The Hashin-Shtrikman bound was presented by Hashin and Shtrikman (1963) where K, µ and f represent the bulk moduli, shear moduli and volume fractions of the individual phases 1 and 2. The upper bound is calculated when the stiffest material is donated by 1 and softest material is denoted by 2. The lower bound requires that the stiffest material is donated by 2 and the softest material is donated by 1. The modified Hashin-Shtrikman upper bound is a modified assumption of Hashin-Shtrikman bound that separates consolidated and unconsolidated elastic moduli trends according to rock critical porosity, φ c . This bound is constructed replacing the volume fractions of the individual phases (1 and 2) in Eqs. (8) with c (Nur et al. 1991(Nur et al. , 1995. According to Avseth et al. (2005) and Mavko et al. (2009), this bound describes trends of sediments as it becomes progressively compacted and cemented away from the lower (Reuss) bound due to the rock critical porosity and sorting (Fig. 6a). The contact cement model assumes that porosity reduces owing to the uniform deposition of cement on the surface of the sand grains and only a small amount of cement deposited at grain contacts is required to significantly increase the stiffness of the rock (Simm and Bacon 2014). This model works in line with the critical porosity model where it describes the trend of a rapid increase in elastic stiffness of sand without a great change in porosity as cement is introduced in the sediments as shown in Fig. 6b (Dvorkin et al. 1994, Avseth et al. 2005. The critical porosity, φ c , values in the Niger Delta Basin range between 0.36 and 0.40.

Acoustic impedance, AI and porosity crossplot analysis
Acoustic impedance (AI) is known as the product of density ( ) and compressional velocity (Vp). This impedance is often used to discriminate lithologies, fluid types and Schematic diagram of the rock physics models for sands at critical porosity (Marion 1990) to predict porosities at well locations using log and seismic volume. Omudu et al. (2007) and Abe et al. (2018) used acoustic impedance, AI crossplot and seismic inversion to discriminate brine filled reservoirs from hydrocarbon-filled reservoirs in an oil field of the Niger Delta Basin. Despite its suitability for discriminating lithologies and fluid types, gamma ray log is often used to guide lithology classifications. In this study, acoustic impedance crossplot with porosity was utilized to discriminate lithologies and fluid types in space. This was done using the crossplot/polygon tool in RokDoc 6.1.4 to carefully zone and separate lithologic/fluid type boundaries from the constrained rock physics model template to better understand the rock-fluid relations. The zoned polygons serve as fluid and lithology separator in the studied Wells D, A and C.

Petrophysical analysis
The reservoir petrophysical properties and their average values are presented in Fig. 7 and Table 2. Results showed that the average porosity values for Wells A, C and D are 22.4%, 23.7% and 23.3% with water saturation values of 8.7%, 3.2% and 11.2%, respectively. The result of average volume of shale indicate that 0.02 v/v and 0.052 v/v are calculated for Wells A and D. This value according to Hilchie (1978) is within the limiting value of 0.15v/v for a very good reservoir, thus will not affect the reservoir effective porosity negatively. In Well C, the average volume of shale was calculated as 0.16 v/v which was slightly above the limiting value of 0.15v/v. Observation shows that the average porosity and effective porosity values in Wells A and D did not vary much as both wells are located in the same fault block. In Well C, average porosity and average effective porosity varied due to the significant increase in the volume of shale and change in depositional environment caused by faulting activities as the well is sited in different fault blocks (Table 2 and Fig. 2b).

Well log rock physics template and crossplot analysis
Rock physics template crossplots and analysis of Wells A, D and C were done using mineral and fluid properties as listed in Table 3 as an input parameter. Two lithofacies (shale and sandstone) including infill fluids (brine and gas) were zoned on the log and crossplot space following confirmation from gamma ray, resistivities and neutron-density logs. Figure 8 highlights the lithofacies and fluid types distribution at various depths of Wells A, D and C. The highlighted grey, red and blue colour polygons in Fig. 8a, 8b and 8c are lithofacies clusters containing shale and sandstone filled with hydrocarbon gas and brine. Results show that the red polygons are characterized by low gamma ray values which represent a clean and unconsolidated reservoir with less cement on the grains. Analysis indicated that the modified Hashin-Shtrikman

Reservoir deposits
Fluid type   upper bound line successfully marked the gas-brine boundary on the clean and unconsolidated reservoir due to its sensitivity to fluid densities. These were confirmed using neutron-density crossover and high resistivity values from logs to match the gas and brine-filled zones ( Fig. 8a  and 8b). Observably in Well C (Fig. 8c), the modified Hashin-Shtrikman upper bound line was successful in confirming the gas-brine boundary even in the absence of the neutron porosity log. Crossplot results indicate that high gamma ray sediment clusters on the low-porosity end (Reuss lower bound line) as shown in Fig. 8a and 8b, emphasizing the predominance of clay particles which is an indication that the rock materials bounding the hydrocarbon zones are stiff as a result of compaction. In Fig. 8c, the high gamma ray sediment clusters on the Reuss lower bound line were caused by the compaction process and presence of feldspathic sands associated with shelf-edge deltas. This effect alters the elastic properties of the rock materials, net to gross, reduces sorting and porosity of the reservoir. In Well C, the gas-filled zone revealed fair reservoir quality due to poor sorting with cement in the grains of the reservoir. This is evidenced in the clusters away from the high porosity end (above the modified Hashin-Shtrikman upper bound line). In general, results revealed that the quality (porosity) of the reservoir which is a function of sorting, grain contact and compaction gradually increases above critical porosities (0.36-0.40) in Wells A and D as shown in Fig. 8a and 8b. Although their average porosities (Table 2) indicate very good quality, observations showed that the Ozifa reservoir qualities in Wells A, D and C decrease gradually as cementation increases towards the Reuss lower bound which mark the boundary for stiffest rock. Results from Fig. 8a and 8b displayed that above critical porosity, clean and unconsolidated sands dominate the reservoir. This is because newly deposited fine-grained sediments that are free of cement act as suspension fallout in the range of critical porosity. This study centres on rock physics template bounds, models application on petrophysical studies to understand the behaviour and interactions of geological variables such as sorting, net-to-gross, grain contacts and saturation on reservoir quality and fluid types discrimination. In summary, Fig. 9 represents the locally derived conceptual model of the rock physics template adopted for successful characterization and discrimination of lithology and fluid types in the study area.

Conclusion
The petrophysical analysis indicated that the Ozifa reservoir sands possess very good pore connectivity, sorting with sufficient hydrocarbon (gas) across the studied wells. The rock physics models (modified Hashin-Shtrikman upper bound, Voigt upper model and Reuss lower bound) used to predict and characterize the litho-fluid relationship were sensitive to factors such as pore fluids, grain contacts, sorting and clay distribution. The pore geometry, net to gross and sorting increase as the quality of the reservoirs (porosity) increases above the modified Hashin-Shtrikman upper bound line. This bound line serves as the key boundary line for brine and gas-filled sands. In conclusion, this study will serve as a model to understand better the effects of moduli generated bounds such as modified Hashin-Shtrikman upper bound on geological variables such as sorting, clay distributions on the reservoir porosity and fluid distribution which in turn affects the reservoir quality. In order words, incorporating this approach during fluid replacement modelling (FRM) as a moduli parameter will effectively reveal more understanding on litho-fluid interactions of clastic reservoirs.