An Improved Method of Clay-Induced Rock Typing Derived from Log Data in Modelling Low Salinity Water Injection: A Case Study on an Oil Field in Indonesia

Abstract

Low salinity water injection (LSWI) is an emerging way to improve waterflood performance through chemical processes. The presence of clay minerals is one of the required parameters to successfully implement LSWI in sandstone formations. The ability of clays to exchange the cations, represented by cation exchange capacity (CEC), leads to oil detachment from the rock surface and changes the formation wettability toward water-wet. There are still limited studies that discuss the implementation of specific CEC models in the field-scale LSWI reservoir simulation. This paper attempts to propose an improved method of clay-induced rock typing that can be representatively implemented for field-scale reservoir simulation. The scope of this study is limited to a sandstone reservoir from an oil field in Indonesia. The oil is considered light, and the reservoir contains main clay minerals, including kaolinite and illite, and a trace of chlorite was also found from the XRD evaluation. CEC can be derived from log data, while rock type can also be estimated from log data by using the artificial neural network method. The main finding is that the combination of those variables, i.e., log data, rock properties, and CEC, results in an improved method to characterize and classify the clay into three types associated with conventional rock types. The classification obtained by the clay typing method can be utilized as an input for advanced LSWI modeling, which is expected to provide more robust results. Furthermore, dispersed clay has a strong influence on the magnitude of cation exchange capacity rather than laminar and structural clays.

1. Introduction

Water injection has been viewed as a proven method to maintain reservoir pressure and displace oil towards production wells, while less attention has been given to the role of the chemistry of the injected water and its impact on oil recovery. Chemically, the salt composition of injected water in conventional waterflood might lead to salt deposition, which potentially reduces the oil production rate [1,2]. Low salinity water injection (LSWI) is an emerging way to improve waterflood performance that exploits the involved chemical processes to the advantage of improving oil recovery. LSWI is an injection of chemically tuned brines into sandstone or carbonate reservoirs [3]. The oil recovery increases due to injected low salinity water through specific mechanisms, which include those chemical processes. Several underlying mechanisms that were published are fines migration [4], pH variation or in situ surfactant generation [5,6], multicomponent ion exchange (MIE) [7], electrical double layer expansion [8], and many other proposed mechanisms. The different proposal regarding the basic mechanisms of LSWI leads to the unique interaction of the crude oil–brine–rock (COBR) system, which cannot be separated from each other during the observation.

The presence of clay minerals is one of the required parameters to implement LSWI in sandstone formation successfully. The ability of clays to exchange the cations composing the brine and inject low salinity water leads to oil detachment from the rock surface and changes the formation wettability toward water-wet [7]. That ability is then represented by the parameter of cation exchange capacity (CEC) [9]. The CEC was previously found to not only depend on the type of clay mineral but also on how those clays are distributed in the formation [10]. Clay minerals can be distributed in the forms of laminar, dispersed, and structural distributions. Those distributions affect the effective porosity of the reservoir [11]. As a result, the form of distributions determines the degree of clay surface exposure that affect the efficiency of LSWI.

The CEC is mainly measured in the laboratory by using the wet chemistry method in the unit of milliequivalents (meq) per formation bulk mass. The sample is crushed into fine grain, and the reactivity of the fines to the chemical (mainly methylene blue) is identified to represent the CEC value. This method basically leads to erroneous interpretation while clay minerals are not fully exposed to fluids within the porous media. The phenomenon is controlled by the distribution of clay in the rock formation. Waxman and Smits [12] observed the deviation of the Archie equation on the shaly sandstone in terms of excess conductance and then determined the CEC as a source of that phenomenon. As the CEC is related to the conductance of the rock sample, they improved the Archie equation, which can overcome the presence of clay in the shaly sandstone sample. Due to limited data from laboratory analysis, other researchers tried to utilize resistivity data from well logs to determine the CEC, then namely log-derived CEC. Juhasz [13] used normalized CEC as represented per pore volume (Q_vn), which is a function of shale volume fraction (V_sh). The Q_vn was then denormalized by multiplying it with that of full shale (V_sh ≈ 1). Researchers from Louisiana State University (LSU) developed an outstanding method to determine the CEC by using the log data from resistivity, spontaneous potential, neutron, and density logs, referred to as the Modified LSU Model [14]. Log-derived CEC determination becomes very useful in its development for not only the petrophysical aspect but also other aspects of petroleum engineering. Demircan et al. [15] utilized this method in drilling optimization. The existence of shale was related to CEC, by which the knowledge of formation characteristics helps prevent ineffective time during either drilling activity or future offset wells.

