Abstract

The optimization method to evaluate the tight reservoir porosity is a difficult technique to use due to its complexity and instability. This paper proposes an improved optimization method to calculate the porosity of tight reservoirs. First, we applied the matrix model which modified the multicomponent model to the problem and it improved the results by deducing a mathematical model. Second, we used the Simulated Annealing Algorithm to calculate the incoherence function, and then based on statistical theory, we obtained the most optimal results. Examples show that the method is effective, and despite the lack of the local experience parameters, its application is valuable in order to evaluate the porosity.

1. Introduction

The optimization log interpretation method provides a variety of logging information which is an effective way to evaluate complex oil and gas reservoirs. After it was proposed in 1980 [1], optimization log interpretation method technology has developed and improved. Although the firm Schlumberger Ltd. boasted about how the technology could be suitable for all types of reservoirs around the world, it was difficult to obtain a satisfactory result when evaluating tight reservoirs. Because the pore system in tight reservoirs presents a much greater variety of geometrical characteristics and structures compared with clastic reservoir, then the complexity of the pore system affects the results because the parameter is difficult to ascertain. Therefore, the development and improvement of the optimization log interpretation method to solve the above problem are required.

The traditional optimization log interpretation method contains two important processes: () the establishment of the pore model and the response equations and () the calculation of the response equations.

The traditional pore model is a multimineral model or multicomponent model [2, 3]. In the model, the reservoirs act as a unit of a space geological body and are made up of partially homogeneous components, such as minerals and porosity. The purpose of this division is to simplify the calculation process and programming design [4]. However, arguments presuppose that all types of minerals and geologic parameters in multimineral or multicomponent models must be known. Those models cannot be used in geographical areas where the mineral components are difficult to distinguish, especially at tight reservoirs. In this paper, we put forward the mixed-matrix model to improve traditional model following the Fuzzy Clustering analysis. An important characteristic of this new model is synthesizing multimineral as a synthesis. For instance, with this model, the user does not need to know the physical parameters of all the types of mineral present in the system. Therefore, the response equations are calculated for the optimization log interpretation mathematical model.

The calculation of the response equations is the optimization algorithm. In the process of algorithm research, a large number of scholars have used different algorithms to propose mathematical models as it is shown in Table 1 [5]. For example, Global Program started using the steepest descent algorithm [1], and Optima Program used the conjugate gradient algorithm [6]. These programs were a good approach evaluating conventional reservoir wells. In China, Ouyang et al. used the simplex algorithm to evaluate clastic reservoir [7] and lamprophyre reservoir [8], using BFGS or DFP to evaluate oil-bearing igneous for better results [9]. All of the above optimization algorithms have an initial value and the initial value greatly influenced the results [10]. Stoffa and Sen introduced the genetic algorithm (GA), in the optimization algorithms [11], and GA does not need the initial value of the input which means GA is not restricted by the initial value [12]. Unfortunately, GA was not stable method to use due to the repeated calculation it involves. An important improvement in the algorithm would be to eliminate the influence of the initial value and the computational stability. Kirkpatrick et al. used the simulated annealing to solve a similar problem [13]. In this paper, we apply the Simulated Annealing Algorithm (SAA) to calculate the response equations, because not only is it a global optimization algorithm for functions of continuous variables but also it depends slightly on the initial value [14].

Firstly, we start comparing the correlation between well logging and core porosity to select the well logging types. Secondly, we have an introduction on mixed-matrix model to improve multimineral model. Thirdly, we use SAA for optimization and then we use the Monte Carlo random method to select the initial parameters of the optimal values. The method was then tested with set data from a tight reservoir. Finally, we conclude with a discussion of advantages and limitations of the method.

2. The Improved Optimization Method

2.1. The Selection of Sensitive Logging Curves

The use of more logging curves is likely to prove more reliable results [15]. However, it becomes less reliable with increasing the response equations and increasing the parameters [16]. Thus, it is crucial to select the more correlated logging curves as possible.

The porosity logs available for this study were obtained from ten wells, the parameters cover acoustic (AC), density (DEN), compensated neutron log (CNL), and natural gamma rays (GR), among others. The porosity datum from the ten wells come from 1773 samples. The data was obtained from tight sandstone and shales reservoir. And then correlations were calculated using the core porosity data and AC, DEN, CNL, and GR as shown in Figures 1(a)1(d). The overall porosity ranged from 0% to 12% (average is 6.48%), the overall AC is within the scope of 156.8 us/m to 338.9 us/m (average is 225.5 us/m), the overall DEN was in a scale of 2.17 g/cm3 to 2.99 g/cm3 (average is 2.54 g/cm3), the overall CNL was from 2.9% to 39.3% (average is 15.9%), and the overall GR ranged from 10 API to 330 API (average is 90 API).

In Figures 1(a)1(d), the correlation coefficients between the core porosity and the well log value of AC, DEN, CNL, and GR are 0.1854, 0.2198, 0.2634, and 0.069, respectively. The AC, DEN, and CNL (triporosity loggings) showed a higher correlation coefficient. In this paper, we selected AC, DEN, and CNL to calculate the porosity, similar to the usual optimization interpretation used in conventional reservoir.

