Optimizing and accelerating history matching progress of numerical reservoir simulation by using material balance analysis

In view of the problem of huge computations and multiple solutions, a method for optimizing and accelerating the progress of history matching by using material balance analysis was established. Based on simple measuring data, such as average pressure, daily production rate and PVT test, a material balance model can be constructed for a three phase system of hydrocarbon reservoir. According to the material balance calculation, original oil in place and aquifer parameters can be obtained as optimized parameters for numerical simulation, which contributes to solve the multi-solutions problem of history matching. Besides, by conducting uncertainty analysis using Monte-Carlo simulation method, one can determine the adjustment range of input parameters, which helps the engineer definite a clear direction, reduce the number of computations and consequently accelerate the progress of history matching. The reliability of the proposed method was verified according to two examples.


Introduction
The advantage of reservoir numerical simulation is that it can integrate different types of data, such as geology, geophysics, reservoir fluid and rock properties, into one platform for production prediction. There is a great deal of uncertainty in reservoir model, which needs to be modified by history matching. History matching has a huge amount of computation and the results have multiple solutions, If a reasonable initial condition and adjustable range of parameters can be screened out before history matching, the uncertainty can be greatly reduced, the process of matching can be optimized and more reasonable combination of reservoir parameters can be obtained.
Before the rapid development and wide application of reservoir numerical simulation technology, material balance method was used in reservoir development management as an important tool for reservoir performance analysis. With the continuous development of computer technology, the material balance method based on zero dimensional model is gradually replaced by exquisite three-dimensional reservoir numerical simulation. This is mainly because the principle of material balance is applied to every grid scale in reservoir numerical simulation, and it seems that separate material balance analysis has become unnecessary. However, compared with numerical simulation, the material balance analysis method has the following advantages: (1) There are many uncertainties in 3D reservoir model, which are related to reservoir shape, petrology, geophysics and so on. However, the material balance equation is only based on the conservation of production, pressure and PVT data, and its uncertainty is relatively small.
(2) It is impossible to obtain the contribution rate of different driving energy to production only by 3D reservoir model, and the influence degree of different model parameters on matching indexes can not be judged in advance. Parameter adjustment is often decided by experience and personal judgment. By calculating the driving index, the material balance analysis method can quantitatively give the different driving energy sources in the oil and gas production process.
(3) The material balance analysis method studies the relationship between the geological reserves, the remaining reserves and the recovery of the reservoir. It does not involve the flow direction of the fluid, and does not consider the spatial changes of rock and fluid properties, the hydrodynamic factors of fluid flow in porous media, fluid differentiation, reservoir geometry, well location or production of different fluids.
Through the above analysis, we think that before the numerical simulation, some valuable information can be obtained through the material balance analysis, which can be used to guide the reservoir numerical simulation. Therefore, this paper establishes an accelerated method of historical matching optimization of reservoir numerical simulation based on material balance analysis. Through the research of predicting global sensitivity parameters and uncertainty analysis, it forms the overall macro understanding of the reservoir and accelerates the process of historical matching.

Material balance equation
Since Schilthuis R J [1]first established the material balance equation of reservoir based on the principle of material conservation in 1936, it has been widely applied and developed in reservoir engineering. It is one of the basic methods for reservoir engineers to interpret and predict reservoir dynamics.

General equation of material balance
The linear form of the material balance equation is where : where N is initial oil reserves, sm 3 ;m is gas cap ratio, production, cumulative gas production, cumulative water injection, cumulative gas injection, cumulative water influx in turn, sm 3 ; i p is initial reservoir pressure, MPa; p is average reservoir pressure, MPa; o B is crude oil volume factor,rm 3 /sm 3 ; s R is dissolved gas oil ratio,sm 3 /sm 3 ; g B is gas volume factor,rm 3 /sm 3 ; w C is compressibility of formation water,MPa -1 ; f C is rock compressibility,MPa -1 ; wc S is irreducible water saturation.
In general, the production data can be obtained directly from the field, the average formation pressure can be measured by the downhole pressure gauge, and the pressurevolume-temperature parameters of fluids and rocks can be measured by the laboratory. The original oil in place N ,gas cap ratio m and aquifer influx e W is unknown. Among those variables, the aquifer influx e W is the function of aquifer radius a r ,thickness h ,permeability a k ,total compressibility t c ,Water influx Angle θ and so on. The aquifer influx is the very important part of material balance equation. There are various models and methods to calculate water influx [1][2][3][4][5][6]. For convenience, we use two common analytical aquifers of Eclipse in this paper. They are Carter-Tracy unsteady aquifer [5] and Fetkvoich pseudo-steady aquifer [6].