Several researchers included the CEC model in the LSWI simulation study. Jerauld et al. [16] used kaolinite content to represent clay minerals in simulating LSWI numerically. They measured kaolinite content from different sources, i.e., core, tracer evaluation, etc., as a function of recovery factor to represent relative permeability shift. Dang et al. [9] proposed a new approach to model the clay distribution in the reservoir. Different facies were assigned relating to different grain sizes of the rock. They constructed a clay-dependent grain size in which high clay content belongs to fine grain and vice versa. Moreover, the generalized correlation between CEC and clay content was used in a study conducted by Seilsepour and Rashidi [17]. However, many simulation studies on a laboratory scale focused on geochemistry modeling during LSWI. Omekeh et al. [18] and Fjelde et al. [19] simulated a coreflood experiment that included an ion-exchange model by using a constant CEC value at all times. Dang et al. [20] used a single CEC and ion equivalent fraction as interpolant of relative permeability shift to mimic the coreflood experiment conducted by Fjelde et al. [19].

There are still limited studies that discuss the implementation of specific CEC models in the field-scale LSWI reservoir simulation. Field-scale reservoir simulation inevitably deals with heterogeneous petrophysical properties such as those characterized by rock types. Motivated by the studies regarding the CEC model for LSWI numerical simulation previously elaborated, this paper attempts to propose an improved method of clay-induced rock typing that can be representatively implemented for field-scale reservoir simulation. This work also relates the determined CEC with rock types derived from log data that are mainly used in reservoir simulation. The relationship gives a new extending method of conventional rock typing that can be used for LSWI modeling appropriately in field-scale.

2. Field Overview

S Field is an oil field located in South Sumatera Province, Indonesia. The field was discovered in 1990 with an area of approximately 4510 acres. The field has two producing structures which are BX and TX Structures. This study focuses only on the BX Structure. TBX sand layer, which is the focus of this study, is the biggest and the main producing layer of this structure. The location area is shown in Figure 1. There are three production wells (B-1, B-3, and B-5) currently producing from the structure and one suspended well (B-4).

The active oil production layer in the BX Structure originates from the TBX-3 Layer, which is a sandstone reservoir from the Talang Akar Formation member, as one of the largest oil-bearing formations in the South Sumatera basin, with an average depth of ±1500 m TVDSS with a thickness of ±3–15 m as a deposit product of the fluvio–deltaic depositional environment [21]. It is reported to initially have 1800 psig reservoir pressure and 227 °F reservoir temperature. The reservoir has 0.156 average porosity and 493.7 md average permeability. It contains light oil with an oil gravity of 33.4°API and bubble point pressure of 1205 psi. The drive mechanism of the reservoir is considered to be a strong water drive. The reservoir has formation water with 20,000 ppm of salinity. The reservoir also contains main clay minerals, including kaolinite and illite, and a trace of chlorite was also found from the XRD evaluation (see

Table 1. Reservoir properties, CEC value, and clay minerals evaluation of BX Structure.

Well	Sample Number	Depth (Meter)	Porosity (%)	Permeability (mD)	Grain Density (gr/cc)	CEC (meq/100gr)	Qv (meq/mL)	Total Clay (%)	Mineral Clay Prediction (XRD)
B-4	8	1509.20	13.61	1.92	2.64	1.60	0.27	10.0	Clean
	9	1511.00	17.24	94.87	2.64	1.20	0.15	26.0	Kaolinite, Illite
	13A	1570.00	6.41	0.08	2.65	3.40	1.32	50.0	Kaolinite, Illite
	14A	1571.00	7.64	0.06	2.65	14.40	4.61	32.0	Kaolinite, Illite, Clorite
	19A	1603.40	9.26	1.25	2.66	1.00	0.26	10.0	Clean
	20	1604.20	11.84	75.29	2.64	0.80	0.16	17.0	Clean
B-5	3A	1431.78	15.50	0.38	2.66	2.19	0.32	20.8	Kaolinite, Illite, Clorite
	12C	1522.92	20.00	64.70	2.65	3.91	0.41	14.1	Kaolinite, Illite
	13B	1523.44	22.70	704.00	2.65	1.62	0.15	25.0	Kaolinite, Illite
	14C	1524.73	26.10	2360.00	2.65	1.37	0.10	17.4	Clean

).

