A Developed Optimization Method of Tight Reservoir Porosity

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.


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: (1) the establishment of the pore model and the response equations and (2) 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.

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/cm 3 to 2.99 g/cm 3 (average is 2.54 g/cm 3 ), 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.

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  ma ,  cl , and .
(2)  ma is the mixed-brittle matrix composed of  QUA ,  FEL ,  CAL ,  DOL , and so on.
(3)  cl is the mixed-clay matrix mainly composed of  KAO ,  MON ,  CHL , and so on.
(4)  is total porosity mainly composed of  OG and  W .
Three kinds of component are determined to reduce the volume for the target porosity to build the linear equation ] . (1)

The Establishment of Incoherence
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

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  ma ,  cl , 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 Δ ma 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.

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 V 0 which is randomly generated by Monte Carlo random method, and it generates a succession of points: V 0 , V 1 , . . ., V  , . .., tending to the global minimum of the incoherence function.New candidate points are generated around the current point V  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 V  is accepted or rejected according to the Metropolis criterion [18].
If Δ ≤ 0, then accept the new point: V +1 = V  ; else accept the new point with probability:  Δ =  −Δ/ , where Δ = (V  )−(V  ),  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].

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.

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) Por ave tends to be stable gradually with increasing the iterate index.(b) The stable Por ave tends to Por core with increasing the iterate index.(c) Por ave 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.

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 abscissa is the calculation porosity and the ordinate is Por core , as shown in Figure 7.It is clear that the correlation coefficient is 0.58 and it increases to 0.8261 from Por Y to Por MC .Therefore, we believe that Monte Carlo random method is better than Yong's method.

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 mixedmatrix method is better than multicomponent model and multicomponent model is better than Wyllie formula.

The Physical Parameters.
As it is shown in (1), there are two sorts of the parameters with three kinds of components ( ma ,  cl , and ) and also the physical parameters (Δ ma ,  ma , Φ Nma , Δ cl ,  cl , Φ Ncl , Δ  ,   , and Φ N ).It is also significant regarding the other parameters except for .We   10 presents that the physical parameters V(Δ ma ,  ma , Φ Nma , Δ cl ,  cl , Φ Ncl , Δ  ,   , Φ N ) conform to the actual situation.

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 mixedbrittle 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.

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.

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.

Figure 1 :
Figure 1: The crossplot between core data and well logging.(a) The crossplot of core porosity and acoustic.(b) The crossplot of core porosity and density.(c) The crossplot of core porosity and neutron.(d) The crossplot of core porosity and natural gamma ray.

Figure 2 :
Figure 2: Mixed-matrix model based on the improved multimineral model.
of random moves, each along one coordinate direction.Accept or reject each point according to the Metropolis criterion.Record the optimum point reached so far.Rest cycles = 0 Reduce T Reset adjustments = 0 End

Figure 3 :
Figure 3: The schematic diagram of SSA minimization.

Figure 4 :
Figure 4: The calculation results of the statistical distribution.(a) The statistical distribution of point a.(b) The statistical distribution of point b.(c) The statistical distribution of point c.(d) The statistical distribution of point d.

3 Figure 5 :Figure 6 :
Figure 5: The influence of the cycle index for calculation.

Figure 7 :Figure 8 :Figure 10 :
Figure 7: The 100 times calculations convergence.(a) The convergence of point a.(b) The convergence of point b.(c) The convergence of point c.(d) The convergence of point d.

Nomenclature
ma : Volume fraction of the mixed-brittle matrix  cl : Volume fraction of the mixed-clay matrix : Porosity or the pore space volume  QUA : Volume fraction of quartz  FEL : Volume fraction of feldspar  CAL : Volume fraction of calcite  DOL : Volume fraction of dolomite  KAO : Volume fraction of kaolinite  MON : Volume fraction of montmorillonite  CHL : Volume fraction of chlorite  OG : Volume fraction of oil and gas  W : Volume fraction of water Δ: Acoustic : Density Φ N : Neutron Δ ma : Acoustic of the mixed-brittle matrix Δ cl : Acoustic of the mixed-clay matrix Δ  : Acousticof fluid  ma : Densityofmatrix  cl : Densityofclay   : Densityoffluid

Table 2 :
The physical parameter range.

Table 3 :
The samples of triporosity logging data and calculated results.

Table 4 :
The distribution of calculation results.

Table 5 :
The contrast with other methods.

Table 6 :
The range of calculation results and the core.Table6reveals that  ma ,  cl , and  can be compared with the core, and Figure