Abstract
　　The ensemble-based history matching technique has been successfully applied to simultaneously estimate multiple petrophysical parameters for hydrocarbon reservoirs. The tuning petrophysical properties include horizontal and vertical permeability， porosity and three-phase relative permeability curves. Four scenarios with different combination of the tuning parameters have been evaluated. The ensemble-based history matching technique is found to be capable of estimating multiple petrophysical parameters by conditioning the reservoir geological models to production history. The uncertainty range of production data generated from the updated models is reduced compared to that of initial models. However， the history-matched models may not always provide good production prediction results， especially when absolute permeability and relative permeability are tuned simultaneously. This further illustrates the non-uniqueness of the history matching solutions. In addition， three-phase relative permeability curves are found to be estimated with good accuracy when absolute permeability fields are known.
　　Key words： Petrophysical parameters； Assisted history matching； Ensemble kalman filter （EnKF）； PUNQ-S3 model
　　Li， H.， & Yang， D.Y. （2012）. Estimation of Multiple Petrophysical Parameters for Hydrocarbon Reservoirs with the Ensemble-Based Technique. Advances in Petroleum Exploration and Development， 4（2）， -0. Available from： http：//www.cscanada.net/index.php/aped/article/view/j.aped.1925543820120402.812
　　DOI： http：//dx.doi.org/10.3968/j.aped.1925543820120402.812
　　INTRODUCTION
　　History matching is a process aims to find a model such that the difference between the performance of the model and the production history of a reservoir is minimized （Tavassoli， 2004）. It is a typical inverse problem where input parameters （e.g.， permeability and porosity） of the flow system are estimated from matching system output （e.g.， production） data. This inverse problem becomes more complicated as the parameters to be estimated normally outnumber the output data. Traditionally， history matching is performed by reservoir engineers based on a trial and error approach. Extensive manpower and experience are needed to conduct the manual history matching on a field-wide case. Recently， computer assisted history matching technique has been utilized to assist in automatically adjusting the parameters and accelerating the history matching process. Assisted history matching is a semi-automatic process where reservoir engineers use their engineering knowledge and experience with assistance of a computer program to perform history matching. With the advent of modern computer technology， high computational capacity allows reservoir engineers to use advanced nonlinear optimization methods to solve the inverse history matching problem. 　　Several optimization algorithms have been proposed in the literature. One of the main categories of optimization algorithm is the gradient-based algorithm， which includes the steepest descent， Newton， quasi-Newton， and conjugate gradient methods. For example， the Broyden-Fletcher-Goldfarb-Shanno （BFGS） method （Yang & Watson， 1991）， limited memory BFGS （LBFGS） method （Eydinov et al.， 2009） and Levenberg-Marquardt algorithm （Reynolds et al.， 2004； Chen et al.， 2008） have been used to estimate relative permeability. These methods require the input of the gradient （first derivative） and/or Hessian matrix （second derivative） of the objective function. Computation of the gradient requires generating sensitivity coefficients， which is the core of the above-mentioned algorithms. Such generated sensitivity coefficients are dependent on the number of model parameters， number of observed data， and efficiency of the linear system solver. Another category of algorithm， which has been used in assisted history matching， is the global optimization algorithm. Global algorithms， such as genetic or evolutionary algorithm （Schulze-Riegert， 2002； Choudhary et al.， 2007； Yang et al.， 2009） and simulated annealing （Datta-Gupta et al.， 1995； Vasco et al.， 1997）， have been used to overcome some of the pitfalls of the gradient-based methods because they do not require derivative computation and use only the objective function value.
　　In the past decade， the ensemble Kalman filter （EnKF） technique has been proved to be an efficient data assimilation method and successfully used in assisted history matching for estimating reservoir petrophysical parameters， such as porosity， absolute and relative permeability， fluid-contact depth. Numerous efforts have been made not only to investigate the characteristics of EnKF method， but also to improve the algorithm itself. Few attempts， however， have been made to study impacts of simultaneously tuning multiple parameters （e.g. absolute permeability and relative permeabilities） on the estimation and prediction results. More and more permanent downhole sensors are being deployed for monitoring and recording pressure， temperature， and/or flow rates. In this case， since the data output frequency is very high， it is important to incorporate these recorded parameters as soon as they are available to keep the reservoir model always up-to-date. Traditional history matching is not suitable for such a purpose because of the heavy computational burden and the high data sampling frequency. 　　In this paper， an ensemble-based history matching technique has been successfully applied to simultaneously estimate multiple petrophysical parameters for a hydrocarbon reservoir when dealing with different tuning scenarios. The well-known UN certainty Quantification （PUNQ）-S3 reservoir model， which is a standard synthetic test case that is based on a real field， was used as the testing model. The tuning petrophysical properties include horizontal and vertical permeability， porosity and three-phase relative permeability curves. Efforts have been made to tune and estimate petrophysical parameters for four different scenarios， and the respectively estimation and prediction results are subsequently presented and discussed. Finally， conclusions are presented.
　　1. ENSEMBLE KALMAN FILTER
　　The ensemble Kalman filter （EnKF） was introduced to overcome some of the problems of the extended Kalman filter. In particular， instead of directly estimating the necessary statistics based on linear assumptions， the EnKF method uses an ensemble of model representations from which all necessary statistics can be directly computed approximately under the Bayesian framework. The EnKF method has been widely applied in the areas of weather forecasting， oceanography and hydrology （Hamill， 2006）. In these applications， only the dynamic variables need to be tuned. It has recently been introduced to the petroleum industry for continuously updating reservoir geological models， within which both static and dynamic variables are simultaneously tuned to assimilate new measurements. The EnKF technique is a Monte-Carlo method where the covariance matrix is updated from a limited number of ensembles， mathematical description of the method can be found in detail elsewhere （Evensen et al.， 2006）.
　　Compared to the gradient-based optimization methods， the EnKF method does not require the gradient of the objective function. In addition， the EnKF method only uses the input and the output of a reservoir simulator so that it can be integrated with any reservoir simulator. Due to these features， EnKF method is easier to implement for assisted history matching compared to the gradient-based optimization methods. N?vdal et al. （2003） used the EnKF method to update the static parameters in near-well reservoir models by tuning the permeability field. Brouwer et al. （2004） used the combination of the EnKF method for continuous model updating with an automated adjoint-based technique to optimize the waterflooding strategy. Results from the previous studies showed that the EnKF method was very efficient and robust. However， the estimation quality of permeability field may deteriorate at a later time of the history matching period， which was often referred to as the “filter divergence”. Gu and Oliver （2005； 2007） used the EnKF method to update the porosity and permeability fields as well as the saturation and pressure fields， and then applied it to match the three-phase production data at wells from the PUNQ-S3 model. Gao et al. （2006） found that the randomized maximum likelihood （RML） method and the EnKF method provide comparable uncertain characterization results for the PUNQ-S3 model. Heidari et al. （2011） updated the porosity and permeability fields of PUNQ-S3 model using different ensemble size and assimilation step and found that a larger ensemble size and a smaller time step result in better estimation results. 　　It was found that accuracy of the traditional EnKF methods for nonlinear and multimodal problems may not be guaranteed （Zafari & Reynolds， 2007； Li & Reynolds， 2009）. Li and Reynolds （2009） presented two iterative EnKF （IEnKF） procedures to solve the strongly nonlinear systems for history matching. The application of this method on a synthetic study showed that in some cases an adjoint solution back to time zero can improve the match and prediction quality； however， the iterations are needed only at time steps when data cannot be matched. Chen et al. （2009） proposed a reparameterization method to handle strong non-Gaussian properties in history matching. In addition， different forms of localization functions （Arroyo-Negrete et al.， 2008； Chang et al.， 2010； Chen & Oliver， 2010a， b） are introduced to integrate with the EnKF method to avoid unnecessary updates far away from observation data. More comprehensive reviews of the application of the EnKF technique in updating reservoir simulation models have been presented by Aanonsen et al. （2009） and Oliver and Chen （2010）， respectively.
　　One of the limitations pertaining to the EnKF algorithm is that significant inconsistency may occur between the updated static and dynamic variables after the EnKF updating step. This is because the dynamic parameters and the static parameters are updated sequentially. To overcome such a limitation， an iterative EnKF algorithm （Wen & Chen， 2006， 2007） is used in this study to ensure that the updated static and dynamic parameters are consistent at each time step. More specifically， at each assimilation time point， a confirmation procedure is performed after the EnKF forecast step. The confirmation procedure is basically to run the simulation again from the last time step to the assimilation time step using the updated static variable to obtain the updated dynamic variable. In this study， the confirmation procedure is performed at every time step， though this may not be necessary for other cases as shown by Wen and Chen （2007）.
　　2. RELATIVE PERMEABILITY REPRESENTATION MODEL
　　In the context of history matching， static vectors， m， of a model mainly include the local grid discrete parameters （e.g.， horizontal and vertical permeability） and global parameters （e.g.， relative permeability data）. A representation model is required to generate relative permeability curves used for reservoir simulation， while it should have limited numbers of controlling parameter， which can be included as part of the static vector. There are two categories of relative permeability representation models： the parametric model and the non-parametric model. The parametric model uses explicit equations to generate relative permeability curves， assuming relative permeability curves fit into the shape of a certain type of parametric model （e.g.， power law model）. Due to its simplicity， the power law model has been widely used to represent relative permeability curves （Lee et al.， 1987； Reynolds et al.， 2004； Li et al.， 2009）. The non-parametric model is more general and flexible as there is no assumption regarding the shape of relative permeability curves. For example， the B-spline model has been utilized in various studies to represent relative permeability （Watson et al.， 1988； Yang and Watson， 1991； Kulkarni & Datta-Gupta， 2000； Okano et al.， 2005； 2006； Eydinov et al.， 2009； Li et al.， 2010）. 　　In this study， the power law model was used to represent the three-phase relative permeability curves. For the oil-water system，
　　where krw and kro are the relative permeabilities of water and oil， respectively； aw and ao are relative permeabilities of water at Sw=1-Sro and oil at Sw=1-Swc， respectively； bw and bo are exponential or shape factors for determining the shape of relative permeability curves； and the dimensionless water saturation， SwD ， is defined as，
　　where krg and krog are relative permeabilities of gas and oil， respectively； ag and aog are gas relative permeability at Sg=1-Sorg-Swc and the relative permeability to oil at Sg=Sgc， respectively； Sg is gas saturation； Sorg is residual oil saturation； and bg and bog are exponential factors.
　　As for three-phase flow in reservoirs， water and gas relative permeabilities can still be calculated by using Equations （1） and （4）； however， the relative permeability to oil， kro， is calculated by using Stone’s Model II method （Stone， 1973） as follows.
　　From Equation （1） to Equation （6）， it can be seen that a total of twelve parameters define relative permeability curves in a three-phase flow system. If the phase relative permeability is defined by normalizing the effective permeability of each phase by the absolute oil permeability at irreducible water saturation （Eydinov et al.， 2007； Chen et al.， 2008）， then ao=aog=1. So relative permeability vector， mkr， which is part of the static vector， m， can be expressed as，
　　These parameters can be classified into two groups， i.e.， the first six parameters are endpoints and the remaining four parameters are shape factors.
　　The relative permeability vector， mkr ， must be obtained to evaluate relative permeability of a three-phase system. It remains a challenging task to accurately determine the endpoints of relative permeability curve by using the assisted history matching （Chen et al.， 2008）. To simplify the model， it is often assumed that endpoints are known and thus the number of parameters in Equation （7） can be further reduced. In this study， the endpoints of the relative permeability curves were assumed to be known， thus only shape factors were estimated using history matching.
　　3. PUNQ-S3 MODEL
　　The performance of the ensemble-based history matching technique was investigated using the well-known Production forecasting with the PUNQ-S3 reservoir model， which has been used by many research groups to test the ability of different methods in terms of history matching and uncertainty qualification （e.g.， Floris et al.， 2001； Gu & Oliver， 2005； Gao et al.， 2006； Liu & Mcvay， 2010； Wang et al.， 2010）. The PUNQ-S3 model has been taken from a reservoir engineering study on a real field operated by Elf Exploitation （Floris et al.， 2001）. It was qualified as a small-size industrial reservoir engineering model. The data available for the PUNQ-S3 model included the porosities and permeabilities at well sites， and the synthetic production history of the first 8 years. 　　The task is to predict the total oil production after 16.5 years including uncertainty quantification. The model contains 19×28×5 uniform grid blocks with an areal dimension of 180×180 m2， among which 1761 blocks are active. The oil， water and gas density at surface conditions is 911.93， 1032.04 and 0.83 kg/m3， respectively. The gas-oil contact and the oil-water contact are located at 2355.0 and 2394.7 m， respectively. As shown in the top structure map （see Figure 1）， the field is bounded to the east and south by a fault， and linked to the north and west to a fairly strong aquifer. A small gas cap is located in the center of the dome-shaped structure. The field initially contains 6 production wells located around the gas-oil contact. Due to the presence of a strong aquifer， there is no injection well in the reservoir. The geometry of the field has been modeled using corner-point geometry.
　　The original PUNQ-S3 reservoir model is generated with a simulation period of 16.5 years. The first year is for extended well testing， which consists of four three-month production periods， each having its own production rate. After the extended well test， all wells are shut-in for three years. During field production， two weeks of each year are needed for each well to do a shut-in test to collect the corresponding pressure data. The wells are operated under production constraint， i.e.， after falling below a limiting bottomhole pressure （BHP）， they will be switched to the preset BHP constraint， which is 12000 kPa for all 6 wells. The “true” total oil recovery after the 16.5 year is 3.92×106 Sm3. The reference production data are generated by using a reservoir simulator.
　　The reference production data types used in history matching process include BHP， oil production rate （OPR）， water-cut （WC） and gas-oil ratio （GOR）. Production data are available at 15 time steps， the number of the available reference production data at each time step is presented in Table 1. Note that the number of reference production data is different in different time steps. The observation data is generated by adding Gaussian noise where zero mean has been added to the reference production data to mimic the errors in such data. The noise level on BHP， GOR， WC and OPR are 100 kPa， 10%， 1% and 0.0001 Sm3/d， respectively. It should be noted that observation data noise level is not correlated to time.
　　The heterogeneity of porosity and permeability fields is modeled using the Gaussian random fields. The properties generated in each layer are performed independently. In each layer， the correlation coefficient between porosity and horizontal permeability is 0.8， correlation between horizontal and vertical permeability is also 0.8. The detailed assumptions and procedures that are used to generate the reference case can be found elsewhere （Floris et al.， 2001）. The reference horizontal and vertical absolute permeability and porosity fields for Layer #1 are shown in Figure 2. Parameters of other layers can be found elsewhere （Li， 2010）. The reference relative permeability （see Figures 3a and b） is generated using the power law model together with the parameters listed in Table 2. 　　4. TESTING SCENARIOS
　　In petroleum industry， scarce information is available at well sites， thus there normally exists a high level of uncertainty in various subsurface petrophysical properties like porosity （φ）， horizontal permeability （kh）， vertical permeability （kv） and relative permeability （kr）. In this study， four different testing scenarios have been set up in order to evaluate performance of estimating different types of reservoir petrophysical properties using the newly developed ensemble-based history matching technique.
　　Scenario #1： Porosity together with horizontal and vertical absolute permeability fields are estimated by history matching the production history， assuming that relative permeability for oil， water， and gas phases are known without uncertainty.
　　Scenario #2： Only relative permeability for oil， water， and gas phases are assumed to be the unknown parameters， while porosity and horizontal and vertical absolute permeability fields are known without uncertainty.
　　Scenario #3： Horizontal and vertical absolute permeability， relative permeability and porosity are all tuned simultaneously during the history matching process.
　　Scenario #4： Relative permeability for oil， water， and gas phases together with porosity are estimated by history matching， while horizontal and vertical absolute permeability fields are assumed to be known without uncertainty.
　　The detailed configurations of four testing scenarios are shown in Table 3. Note that cells with the symbol “+” mean the corresponding petrophysical properties are tuned in the history matching process， whereas the petrophysical properties with “-” marks are assumed to be known without uncertainty and thus will not be tuned. In each scenario， as listed in Table 3， the total number of parameters in each member varies accordingly when different petrophysical parameters are tuned. For example， in Scenario #1， three types of petrophysical parameters are tuned in each grid block， thus a total of 7980 parameters are tuned in each model. It is worthwhile mentioning that， the same reference production data， as listed in Table 1， are used for all scenarios. For Scenarios #1， #3 and #4， 80 ensemble members are used. An ensemble size of this order has been applied to similar applications （Gu & Oliver， 2005； Lorentzen et al.， 2005）. For Scenario #2， a smaller ensemble size of 50 is used because only four relative permeability shape factors are tuned in this scenario. 　　5. RESULTS AND DISCUSSION
　　5.1 Scenario #1
　　For Scenario #1， both initial permeability and porosity fields need to be generated. Relative permeability curves are not tuned in Scenario #1， and thus do not need to be initialized. The initial ensembles were generated， assuming there is no correlation between porosity and permeability. The absolute permeability and porosity fields of the initial ensembles are generated using the Gaussian geostatistical simulation tools embedded in the CMG Builder （Version 2010.10）. These parameter fields are generated using geostatistical parameters that are used to generate the true case （Gao et al.， 2006）. The porosity， horizontal and vertical permeability in each layer is generated by constraining to the hard data that are observed at well sites.
　　The mean value of the porosity， horizontal permeability and vertical permeability of Layer #1 of the initial ensembles are shown in Figure 2. These fields are average values for 80 conditional geostatistical simulations. It can be seen from these figures that， the mean property fields of initial ensemble are much smoother than the reference cases due to the fact that when the ensemble size is large， the mean value of the Gaussian sequential simulation converges to the Kriged solution （Yarus， 1994）. As a result of taking average among all ensemble members， the extreme values of the properties in Layers #1 and #3 are not found in the mean of the initial ensemble. For example， comparing Figure 2a with Figure 2d， it can be seen that the high porosity zone disappeared in the average initial porosity field. To meet the Gaussian assumption of the EnKF method， a natural logarithm transformation is used for the permeability data in history matching process， i.e.， the log values of horizontal and vertical permeability are tuned using the ensemble-based technique in the history matching process. The porosity data is not transformed.
　　The mean value of the updated porosity， horizontal permeability and vertical permeability fields of Layer #1 from Scenario #1 are also shown Figure 2. It can be seen that， the updated values are closer to the true case. For example， even though the updated porosity in Layer #1 （Figure 2g） is not identical to the reference case （see Figure 2a）， the region of high porosity of the updated case is much closer to the reference case than that of the initial case （Figure 2d）. The updated horizontal permeability in Layer #1 （Figure 2h） is higher than that of the initial case （Figure 2e）， and the high/low streak characteristic featured in the true case （Figure 2b） is captured in the updated models. 　　Note that the values in each layer are truncated according to the corresponding upper and lower limits of true values. Because of truncation， the overshooting problem of permeability and porosity is not exhibited in these figures； however， it does occur in the history matching process. This has also been reported by other researchers （Gu and Oliver， 2005； Naevdal et al.， 2005； Dong et al.， 2006； Gao et al.， 2006； Skjervheim et al.， 2007）. The overshooting is caused mainly by the strong non-linearity inherent in multiphase flow simulations. It is necessary to perform further investigation on this issue in the future.
　　Prior to the history matching， the initial ensembles are used in the simulation to model and evaluate reservoir performance of 16.5 years in order to investigate the uncertainty involved in the initial ensemble members. The cumulative oil and water production for the whole field is shown in Figure 4. Data for both the history matching period （first 8 years） and subsequent 8.5 years of prediction are presented. In all these figures， the red lines denote the reference case， while the grey lines represent results from different ensemble members. As shown in Figure 4， field cumulative water production generated from initial models are much higher than that of the reference case. However， field cumulative oil production is well matched. This is due to the fact that OPRs are main constraint for all producers， and the targeted OPRs can normally be met in the production period.
　　The cumulative oil and water production for the whole field generated using the updated reservoir models in Scenario #1 are also presented in Figure 4 for comparison. It can be seen that the uncertainty range of the predicted cumulative oil production （Figure 4c） using updated models is smaller than that of initial ones （Figure 4a）. The cumulative water production predictions generated using the updated models （Figure 4d） are much closer to the reference case than that of the initial models （Figure 4b）. In addition， the reference production line locates near the center of the forecasted performance using the updated model. This indicates that the updated models provide a reasonable smaller uncertainty range of the forecasted cumulative water production.
　　On a single well scale， the BHP， WC and GOR of Well Pro-11 using the initial and updated reservoir models in Scenario #1 are shown in Figure 5. Note that all these values are generated by conducting reservoir simulation from time zero using the updated reservoir models. Compared to the production results generated using the initial ensembles， the updated models provides better history matching and prediction results. The reference case， denoted by red lines， located in the middle of the space that is spread by the production curves generated using the updated models， indicating that the updated models generate an unbiased estimation of the petrophysical properties of the reservoir. 　　5.2 Scenario #2
　　In this scenario， attempts have been made to use the ensemble-based history matching technique to estimate three-phase relative permeability data. Only relative permeability data are tuned， while all other petrophysical data are assumed to be known without uncertainty. Consequently， the initial models are only different from each other in terms of relative permeability curves.
　　The endpoints are assumed to be known， so only four shape factors of relative permeability are tuned. Relative permeability curves are initialized by assigning random numbers to the shape factors. The upper and lower boundaries of shape factors are listed in Table 4. The mean and updated value of the initial relative permeability curves are shown in Figures 3c and d， respectively. Considering that only four shape factors need to be tuned， a small ensemble size of 50 is used in the history matching process. For the oil-water relative permeabilities， the updated relative permeability curve of oil is higher than the reference one， while the water relative permeability curve is slightly lower than the reference curve. The mismatch mainly occurs at high water saturation zones， which have not appeared in the first 8 years of production. As for the oil-gas relative permeabilities， the gas relative permeability curve is estimated with good accuracy， whereas the relative permeability curve of oil is higher than the reference one. The mismatch mainly occurs at high gas saturation zones， whose information has also not been revealed in the production history.
　　To test the uncertainty involved in the 50 initial ensemble members， simulation is first conducted with the initial ensemble members， while the generated cumulative oil and water production are presented in Figures 6a and b， respectively. It is interesting to find out that the predicted oil and water production are close to even distribution. This may be ascribed to the fact that shape factors of relative permeability are initialized using random numbers with even distribution. Furthermore， a huge uncertainty involved in the initial ensemble is clearly demonstrated by the widespread water production forecast shown in Figure 6b.
　　The reservoir models are updated by history matching the production data collected for the first 8 years. After history matching， the updated reservoir models are used to predicate the reservoir performance in the subsequent 8.5 years. The cumulative oil and water production for the whole field generated using the updated reservoir models in Scenario #2 are shown in Figures 6c and d， respectively. Compared to the results generated by the initial ensemble members， the updated models provide a much better prediction of reservoir performance. It can be seen that the predicted water production rate using the updated models are close to the true case， and the uncertainty range is much smaller than that resulting from the initial ensemble members. This indicates that a good estimation of relative permeability is obtained via history matching. 　　Figure 7 depicts the relationship between the estimated shape parameters and time index. As time index increases， shape factor changes gradually towards its corresponding true values， which is shown as a horizontal line in each individual case. In addition， the deviation of each shape factor， shown as error bars in Figure 7， is reduced as time index increases. At end of the history matching process， the estimated shape factors are all close to their corresponding true values. It should be note that this result is obtained when a high uncertainty range is set for the shape factors， and better estimation results would be expected if a lower range of uncertainty was given.
　　5.3 Scenario #3
　　In Scenario #1， demonstration has been made to estimate the local petrophysical properties such as porosity and absolute permeability by using the ensemble-based history matching technique. In Scenario #2， relative permeability， a global petrophysical property， has been estimated with satisfactory results. Subsequently， in Scenario #3， attempts will be made to simultaneously estimate the global and local petrophysical properties in the PUNQ-S3 reservoir. The initial porosity and permeability fields are the same as those used in Scenario #1， while relative permeability curves are initialized using the same boundary conditions as those used in Scenario #2.
　　The cumulative oil and water production generated using 80 initial ensemble members are shown in Figures 8a and b， respectively. Compared to the initial case in Scenario #1 and #2， the uncertainty range of Scenario #3 is the largest， this is because more unknown parameters are involved in Scenario #3. The cumulative oil and water production generated using 80 updated ensemble members are shown in Figures 8c and d， respectively. For the cumulative oil production， the updated models generate much better predicting results compared to that of the initial model. The oil production has been well matched for the first 8 years， and subsequently the predicted OPR is very close to that of the reference case. As for cumulative water production， the water breakthrough time is improved compared to that of initial case （Figure 8b）. However， many of the updated ensemble members predict a higher water production than that of the reference case during the predicting period of 8.5 years. This is attributed to the biased parameter estimation results as discussed subsequently.
　　The mean value of the updated porosity， horizontal permeability and vertical permeability fields of Layer #1 from Scenario #3 are shown in Figures 2j， 2k， and 2l， respectively. It can be seen from these figures that the updated values are much closer to the true ones compared to the initial values shown in Figures 2d， 2e， and 2f， respectively. The high permeability and porosity streaks are well captured in the updated models. The estimated permeability fields， nevertheless， are in general lower than those of the reference case. This is mainly due to the fact that both absolute permeability and relative permeability are controlling the mobility of the fluids. The product of these two properties other than the individual one is used in the reservoir simulation process. For example， in order to achieve a high OPR， either the absolute permeability or relative permeability to oil can be increased. As shown subsequently， the estimated relative permeability to oil is greatly increased， leading to the estimated low permeability values. 　　The initial and updated relative permeability curves are shown in Figures 3e and f， respectively. For the oil/water relative permeabilities， it is clearly shown that relative permeability to oil and water are not matched. The estimated water relative permeability is much lower than the reference and initial case. The reason for the downward tuning of the water relative permeability in the history matching is to delay the early water breakthrough in initial ensemble members， as shown in Figure 8b. As for oil， a much higher estimation of relative permeability to oil is reached at the end of history matching. This is due to the overshooting problem of the EnKF method as well as simultaneously tuning the absolute and relative permeabilities. In the updating process， the estimated bo of many ensemble members fall below the preset lower boundary； consequently， the ultimate mean value of bo is very low， resulting in high estimated oil relative permeability values. For the gas/oil relative permeabilities， the estimated shape factors are not changed much compared to their initial values.
　　From the above discussion on history matching， performance prediction and parameter estimation， it is found that simultaneous adjustment of absolute and relative permeabilities may lead to erroneous estimation of relative permeability； however， the updated models provide good history matching and production prediction results. This is due to the fact that the history matching solution is generally non-unique （Oliver et al.， 2008； Wang et al.， 2010）.
　　5.4 Scenario #4
　　To further investigate the feasible conditions of accurately estimate relative permeability， Scenario #4 is constructed to update porosity and relative permeability simultaneously with known absolute permeability. The shape factors of initial relative permeability are generated using the boundary data listed in Table 4. The endpoints are assumed to be known. The initial porosity fields are the same as that of the initial models in Scenario #1. The permeability fields are assumed to be known. The ensemble size is 80， which is the same as that in Scenarios #1 and #3.
　　To examine the uncertainty of the initial ensemble members， simulation of 16.5 years using 80 initial models has been conducted. The cumulative oil and water production obtained by initial ensemble models are shown in Figures 9a and b， respectively. Similar to Scenarios #2 and #3， when relative permeability curves are tuned in history matching， a much wider range of cumulative water production prediction has been generated. This means that relative permeability data have a profound impact on the production prediction results. After history matching， the cumulative oil and water production generated by the updated ensemble models are shown in Figures 9c and d， respectively. Compared to the initial case， the uncertainty range of reservoir performance of the updated case is much smaller. Both the true oil and water production are preserved in the updated models， though the mean value of the predicted water production is still higher than the true value because the estimated relative permeability to water is high （See Figure 3g）. 　　For relative permeability estimation results， the initial and updated relative permeability curves are presented in Figures 3g and h， respectively. Compared to the estimation results in Scenario #3 （see Figures 3e and f）， the overshooting of bo does not occur， while a good estimation of relative permeability to oil is achieved after history matching. Relative permeability to water does not change much compared to its initial one. For the gas-oil relative permeabilities， the estimated gas relative permeability is improved slightly after the history matching. As for the oil relative permeability， there still exists a mismatch between the updated oil relative permeability curve and the true one. Mismatch mainly appears at high gas saturation zones， due to the fact that the first eight years of production mainly reveals little information about relative permeabilities in high gas saturation ranges. The mean value of updated porosity fields can be found elsewhere （Li， 2010）.
　　5.5 Comparison with PUNQ-S3 Results in the Literature
　　The PUNQ-S3 model is built to test the ability of different methods in terms of history matching and uncertainty qualification. Many research groups have published their results， among which Floris et al. （2001） summarized some of the results as listed in Table 5. To minimize the effect caused by using different simulators， the results are normalized by dividing the cumulative oil production at 16.5 years by the corresponding reference values.
　　Figure 10 is the box plots accounting for the normalized prediction results， which are presented in terms of normalized “P10”， “P50” and “P90” values. The upper boundaries of boxes represent the normalized “P10” value， while the horizontal lines in boxes show the normalized “P50” values and the lower boundaries of boxes show the normalized “P90” values. Five out of eleven scenarios show that the estimated uncertainty ranges do not include the reference case. Since the first three scenarios are generated using homogeneous layer models， a much larger uncertainty range in prediction is obtained. All the remaining scenarios assume that the models are heterogeneous. Scenarios are different in tuning parameters， uncertainty qualification methods， and optimization methods and reservoir simulators that are used. More details of the first eleven scenarios are summarized by Floris et al. （2001）. It can be seen that the uncertainty range generated using different history matching methods are different from each other clearly demonstrating that parameterization， spatial technique， optimization and uncertainty qualification techniques contribute to the final results. 　　The boxes filled in grey are results from four scenarios of this study. Compared to previous results， good production performance predictions have been achieved using the technique developed in this study. The reference case was located in the uncertainty ranges of results generated from all four data assimilation scenarios. It is found from Scenario #1 that， when relative permeability data are known， a good production prediction is achieved and that the normalized “P50” value is close to 1. This means the mean value of the prediction is very close to the reference value. The range between the “P10” and “P90” values reflects a reasonable uncertainty range of the prediction. For Scenario #2， a higher prediction result is observed when relative permeability is tuned， though the uncertainty range of the prediction is comparable to that of Scenario #1. This further demonstrates that relative permeability curves， as global parameters， affect prediction results significantly. However， the uncertainty range of Scenario #1 and #2 is narrow. This indicates that， when the number of parameters is less during history matching， the uncertainty of its predictions tends to be underestimated （Celaya & Wahr， 1996）. The uncertainty range of Scenario #3 is much larger than that of Scenarios #1 and #2 due to the fact that more parameters are tuned in Scenario #3. The predicted “P50” value of the oil production is almost identical to the reference value， even though relative permeability curves are not well matched. Scenario #4 generates a similar uncertainty range as that of Scenario #3. This may be caused by the truncation of the out-of-boundary porosity in the updating process.
　　CONCLUSIONS
　　The ensemble-based history matching technique has been successfully applied to the well-known standard PUNQ-S3 reservoir model for estimating the multiple petrophysical parameters. Four different testing scenarios， with different combination of tuning petrophysical parameters， have been examined. It is found that the ensemble-based technique is capable of estimating petrophysical parameters by conditioning the reservoir geological models to production history. Compared to the results presented in the literature， good production performance predictions have been achieved by using the technique developed in this study. The reference case is located in the uncertainty ranges generated from all four data assimilation scenarios in this study.
　　The selection of parameters will affect the parameter estimation and production prediction results. It is shown from this study that relative permeabilities have profound impacts on the history matching and prediction results. Relative permeability can be estimated with good accuracy when absolute permeability fields are known. Parameters can be estimated with good accuracy if the types of parameter are less， though the uncertainty range of the prediction is underestimated. Also， a good history matching together with a large range of prediction uncertainty is observed when relative permeability is tuned simultaneously with porosity and absolute permeability， though accuracy of the estimated relative permeability is poor. In some cases， history-matched models may not provide good parameter estimation. This further illustrates the non-uniqueness of the history matching. The overshooting of the EnKF method may cause the updated model to be physically unrealistic， and further investigation is necessary. 　　ACKNOWLEDGEMENTS
　　The authors acknowledge a Discovery Grant awarded to D. Yang from the Natural Sciences and Engineering Research Council （NSERC） of Canada.
　　REFERENCES
　　Aanonsen， S. I.， N?vdal， G.， Oliver， D. S.， Reynolds， A. C.， & Vallès， B. （2009）. The Ensemble Kalman Filter in Reservoir Engineering-A Review. SPE Journal， 14（3）， 393-412.
　　Arroyo-Negrete， E.， Devegowda， D.， Datta-Gupta， A.， & Choe， J. （2008）. Streamline-Assisted Ensemble Kalman Filter for Rapid and Continuous Reservoir Model Updating. SPE Reservoir Evaluation & Engineering， 11（6）， 1046-1060.
　　Brouwer， D. R.， & Jansen， J. D. （2004）. Dynamic Optimization of Water Flooding with Smart Wells Using Optimal Control Theory. SPE Journal， 9（4）， 391-402.
　　Celaya， M.， & Wahr， J. （1996）. Aliasing and Noise in Core-surface Flow Inversions. Geophysical Journal International， 126（2）， 447-469.
　　Chang， H.， Chen， Y.， & Zhang， D. （2010）. Data Assimilation of Coupled Fluid Flow and Geomechanics Using the Ensemble Kalman Filter. SPE Journal， 15（2）， 382-394.
　　Chen， S.， Li， G.， Peres， A. M. M.， & Reynolds， A. C. （2008）. A Well Test for In-Situ Determination of Relative Permeability Curves. SPE Reservoir Evaluation & Engineering， 11（1）， 95-107.
　　Chen， Y.， & Oliver， D.S. （2010a）： Cross-Covariances and Localization for EnKF in Multiphase Flow Data Assimilation. Computers & Geosciences， 14（4）， 579-601.
　　Chen， Y.， & Oliver， D. S. （2010b）. Ensemble-based Closed-Loop Optimization Applied to Brugge Field. SPE Reservoir Evaluation & Engineering， 13（1）， 56-71.
　　Chen， Y.， Oliver， D. S.， & Zhang， D. （2009）. Data Assimilation for Nonlinear Problems by Ensemble Kalman Filter with Reparameterization. Journal of Petroleum Science and Engineering， 66（1-2）， 1-14.
　　Choudhary， M. K.， Yoon， S.， & Ludvigsen， B. E. （2007）. Application of Global Optimization Methods for History Matching and Probabilistic Forecasting-Case Studies. Paper SPE 105208 presented at the 15th SPE Middle East Oil & Gas Show and Conference， Kingdom of Bahrain， March 11-14.
　　Datta-Gupta， A.， Vasco， D. W.， Long， J. C. S.， D’Onfro， P. S.， & Rizer， W. D. （1995）. Detailed Characterization of a Fractured Limestone Formation by Use of Stochastic Inverse Approach. SPE Formation Evaluation， 10（3）， 133-140.
　　Dong， Y.， Gu， Y.， & Oliver， D.S. （2006）. Sequential Assimilation of 4D Seismic Data for Reservoir Description Using the Ensemble Kalman Filter. Journal of Petroleum Science and Engineering， 53（3）， 83-99. 　　Evensen， G. （2006）. Data Assimilation： The Ensemble Kalman Filter. Berlin， Germany： Springer.
　　Eydinov， D.， Gao， G.， Li， G.， & Reynolds， A.C. （2009）. Simultaneous Estimation of Relative Permeability and Porosity/Permeability Fields by History Matching Production Data. Journal of Canadian Petroleum Technology， 48（12）， 13-25.
　　Floris， F.J.T.， Bush， M.D.， Cuypers， M.， Roggero， F.， & Syversveen， A.R. （2001）. Methods for Quantifying the Uncertainty of Production Forecasts： A Comparative Study. Petroleum Geoscience， 7（SUPP）， 87-96.
　　Gao， G.， Zafari， M.， & Reynolds， A.C. （2006）. Quantifying Uncertainty for the PUNQ-S3 Problem in a Bayesian Setting with RML and EnKF. SPE Journal， 11（4）， 506-515.
　　Gu， Y. & Oliver， D.S. （2005）. History Matching of the PUNQ-S3 Reservoir Model Using the Ensemble Kalman Filter. SPE Journal， 10（2）， 217-224.
　　Gu， Y.， & Oliver， D.S. （2007）. An Iterative Ensemble Kalman Filter for Multiphase Fluid Flow Data Assimilation. SPE Journal， 12（12）， 438-446.
　　Hamill， T. M. （2006）. Predictability of Weather and Climate （pp. 124-156）. Cambridge， UK： Cambridge Press.
　　Heidari， L.， Veronique， G.， Mickaele， L.R.， & Hans， W. （2011）. History Matching of Reservoir Models by Ensemble Kalman Filtering： the State of the Art and a Sensitivity Study. In Y.Z. Ma & P.R. La Pointe （Eds.）， Uncertainty Analysis and Reservoir Modeling： AAPG Memoir 96 （pp. 249-264）.
　　Kulkarni， K.N.， & Datta-Gupta， A. （2000）. Estimating Relative Permeability from Production Data： A Streamline Approach. SPE Journal， 5（4）， 402-411.
　　Lee， T.Y.， & Seinfeld， J.H. （1987）. Estimation of Absolute and Relative Permeabilities in Petroleum Reservoirs. Inverse Problem， 3（4）， 711-728.
　　Li， G.， & Reynolds， A.C. （2009）. Iterative Ensemble Kalman Filters for Data Assimilation. SPE Journal， 14（3）， 496-505.
　　Li， H. （2010）. Ensemble-Based History Matching and Its Application in Estimating Reservoir Petrophysical Properties （Doctoral dissertation）. University of Regina， Regina， Saskatchewan.
　　Li， H.， Chen， S.， Yang， D.， & Tontiwachwuthikul， P. （2010）. Ensemble-Based Relative Permeability Estimation Using B-Spline Model. Transport in Porous Media， 85（3）， 703-721.
　　Liu， C.， & McVay， D.A. （2010）. Continuous Reservoir-Simulation-Model Updating and Forecasting Improves Uncertainty Quantification. SPE Reservoir Evaluation & Engineering， 13（4）， 626-637.
　　Lorentzen， R.J.， N?vdal， G.， Vallès， B.， Berg， A.M.， & Grimstad， A.A. （2005）. Analysis of the Ensemble Kalman Filter for Estimation of Permeability and Porosity in Reservoir Models. Paper SPE 96375， presented at SPE Annual Technical Conference and Exhibition， Dallas， TX， October 9-12. 　　N?vdal， G.， Johnsen， L.M.， Aanonsen， S.I.， & Vefring， E. （2005）. Reservoir Monitoring and Continuous Model Updating Using the Ensemble Kalman Filter. SPE Journal， 10（1）， 66-74.
　　Okano， H.， Pickup， G.E.， Christie， M.A.， Subbey， S.， & Sambridge， M. （2006）. Quantification of Uncertainty Due to Subgrid Heterogeneity in Reservoir Models. Paper SPE 100223， presented at SPE EUROPEC/EAGE Annual Conference annd Exhibition， Vienna， Austria， June 12-15.
　　Okano， H.， Pickup， G.E.， Christie， M.A.， Subbey， S.， Sambridge， M.， & Monfared， H. （2005）. Quantification of Uncertainty in Relative Permeability for Coarse-Scale Reservoir Simulation. Paper SPE 94140， presented at SPE EUROPEC/EAGE Annual Conference and Exhibition， Madrid， Spain， June 13-16.
　　Oliver， D. & Chen， Y. （2010）. Recent Progress on Reservoir History Matching： A Review. Computers & Geosciences， 15（1）， 185-221.
　　Reynolds， A.C.， Li， R.， & Oliver， D.D. （2004）. Simultaneous Estimation of Absolute and Relative Permeability by Automatic History Matching of Three-Phase Flow Production Data. Journal of Canadian Petroleum Technology， 43（3）， 37-46.
　　Schulze-Riegert， R.W.， Axmann， J.K.， Hasse， O.， Rian， D.T.， & You， Y.L. （2002）. Evolutionary Algorithms Applied to History Matching of Complex Reservoirs. SPE Reservoir Evaluation & Engineering， 5（2）， 163-173.
　　Skjervheim， J.A.， Evensen， G.， Aanonsen， S.I.， Ruud， B.O.， & Johansen， T.A. （2007）. Incorporating 4D Seismic Data in Reservoir Simulation Models Using Ensemble Kalman Filter. SPE Journal， 12（3）， 282-292.
　　Stone， H.L. （1973）. Estimation of Three-phase Relative Permeability and Residual Oil Data. Journal of Canadian Petroleum Technology， 12（4）， 53-61.
　　Tavassoli， Z.， Carter， J.N.， & King， P.R. （2004）. Errors in History Matching. SPE Journal， 9（3）， 352-361.
　　Vasco， D.W.， & Datta-Gupta， A. （1997）. Integrating Multiphase Production History in Stochastic Reservoir Characterization. SPE Formation Evaluation， 12（3）， 149-156.
　　Wang， Y.， Li， G.， & Reynolds， A.C （2010）. Estimation of Depths of Fluid Contacts by History Matching Using Iterative Ensemble-Kalman Smoothers. SPE Journal， 15（2）， 509-525.
　　Watson， A.T.， Richmond， P.C.， Kerig， P.D.， & Tao， T.M. （1998）. A Regression-Based Method for Estimating Relative Permeabilities from Displacement Experiments. SPE Reservoir Engineering， 3（3）， 953-958.
　　Wen， X.H.， & Chen， W.H. （2006）. Real-Time Reservoir Model Updating Using Ensemble Kalman Filter with Confirming Option. SPE Journal， 11（4）， 431-442.
　　Wen， X.H.， & Chen， W.H. （2007）. Some Practical Issues on Real-Time Reservoir Model Updating Using Ensemble Kalman Filter. SPE Journal， 12（2）， 156-166.
　　Yang， C.， Card， C.， & Nghiem， L. （2009）. Economic Optimization and Uncertainty Assessment of Commercial SAGD Operations. Journal of Canadian Petroleum Technology， 48（9）， 33-40.
　　Yang， P. H.， & Watson， A. T. （1991）. A Bayesian Methodology for Estimating Relative Permeability Curves. SPE Reservoir Engineering， 6（2）， 259-265.
　　Yarus， J. M.， & Chambers， R. L. （1994）. Stochastic Modeling and Geostatistics： Principles， Methods， and Case Studies. Tulsa， OK： AAPG.
　　Zafari， M.， & Reynolds， A. C. （2007）. Assessing the Uncertainty in Reservoir Description and Performance Predictions with the Ensemble Kalman Filter. SPE Journal， 12（3）， 382-391.