3. Methodology

The log-derived CEC method for shaly sandstone proposed by Juhasz [13] was used to determine the CEC value for each well in the BX Structure. The author introduced the concept of normalized CEC per unit total pore volume (Q_vn). Normalized Q_v fills in the gap where there is no specific logging tool to measure the amount of Q_v in the reservoir, especially in the shaly sand formation. Q_v is normally measured in the lab, which requires the acquisition of core. In the absence of or limited core data, Q_v can be determined using the concept of Q_vn proposed by Juhasz [13] as expressed by Equations (1) and (2). The CEC of the shale interval was determined as the maximum Q_v value, while other intervals were denormalized by using Q_vn generated from V_sh data. Q_vn is determined by Equation (1), while the relationship between Q_v and CEC is represented by Equation (2).

$${\mathrm{Q}}_{\mathrm{vn}}=\frac{{\mathrm{V}}_{\mathrm{sh}}/{\mathsf{\varphi }}_{\mathrm{T}}}{{\left({\mathrm{V}}_{\mathrm{sh}}/{\mathsf{\varphi }}_{\mathrm{T}}\right)}_{\mathrm{sh}}}=\frac{{\mathrm{V}}_{\mathrm{sh}}{\mathsf{\varphi }}_{{\mathrm{T}}_{\mathrm{sh}}}}{{\mathsf{\varphi }}_{\mathrm{T}}}=\frac{{\mathrm{Q}}_{\mathrm{v}}}{{\mathrm{Q}}_{{\mathrm{v}}_{\mathrm{sh}}}}$$

(1)

$${\mathrm{Q}}_{\mathrm{v}}=\frac{\mathrm{CEC}\left(1–{\mathsf{\varphi }}_{\mathrm{e}}\right)\mathsf{\rho }}{100\mathsf{\rho }}$$

(2)

ϕ_T and ϕ_e stand for total and effective porosity, respectively, and the subscript “_sh” represents the property of ~100% shale interval.

The concept of clay distribution was used to comprehend the understanding of this study. The basic clay type volume distribution was determined using the rhombus plot proposed by Juhasz [22] as a modification of the previous plot generated by Thomas and Stieber [11], as shown in Figure 2. The generated clay distribution from log data was then validated by the core samples and sediment interpretation. The effect of each clay distribution (Figure 3), i.e., dispersed, structural, and laminar clays, on CEC was analyzed by constructing a curve including V_sh and the distribution type, which is discussed in Section 4.3.

Routine core analysis (RCAL) of the corresponding structure was used to characterize the rock type of the reservoir by using the hydraulic flow unit (HFU) method proposed by Amaefule et al. [23]. Five rock types were classified for the reservoir. Several log data, i.e., resistivity, gamma ray, density, neutron, spontaneous potential, and sonic logs, taken at a depth of each core sample were utilized to generate a correlation for estimating the HFU from un-cored intervals. Guo et al. [24] used the multi-linear regression approach to estimate the HFU from well log data; however, the result may not be quite satisfactory. Thus, an additional approach using a supervised machine learning (ML), the backpropagation of the artificial neural network (ANN) method, possibly yields a better correlation result with cross-validation from testing and validation data from 35 core samples data as a target. The evaluation was performed using MATLAB^TM.

The CEC/Q_v and HFU derived from well log data were then plotted in the same curve. While HFU only considers dynamic properties of porous media, the CEC value characterized at the same HFU should not give a good agreement trendline. According to the identical approach by da Silva et al. [25], V_sh has a good agreement with the CEC empirically. Thus, CEC was then regressed based on the V_sh and the porosity giving the internal CEC model at the same HFU. This extending rock typing is useful for field-scale LSWI simulation. Figure 4 shows the workflow of the proposed method, which was included in that of the conventional reservoir simulation.