2.2. The Establishment of the Mixed-Matrix Model

In order to solve the multimineral model or the multicomponent model, which cannot be used in geographical areas where the mineral components are difficult to distinguish especially the tight reservoir, we modified multicomponent model to the mixed-matrix model in this research. Figure 2 shows the following:(1)In the model, the complex formation as a unit of a space geological body is made up of three components, which are , , and .(2) is the mixed-brittle matrix composed of , , , , and so on.(3) is the mixed-clay matrix mainly composed of , , , and so on.(4) is total porosity mainly composed of and .

Three kinds of component are determined to reduce the volume for the target porosity to build the linear equation

2.3. The Establishment of Incoherence Function and Constraints

The building of optimization log interpretation mathematical model is the basis of the optimization log interpretation, including the establishment of incoherence function and constraints.

According to the principle of least squares methods, the incoherence function of optimization interpretation is established as

Extraordinarily, each logging data has different order of magnitude, so the data should be processed. The procession is prior to normalization, which is for the different measuring dimension. The incoherence equation (2) turns into

After establishing the incoherence function, the constraints also need to be defined to ensure the practical significance. Its constrains are shown on

2.4. The Selection of the Geologic Parameters in Mixed-Matrix Model

The conventional reservoir matrix is mainly composed of quartz, feldspar, calcite, and dolomite; however, the tight reservoir matrix composition is hard to be determined. In mixed-matrix model, we modified the components as , , and . The logging parameters of each component are no longer a certain value but an interval, as shown in Table 2. For example, the value of is determined by the relative content of quartz, feldspar, calcite, dolomite, and so on. At the same time, the minimum incoherence function was obtained like all the other parameters; therefore, the components and physical parameters can be calculated.

2.5. The Simulated Annealing Algorithm and Monte Carlo Initial Value

Szucs and Civan used the simulated annealing method to interpret the multilayer well log [17]; however, the method required layering of reservoirs in advance. In order to find approximate minimum of the incoherence function, we present the Simulated Annealing Algorithm (SSA) combined with Monte Carlo random method to find approximate results. SAA is essentially an iterative random search procedure with adaptive moves along the coordinate directions, depending only slightly on the starting point. It allows uphill moves under the control of a probabilistic criterion, thus tending to avoid the first local minima encountered.

As is shown in Figure 3, SSA proceeds iteratively: starting from a given point which is randomly generated by Monte Carlo random method, and it generates a succession of points: , tending to the global minimum of the incoherence function. New candidate points are generated around the current point applying random moves along each coordinate direction, in turn. If the point falls outside the definition domain, a new point is randomly generated by Monte Carlo random method until a point belonging to the definition domain is found. A candidate point is accepted or rejected according to the Metropolis criterion [18].If , then accept the new point: ;else accept the new point with probability:,where , is a parameter called temperature and the value must be higher than 0, is the termination cycle, is temperature reduction, and is the error tolerance for the function.

We wrote the optimization algorithm using MATLAB software according to the schematic diagram on Figure 3 [19].

2.6. The Calculation Results of Normal Distribution

Four samples were chosen from the 1773 samples (Table 3) to illustrate and eliminate the differences. In Table 3, each of them was repeated 100 times as shown in Figure 4.

In Figure 4, the abscissa is the calculation porosity and the ordinate is the frequency of the results. Figures 4(a)4(d) present that are 2.77%, 3.99%, 4.38%, and 5.74%, respectively, are 2.78%, 3.95%, 4.34%, and 5.75% respectively, and are 1.47%, 1.86%, 1.88%, and 2.11%, respectively, which means that is very close to ; therefore, is almost half of . The calculated results trend is a normal-like distribution. According to the Pauta criterion, theoretically, the probability of one to three standard deviations is 68.3%, 95.4%, and 99.7% in Table 4.

2.7. The Iterate Index

The normal-like distribution will help eliminate the differences or make computational stability. We can calculate many times and average the results to make computational stability. Apparently, the iterate index can lead to prolonging computation time. Therefore, it is crucial to determine the iterate index of calculation. We chose X5 well 3671 m–3679 m for the test.

Figure 5 presents that the iterate indexes are 5, 30, 50, and 100. It can be noted that (a) tends to be stable gradually with increasing the iterate index. (b) The stable tends to with increasing the iterate index. (c) does not change much as the iterate index is up to approximately 30.

Figure 6 shows that it has been calculated 5 times when the iterate index is 30. Figure 6 shows the repeated calculation using the previously mentioned approach. In this case, each calculation is not the same but the overall trend is resemblance.

In Figure 7, the four samples were designed for the convergence which were calculated 100 times. Figures 7(a)7(d) present that convergences are 2.78%, 3.95%, 4.34%, and 5.75%, respectively.

3. Applications and Contrast

3.1. The Influence of Initial Value