Water influx calculate model
The linear form of the material balance equation (1) Carter-Tracy unsteady aquifer In 1960, Carter-Tracy put forward the approximate calculation formula of water influx based on the Hurst-van Everdingen aquifer model. the cumulative water influx in n t is : where B is influx constant, rm 3 /day; D t is dimensionless time; n p ∆ is overall pressure drop in n t ,MPa; is dimensionless pressure; is dimensionless derivative of pressure; which is calculated by the follow formula : (2) Fetkvoich pseudo-steady aquifer In 1971, Fetkovich put forward a new approximate calculation formula of water influx. the water influx formula in the n th time step is: So, the cumulative water influx in n t is : In this formula, ( ) 1 a n p − is the average aquifer pressure in is the reservoir pressure in n t ,MPa; J is influx index,rm 3 /day/MPa; ei W is the maximum water influx,sm 3 ;they can be calculated by the follow formula:

Application of material balance analysis
Based on the research of domestic and foreign scholars on material balance analysis, this paper proposes an idea and method that combines the classical material balance analysis method with reservoir numerical simulation technology: (1) global sensitivity parameters such as reserves and water parameters are obtained through material balance calculation, which play a decisive role in the overall dynamic characteristics of the reservoir, and can be used as the initial optimization parameter combination of reservoir numerical simulation (2) analyze the uncertainty of the parameters, reduce the adjustment range of some adjustable parameters in the history matching, and determine the reasonable adjustment range.

Global sensitivity parameter estimation
As for the calculation of reserves, some experts have proposed some effective methods since the middle of last century [7][8][9], the most widely used methods are Havlena-Odeh method [10], Campbell method [11] and Tehrani method [12].

Havlena-Odeh method
Havlena and Odeh recombined the linear form of the material balance equation as follows: If e W can be expressed as: Then: whereU is water influx constant,rm 3 /day; S is the water body function, and its expression is different for different types of water bodies [13]. If the water model is correct, should be a straight line with slope U and intercept N as shown in Figure 1. When using this method to determine the reserves and water influx constant, we need to use the trial and error method, and constantly assume different water body parameters until t F E and t S E meet the linear relationship.

Tehrani method
Tehrani pointed out that equation (18) is sensitive to the error of average pressure, and proposed that the least square regression should be directly applied to equation (1) to determine reserves and water influx constant [12].
The measured and calculated deviations of F is defined as follows: Then the sum of squares of errors is: Reserve N and water influx constant U can be obtained by minimizing SSD . Similar to the Havlena-Odeh method, equation (19) needs to determine the water body function S first.

Campbell method
Campbell proposed a linear form of the material balance equation: When there is no water influx, t F E with respect to F should be a horizontal line with the value of N ; when the amount of water influx increases, the calculated value of N increases continuously. When the line of Generally, Campbel curve is used to judge the intensity of water influx, as shown in Fig.  2:

The solution method used in this paper
There are many unknown parameters involved in the mass balance equation, especially the calculation of water influx. The traditional calculation method needs many times of manual calculation, so the calculation efficiency is low. In order to make the calculation results more accurate and reliable, this paper uses the method of combining nonlinear regression and linear regression to comprehensively determine the reserves and water body parameters. The technical route is shown in Fig. 3.
The process of nonlinear regression is: the square of the difference between the measured value and the calculated value of cumulative oil production is taken as the objective function, and the reserves and water body parameters can be obtained by calculating the minimum value of the objective function. Therefore, the problem is transformed into the following nonlinear regression problems: That is, solving a set of parameters c makes ( ) f c minimum. The Levenberg-Marquardt method is used in the nonlinear regression program [14], and the iterative formula is shown in equations (24) and (25) ( ) Where J is Jacobian matrix, ; W is weight matrix; λ is damping coefficient; D is difference between calculated value and measured value.
For the linear regression method, since the calculation model of Fetkvoic quasi steady water body and Carter-Tracy unsteady water body can not be expressed as the form shown in equation (17), the Havlena-Odeh linear equation and the Tehrani method are no longer applicable. In this paper, equation (1) is written as follows: Then e F W − with respect to t E should be a straight line with slope N and passing through the origin. According to this feature, the correctness of the water parameters obtained by the nonlinear regression method can be verified and the reserves N can be obtained.

Uncertainty analysis
Although reservoir numerical simulation has become a mature technology, there is still a problem that has not been solved well, that is, the uncertainty of the model. The parameters describing reservoir geology and fluid properties, as well as the recorded historical data of production performance, have varying degrees of uncertainty. In this paper, Monte-Carlo simulation method is used to analyze the uncertainty of numerical simulation model, which can recognize and reduce the uncertainty of numerical simulation model. The specific process can be divided into the following steps: (1) The source and degree of uncertainty of basic data are analyzed.
(2) Set the uncertainty error range of basic data, and use random sampling method to generate N groups of basic data with a certain error range, they obey the random error of certain distribution, or have the deterministic error set manually (3) Based on the material balance calculation, Monte-Carlo simulation is carried out. The uncertainty of the simulation prediction results is analyzed by using the statistical principle, and compared with the measured data, the adjustable range of the uncertain parameters in the numerical simulation is determined.

Example analysis 4.1 Theoretical calculation example
Based on the PUNQ-S3 model [15], a group of production data is built for material balance theoretical analysis by numerical simulation. The data is as Table 1. The model information is as follows. The grid of the three dimensional three phase black oil model is 19×28×5, and there are two faults in the west and south side, and water body is connected in the west and north. There are 5 producers and no injector. Production history is 16.5 years including 1 year for test and 3 years for shut and production for 12.5 years. The reserve is calculated about 17.373×10 6 sm 3 by volumetric method. The reservoir and water parameters in the numerical simulation are list in table 2. The PUNQ-S3 model is used because the complexity of this theoretical model is fit for the realistic reservoir. The target of this example is to show the reasonability to get the global sensitivity parameters of numerical simulation using material balance method and verify the calculation result.

Reserve and water body calculation
Calculating the F and Et by the data in table 1, and Campbell curve is drawn, as Figure 4.   From the figure4 we can judge that strong water influx exists in the reservoir. So the water influx volume is needed to consider when the material balance is analyzed.
The reserve and the water body parameter are analyzed by linear regression and the result is as figure 5. From figure 5 the result difference is obvious with water body and without water body. It shows that water influx affects the material balance calculation much. The calculated oil reserve is 17.159×10 6 sm 3 , similar with the volume method result. The water influx constant B and time dimensionless factor c t are calculated by the water body parameter from the regression analysis..
Based on the water parameter from the non-linear regression, the relationship curve about e F W − and t E is drawn as figure 6. The linear relationship of e F W − and t E is confirmed, and the reserve is calculated by the linear least squares analysis and the result is 17.365×10 6 sm 3 , which is fit with the non-linear regression result.  Finally the pressure is calculated by the parameters got from the regression analysis, and the test pressure is compared with the calculated data in figure 7. Two pressure is in consistent law of change, which shows the calculation result is correct. Now we finish the calculation of reserve and water body parameters as the process of figure 3, and next we will confirm the reasonability of the result.