The workflow within the red dashed line in Figure 4 shows the procedure of this current work. The result of this work was then combined with a geochemical reaction model that was built in reservoir simulation software to represent the wettability alteration. Basically, the alteration is modeled through relative permeability shift in the account of released cations within the MIE mechanism. This affects the mobility of oil and water when the grid in the reservoir model is invaded by the low salinity water. The proposed method is attached to the dynamic model only through the fluids’ relative permeability shifting, while there is no effect on pressure–volume–temperature (PVT) fluids and reservoir static models.

4. Results and Discussion

4.1. Log-Derived CEC Determination

Determination of log-derived CEC using Juhasz’s [13] method requires log data, i.e., density, neutron, and gamma ray. Gamma ray is used to obtain shale volume fraction, while density and neutron logs are utilized to acquire porosity. Normalized Q_v is initially determined for each depth based on the porosity and the shale volume fraction at the corresponding depth (Equation (1)). The porosity of 100% shale used in this work was 0.07 based on log data at the shale zone. Q_v of the 100% shale was then determined at a number where the calculated Q_v and also the CEC are very well matched with those measured in the laboratory. An amount of 25 meq/mL of shale Q_v was chosen as it yields well estimated Q_v compared with those measured, as shown in Figure 5.

The measured CEC basically represents the bulk property due to the method of how those CEC are measured. Furthermore, log-derived CEC is also determined through normalization based on V_sh. In order to clearly understand this property, clay distribution needs to be identified. The distribution should explain how clay affects the CEC physically through how it is distributed within the pore space. The concept and how to obtain clay distribution are elaborated in the following section.

4.2. Clay Distribution

Porosity and shale volume fractions obtained from log data at each reservoir layer were used in applying Juhasz’s method [22] for the determination of clay distribution (Figure 6). According to the curve, most data lie on the dispersed-laminar clay. These kinds of distribution affect the effective porosity, as shown in Figure 3. The effective porosity is determined based on the portion of each distribution. Qualitatively, the clay distributed in the form of dispersed clay contributes the most to the total porosity reduction, and it is followed by laminar and structural clays, respectively. Outliers are also identified in the curve, which may be caused by the presence of gas during log data acquisition.

Clay distribution results were then assessed with log and core data. The cores were taken from three different zones those are shale dominated, transition, and shaly sand. Core samples, micro-resistivity image log, and gamma ray log data are set on the corresponding depth and basically show good agreement with each other, as shown in Figure 7. According to the clay distribution results, laminar clay distribution is dominant in the zone in which clay is major. These results are confirmed by the core sample in which the clay lays, forming thin layers (lamination) within the formation.

4.3. Relationship between Log-Derived Q_v and Clay Distribution

Determination of clay distribution was performed using the method depicted by Juhasz [22], as shown in Figure 6. The position of each point inside the rhombical area of the plots determines the percentage of clay distribution-type existence, namely dispersed, laminar, structural, or combination. The combination meant by this method is either laminar-dispersed or laminar-structural clays. The boundary lines of the rhombical area declare either fully dispersed or fully laminar type of clay, while any points located within means the combination with a certain percentage of each distribution type. Those plots represent depth points within the interval of the depth of observation. Furthermore, in order to evaluate the relationship between CEC and clay distribution, Figure 8 is constructed, which expresses the shale volume in different distribution type compositions and its corresponding Q_v.

The value of Q_v does not necessarily rely on the shale volume contained in the interval. As shown in Figure 8, there are data points that have relatively low shale volume (20–40%) but have a considerable magnitude of Q_v. Those are symbolized by the red color in Figure 8a, in which dispersed clay dominates the distribution. Dispersed clay conceptually fills the pore space and decreases the porosity. Therefore, the portion of dispersed clay is entirely present on the surface of the porous media to contact with brine and hydrocarbon fluids. Furthermore, Figure 8c reveals that structural clay contributes the least towards the magnitude of Q_v. Although structural clay is also present entirely within the pore space, the high shale volume dominated by structural clay yields low Q_v. Additionally, Figure 8b shows the laminar shale somewhat dominates the clay distribution of the interval, and its effect on Q_v magnitude is not clearly identified. However, laminar clay is generally aligned with the V_sh.