Using the last sampling point to be the next point initial value is helpful technique used by different studies, for instance, Yong and Sun [20]. Figure 8 shows the results of using Yong’s initial method for the 1773 samples in this study (Figure 8(a)) and also by Monte Carlo random method (Figure 8(b)); the iterate index used was 30 for the graphics, where is the calculation porosity by Yong’s method and is the calculation porosity by Monte Carlo random method.

The abscissa is the calculation porosity and the ordinate is , as shown in Figure 7. It is clear that the correlation coefficient is 0.58 and it increases to 0.8261 from to . Therefore, we believe that Monte Carlo random method is better than Yong’s method.

3.2. The Contrast with Other Methods

Mixed-matrix model is based on multicomponent model. The contrast of the two types of models has great significance due to their strong connection. We chose JY1 well 2330 m–2352 m which belongs to the shale gas reservoir. We used Wyllie formula, multicomponent model, and our mixed-matrix model to calculate the porosity as shown in Figure 9.

From Table 5, the results differ slightly with the core which is due mainly to the 2340 m interval. Figure 9 reveals the different porosity by different methods, and mixed-matrix method is better than multicomponent model and multicomponent model is better than Wyllie formula.

3.3. The Physical Parameters

As it is shown in (1), there are two sorts of the parameters with three kinds of components (, , and ) and also the physical parameters (,,,, , , , , and ). It is also significant regarding the other parameters except for . We chose JY1 well 2330 m–2352 m for the test and the iterate index is 30.

Table 6 reveals that , , and can be compared with the core, and Figure 10 presents that the physical parameters conform to the actual situation.

4. Discussions and Conclusions

4.1. Mixed-Matrix Model from the Original Model

Tight reservoirs study is complex; they are composed of many kinds of minerals, which are difficult to be determined. Various factors should be considered in the establishment of the reservoir model, so the multicomponent model should be modified to mixed-model, in order to obtain more accurate reservoir parameters.

Mixed-matrix model is based on the ideas of Fuzzy Clustering. Comparing the classic multicomponent or multimineral models, the brittleness minerals are fuzzy as mixed-brittle matrix, the clay minerals are fuzzy as mixed-clay matrix, and the oil and gas and water porosity are fuzzy as the pore space volume. This improvement makes the approach easier and the user did not require all types of mineral’s physical parameters.

4.2. Simulated Annealing Algorithm and the Normal Distribution

The Simulated Annealing Algorithm can provide a high reliability in the minimization of incoherence function. It does not guarantee, of course, finding the global minimum of the incoherence function, but if the function has many good near-optimal solutions, it should find one. It is worth noting that the results of the calculation are slightly different each time, so we should use the normal-like distribution to eliminate those differences and achieve computational stability.

4.3. The Advantages and Limitations

All methods provide advantages as well as limitations which will be discussed next. The advantages of the method are the following: (a) the porosity can be evaluated well even if there is a lack of the regional geological parameters; (b) matrix physical parameters can be accurately obtained; (c) the result of calculation is different for each time, but it has consistent trend in the interval of computation; (d) the user is not required to compute all types of mineral, for instance, physical parameters, so the approach is easier.

On the other hand, the limitations are as follows: for individual well, iterating index of 30 may be not enough to meet the requirement. If the case happened, it will take more time to find the final answer.

Nomenclature

:Volume fraction of the mixed-brittle matrix
:Volume fraction of the mixed-clay matrix
:Porosity or the pore space volume
:Volume fraction of quartz
:Volume fraction of feldspar
:Volume fraction of calcite
:Volume fraction of dolomite
:Volume fraction of kaolinite
:Volume fraction of montmorillonite
:Volume fraction of chlorite
:Volume fraction of oil and gas
:Volume fraction of water
:Acoustic
:Density
:Neutron
:Acoustic of the mixed-brittle matrix
:Acoustic of the mixed-clay matrix
:Acoustic of fluid
:Density of matrix
:Density of clay
:Density of fluid
:Neutron of matrix
:Neutron of clay
:Neutron of fluid
:A vector, , ,
:The incoherence function
:The normalized incoherence function
:th logging curve which is , and
:th tool response function
:The difference between and
:A parameter called temperature except 0
:The termination cycle
:Temperature reduction
:The error tolerance for function
:The core porosity
:The average porosity
:The standard deviation porosity
:The calculation porosity by Yong’s method
:The calculation porosity by Monte Carlo random initial value method
:The minimum of volume fraction of the mixed-brittle matrix
:The maximum of volume fraction of the mixed-brittle matrix
:The minimum of volume fraction of the mixed-clay matrix
:The maximum of volume fraction of the mixed-clay matrix
:The minimum of porosity
:The maximum of porosity.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This study was financially supported by () the National Natural Science Foundation of China, (Grant no. 41302107); () the Ministry of Land and Resources Special Funds for Scientific Research on Public Cause (Grant no. 201311107); () the Fundamental Research Funds for the Central Universities (Grant no. 2652011282); and () the CNPC Innovation Foundation (Grant no. 2012D-5006-0103).