Result confirm
From the above analysis we can see the reserve calculated is fit with the numerical simulation data. The water influx is calculated by water influx model and by numerical simulation in figure 8. The water influx volume calculated by material balance is similar with that got by simulation, which shows the material balance can obtain the input parameters for the numerical simulation. Two water bodies are set in the numerical simulation model. But only one set parameter is got in the material balance calculation. Though the parameters got from the material balance calculation cannot be one to one correspondence relationship with simulation model, the water influx volume is coincidence with the real reservoir. So the water influx volume can be the standard data put in the model and need not adjust often.

Real reservoir example
The real reservoir example is an oil field in the east area of china and developed form 1970. The oil field lies in the middle of the anticlinal structure of the fornix and there are 85 wells totally. The pressure and the water cut test data are list in figure 9 compared with the initially calculated data in the numerical model. From the figure 9 we can see, the difference between the test data and the simulation data. And the pressure calculated increases fast and has great difference with the realistic production data. The last time set in the model is 11748 days and the numerical calculation stop at 6452 days due to the non-convergence. The reason of the difference is analyzed and the injection split is not correct and is big uncertain. The uncertainty analysis method is used to confirm the change of the range of injection volume. And then the injection volume is corrected.
The reservoir has no natural water by the geological research. So the material balance equation can be simplified to the follows.
Then the material balance equation in 1 t time is as follows.
the material balance equation in 2 t time is as follows.
The o q , g q , w q , g i and w i is oil rate, gas rate, water rate and water injection rate, sm 3 /day. Changing the form and getting the water injection rate as follows.
The injection volume can be calculated by equation (31) when the production, formation volume factor, Rs is obtained. It should be noted that equation (31) don't consider the natural water body, and not universal for other model. If the natural water influx is taken into consideration, then the calculated water injection rate w i is the sum of the natural water influx and injection water. Almost all the oil field are put into production with water injection, and natural water influx is little proportion so the equation (31) can be used approximately.
Due to the measurement error the production, formation volume factor, Rs is uncertainty. In this example we are short of the whole pressure testing data. And the uncertainty is added. So we choose Monte-Carlo simulation method to predict the water injection.
First the data uncertainty is analyzed. Because the oil production is the most important data in the whole oil field management. The measurement and statistic is relatively correct. But the water and gas production data is not so accuracy. The measurement in separator has also big error except normal test and human factor. The production of dissolved gas is not measured and the uncertainty is so big. Due to the less of separate measuring, the injection data in different layer is uncertainty too. The error is increased by the measurement condition and the inaccurate correction. The error of pressure data comes from the tools calibration, frequency of measurement and fluid gradient. The PVT data of crude oil and natural gas come from the testing in the laboratory, and then the PVT data has big uncertainty due to the different condition between the laboratory and formation. Usually the oil field lacks the initial data of the earlier time and need to give the reasonable assumption.
Based on the above analysis and combined with the realistic condition, we assume the error of quantity, formation volume factor, Rs is normal distribution, and the standard deviation of every parameter is o q -5%, w q -10%, g q -20%, o B -15%, g B -20%, w B -5%, s R -15%. The simulation is done for 3000 times, and the average of water injection, P10, P50, P90 as figure 10. The formation pressure change is predicted by simulation using the average, P10, P50, P90 injection data. And it is compared with the test data in figure 11. The predicted pressure changes in the range of the test pressure, and the calculated result is get better. Finally we use the calculated average injection data to do the numerical history matching, and the water cut matching curve in figure 12 is getting better and the precision of the simulation is increased by the correction injection data.

Conclusion
(1) This paper presents a method of accelerating the history fitting of reservoir numerical simulation by material balance analysis optimization, and puts forward the specific implementation scheme from the determination of global sensitivity parameters and uncertainty analysis.
(2) The established global sensitivity parameter solving method can obtain more accurate basic reservoir parameters, including reserves and water parameters, and provide better initial conditions for numerical simulation.
(3) The established uncertainty analysis method can analyze the influence of parameter uncertainty and determine the adjustable range of parameters.
(4) Material balance analysis is a very important work before the exquisite threedimensional numerical simulation, which provides a macro overall understanding and reduces uncertainty.