Another remarkable finding is that laminar clay does not associate directly with the decrease in the total porosity, as depicted in Figure 9b. On the other hand, a good agreement is shown by dispersed clay to that of the total porosity (Figure 9a). The higher the portion of dispersed clay, the lower the total porosity. There is no specific correlation that can be identified between structural clay and porosity (Figure 9c). All the above findings confirm the concept of each type of clay distribution (Figure 3). Considering the relationship of V_sh and laminar clay represented in Figure 8b, the increase in V_sh can be described by the arrow sign in Figure 9b. Although V_sh reaches the value of 1, it does not mean that Q_v achieves the highest magnitude. This is possibly due to the different types of clay minerals.

4.4. Log-Derived Rock Type Determination

Log-derived rock typing was performed to determine the rock type at each depth based on core samples processed for RCAL. Conventional RCAL using the HFU method proposed by Amaefule et al. [23] was conducted for native core samples in the objective field resulting in five rock types. Log data such as resistivity, gamma ray, density, neutron, spontaneous potential, and sonic logs were then prepared for the rock type estimation of each depth.

The first attempt was to estimate the discrete rock type (DRT) by utilizing multi-linear regression, as was conducted by Guo et al. [24]. Each log data were initially normalized to give a consistent range of and sensitive value for the entire log data so that the value lies in the range of zero to one. The normalization follows the formula written below:

$$\mathrm{Nx}=\frac{\mathrm{x}/{\mathrm{x}}_{\min }}{{\mathrm{x}}_{\max }/{\mathrm{x}}_{\min }}$$

(3)

where x is any log data that are normalized, x_min and x_max are the lowest and highest value identified from x log, respectively, and Nx is the normalized value of x log at each depth. On the other hand, the basic equation for DRT estimation through multi-linear regression is formulated as follows:

$${\mathrm{DRT}=\mathsf{\lambda }}_0{+\mathsf{\lambda }}_1{\mathrm{Nx}}_{\mathrm{RD}}{+\mathsf{\lambda }}_2{\mathrm{Nx}}_{\mathrm{GR}}{+\mathsf{\lambda }}_3{\mathrm{Nx}}_{\mathrm{RHO}}{+\mathsf{\lambda }}_4{\mathrm{Nx}}_{\mathrm{NPHI}}{+\mathsf{\lambda }}_5{\mathrm{Nx}}_{\mathrm{SP}}{+\mathsf{\lambda }}_6{\mathrm{Nx}}_{\mathrm{DT}}$$

(4)

where $\mathsf{\lambda }$s are regression coefficients, and subscripts of RD, GR, RHO, NPHI, SP, and DT refer to resistivity, gamma ray, density, neutron porosity, spontaneous potential, and sonic logs, respectively.

The results of utilizing the multi-linear regression method do not give a good estimation of DRT, as depicted in Figure 10. The points located out of the colored boxes in Figure 10 indicate missed classification. The coefficients for each log data resulting from this method are listed in

Table 2. Coefficients resulted from multi-linear regression method.

		NxRD	NxGR	NxRHO	NxNPHI	NxSP	NxDT
Coeff. label	${\mathsf{\lambda }}_0$	${\mathsf{\lambda }}_1$	${\mathsf{\lambda }}_2$	${\mathsf{\lambda }}_3$	${\mathsf{\lambda }}_4$	${\mathsf{\lambda }}_5$	${\mathsf{\lambda }}_6$
Coeff. value	9.766	15.894	−6.532	−4.867	11.023	0.676	7.177

. This result triggers the more advanced method required for DRT prediction instead of using a multi-linear regression approach.

Machine learning is promoted for learning the correlation of several bundles of data with the undefined association, as encountered in this case. The ANN method, which is one of the machine learning tools, was used as the second attempt to estimate the DRT. The neural network itself is not an algorithm but rather a framework for many different machine learning algorithms to work together and process complex data inputs.

This study used the backpropagation method of ANN, which is a systematic method of artificial neural network that uses algorithm supervised learning and is usually used by perceptrons with multiple screen layers for changing the weights in the hidden layer. In the backpropagation, the training used is the type of control that uses a weight adjustment pattern to achieve the minimum error value between the output results prediction of real output. This work used seven log data as input layer and normalized them using Equation (3); previously, the feature engineering using heatmap and Pearson correlation were examined. As shown in Figure 11, there are two hidden layers. The first hidden layer, labeled as layer-I, is designated with 30 hidden layers; the second layer, labeled as layer-j, consists of 15 hidden layers; and a single output layer as the target data is the HFU/DRT. In this case, the available core dataset involved is 35 core samples. The whole data were used during the ANN learning process. The proportion of training and testing data is 70:30 in the percentage of total data used.

Based on 35 samples gathered from RCAL, they are divided into three functions, which are train, test, and validation data. In addition to log data, HFU/DRT was also normalized before the ANN was performed. This ANN method yields a good estimation of the core data, as shown in Figure 12a. The result was matched better than that resulting from multi-linear regression. Moreover, there was no missed classification of five rock types as previously performed by the first attempt. Thus, the estimated DRT through the ANN method was used for the following analysis.

4.5. Clay Induced Rock Typing

Clay typing was proposed in this work after generating two key variables, i.e., log-derived Q_v and log-derived DRT. In order to identify the correlation between those variables, a cartesian plot of Q_v vs. $\mathsf{\varphi }$_T was plotted in which the rock type number was represented by a color scale. Based on Figure 13a, rock types 4 and 5 were widely distributed along the fraction of $\mathsf{\varphi }$_T. A reasonable margin of Q_v at a single value of $\mathsf{\varphi }$_T was also identified within those rock types, especially in a range of 0.01 up to 0.125 of porosity. This is caused by the fraction of shale controlling the magnitude of Q_v, as depicted in Figure 13b. According to the log-generated data, none of the data of rock types 1–3 were identified below 0.125 of porosity. That does not mean that these rock types do not comprise the porosity below 0.125 in the reservoir model, as this work only proceeded with log data from several wells.

In order to comply with both considerations, shale volume fraction and generalizing the concept of rock properties were included in the clay typing. The shale volume fraction was divided into three classifications, i.e., 0–0.125 as group 1, 0.126–0.875 as group 2, and 0.876–1 as group 3, to overcome the margin of Q_v in rock types 4 and 5. Each group then represented clay type 1 up to 3, as depicted in Figure 14. The correlations describing the clay types are listed in

Table 3. Correlations and remarks of clay typing for each rock type.

	Correlations	Implementation in RT 1–3	Implementation in RT 4–5	RMSE
Clay type 1	${\mathrm{Q}}_{\mathrm{v}}=0.0443{\mathsf{\varphi }}_{\mathrm{T}}{}^{-1.202}$	-	Vsh: 0–0.125	0.289
Clay type 2	${\mathrm{Q}}_{\mathrm{v}}=0.1130{\mathsf{\varphi }}_{\mathrm{T}}{}^{-1.153}$	All range of Vsh	Vsh: 0.126–0.875	1.356
Clay type 3	${\mathrm{Q}}_{\mathrm{v}}=0.2521{\mathsf{\varphi }}_{\mathrm{T}}{}^{-1.076}$	-	Vsh: 0.876–1	1.289

. Those are the best correlation that can be achieved regardless of the root mean square error (RMSE) values, considering the wide range of available field data. This classification applies only to rock types 4 and 5 as the worst rock types. Additionally, the effect of shale volume fraction is insignificant on the Q_v for rock types 1 to 3. As clay may also be minor in these rock types, a single clay model was preferred instead of using the model that is used for the worst rock types. By considering the statistical data distribution, clay type 2 was chosen for clay typing of rock types 1 to 3 and also for the extended model of those rock types for a porosity below 0.125.

Clay typing in this work is basically a function of two main variables: the formation shale volume and the porosity. Each type is assigned to the conventional rock types based on its specific conditions. The clay typing associated with the rock types simplifies the application of this method in reservoir simulation as rock types are generally well defined. The implementation is by simply defining a new property named “clay types” and assigning formulas with certain conditions corresponding to those that are listed in Table 3. Further procedure and application of the proposed method in reservoir simulations will be published following this article.

5. Conclusions

The results of this study have led to the following conclusions:

• An improved method was developed to characterize and classify the clay types by using the combination of log-derived cation exchange capacity and conventional rock types;
• The classification obtained by the clay typing method developed in this work can be utilized as an input for advanced modeling of low salinity water injection (LSWI) to provide more robust results. Previous studies simply used cation exchange capacity as a single input or were defined by empirical correlation;
• This work identifies that dispersed clay has a strong influence on the magnitude of cation exchange capacity rather than laminar and structural clays. Further investigation is strongly recommended to quantify the effects of clay distribution types and minerals on the clay typing method.