APPLICATION OF MULTILEVEL CONCEPTS FOR UNCERTAINTY QUANTIFICATION IN RESERVOIR SIMULATION by Doaa Mostafa Ali Elsakout Submitted for the degree of Doctor of Philosophy School of Energy, Geoscience, Infrastructure and Society Heriot-Watt University April 2016 The copyright in this thesis is owned by the author. Any quotation from the report or use of any of the information contained in it must acknowledge this report as the source of the quotation or information. Abstract Uncertainty quantification is an important task in reservoir simulation and is an active area of research. The main idea of uncertainty quantification is to compute the distribution of a quantity of interest, for example oil rate. That uncertainty, then feeds into the decision making process. A statistically valid way of quantifying the uncertainty is a Markov Chain Monte Carlo (MCMC) method, such as Random Walk Metropolis (RWM). MCMC is a robust technique for estimating the distribution of the quantity of interest. RWM is can be prohibitively expensive, due to the need to run a huge number of realizations, 45% − 70% of these may be rejected and, even for a simple reservoir model it may take 15 minutes for each realization. Hamiltonian Monte Carlo accelerates the convergence for RWM but may lead to a large increase computational cost because it requires the gradient. In this thesis, we present how to use the multilevel concept to accelerate conver- gence for RWM. The thesis discusses how to apply Multilevel Markov Chain Monte Carlo (MLMCMC) to uncertainty quantification. It proposes two new techniques, one for improving the proxy based on multilevel idea called Multilevel proxy (ML- proxy) and the second one for accelerating the convergence of Hamiltonian Monte Carlo is called Multilevel Hamiltonian Monte Carlo (MLHMC). The idea behind the multilevel concept is a simple telescoping sum: which rep- resents the expensive solution (e.g., estimating the distribution for oil rate on finest grid) in terms of a cheap solution (e.g., estimating the distribution for oil rate on coarse grid) and ‘correction terms’, which are the difference between the high reso- lution solution and a low resolution solution. A small fraction of realizations is then run on the finer grids to compute correction terms. This reduces the computational cost and simulation errors significantly. MLMCMC is a combination between RWM and multilevel concept, it greatly re- duces the computational cost compared to the RWM for uncertainty quantification. It makes Monte Carlo estimation a feasible technique for uncertainty quantification inreservoirsimulationapplications. Inthisthesis, MLMCMChasbeenimplemented on two reservoir models based on real fields in the central Gulf of Mexico and in North Sea. MLproxy is another way for decreasing the computational cost based on con- structing an emulator and then improving it by adding the correction term between the proxy and simulated results. MLHMC is a combination of Multilevel Monte Carlo method with a Hamiltonian MonteCarloalgorithm. ItacceleratesHamiltonianMonteCarlo(HMC)andisfaster than HMC. In the thesis, it has been implemented on a real field called Teal South to assess the uncertainty. Dedicated to the souls of my parents... Acknowledgements I would like to express my deepest gratitude to my supervisors, Prof. Mike Christie and Prof. Gabriel Lord for leaving me alone to study and be independent. ThanksforprovidingfreedomtoexploredifferentideasthroughoutmyPhDjourney. Also, I would like to thank the Uncertainty quantification group for improving my presentation skills and giving feedback. I would like to thank Uncertainty quantifi- cation sponsors for their useful comments on my work especially, Prof. Jonathan Carter. Furthermore, I would like to thank the examiners for this thesis, Prof Pe- ter King and Dr. James Cruise for their useful comments and incredible feedback to improve the thesis. Moreover, thanks to computer support team and the best magician Jack Talbot for fixing software problem issues during my PhD. I would like to thank Ali Danesh Scholarship for funding my PhD study at Heriot-WattUniversity. IwouldliketothanktheAfricanInstituteforMathematical Science for providing funding to conference participation. Thanks to Faculty of Science, Cairo University for giving me a study leave to study PhD. I would like to thank my best friend Samah Alhafian for supporting me during my PhD study. Also, I would like to thank my friends Laila, Maha, Radiha, Razan, Alyaa, and Sohad for their effort to make me enjoying the time here in Edinburgh during my PhD. Moreover, my friends at Uncertainty quantification group Zainab, Alexandra, Junko and Behzad for supporting me during my study. Last but not least, I would first like to thank my mother without her continuous support and encouragement I never would have been able to achieve my goals. I dedicate this PhD for my parents’ souls. i Contents 1 Introduction 1 1.1 Thesis Objectives and Statement . . . . . . . . . . . . . . . . . . . . 4 1.2 Thesis Outline . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5 2 Background Material 7 2.1 Mathematics of a Flow in Porous Media Flow . . . . . . . . . . . . . 7 2.1.1 Darcy’s Law . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8 2.2 Reservoir Simulation . . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.3 Deterministic Versus Stochastic Modelling . . . . . . . . . . . . . . . 10 2.3.1 Fractional Flow in a Stochastic Setting . . . . . . . . . . . . . 11 2.3.2 Uncertainty Sources in Reservoir Simulations . . . . . . . . . . 12 2.4 An Inverse Problem . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 2.4.1 History Matching . . . . . . . . . . . . . . . . . . . . . . . . . 13 2.5 Uncertainty Quantification . . . . . . . . . . . . . . . . . . . . . . . . 15 2.5.1 Why has Probability Theory been Used for Uncertainty quan- tification? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.5.2 Bayesian Inference . . . . . . . . . . . . . . . . . . . . . . . . 16 2.5.2.1 Bayes’ Theorem . . . . . . . . . . . . . . . . . . . . 16 2.6 Objective Function . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.6.0.2 How to Estimate σ . . . . . . . . . . . . . . . . . . 18 tj 3 Stochastic Algorithms–Literature 21 3.1 Optimisation Algorithms . . . . . . . . . . . . . . . . . . . . . . . . . 22 3.1.1 Neighbourhood Algorithm (NA) . . . . . . . . . . . . . . . . . 22 ii 3.1.2 Particle Swarm Optimisation (PSO) . . . . . . . . . . . . . . 23 3.1.2.1 PSO Advantages . . . . . . . . . . . . . . . . . . . . 24 3.2 Neighbourhood Algorithm Bayes (NAB) . . . . . . . . . . . . . . . . 24 3.3 Rejection Sampling (RS) . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.4 Markov Chain Monte Carlo (MCMC) . . . . . . . . . . . . . . . . . . 28 3.4.1 Introduction to Markov Chains . . . . . . . . . . . . . . . . . 28 3.4.2 Metropolis-Hasting Algorithm (MH) . . . . . . . . . . . . . . 31 3.4.2.1 Choosing the Proposal Distribution . . . . . . . . . . 31 3.4.3 Gibbs Sampling . . . . . . . . . . . . . . . . . . . . . . . . . . 33 3.4.4 Advantages and Disadvantage of MCMC Methods . . . . . . . 34 3.4.5 Chain Set-up . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 3.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36 4 Numerical Solution for Conservation Equations–Literature 37 4.1 Finite Difference Schemes . . . . . . . . . . . . . . . . . . . . . . . . 38 4.2 Analytical Solution of the First Order Hyperbolic Equation . . . . . . 39 4.3 Advection Equation . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40 4.3.1 Numerical Solution of the Advection Equation . . . . . . . . . 40 4.3.1.1 Single Point Upstream Weighting (Upwind) Scheme . 41 4.3.1.2 Lax-Wendroff Scheme . . . . . . . . . . . . . . . . . 43 4.4 Buckley-Leverett Equation . . . . . . . . . . . . . . . . . . . . . . . . 44 4.4.1 Derivation of the Buckley-Leverett Equation . . . . . . . . . . 45 4.4.2 Analytical Solution . . . . . . . . . . . . . . . . . . . . . . . . 46 4.4.3 The Numerical Solution of the Buckley-Leverett Equation . . 47 4.5 Pressure Equation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 48 4.5.1 Analytical Solution . . . . . . . . . . . . . . . . . . . . . . . . 48 4.5.2 Similarity Solution . . . . . . . . . . . . . . . . . . . . . . . . 49 4.5.3 The Numerical Solution of the Pressure Equation . . . . . . . 50 4.5.3.1 Stability Condition for Explicit and Implicit Schemes 50 4.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 iii 5 Multilevel Monte Carlo for Porous Media Flow 52 5.1 Monte Carlo Method . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 5.1.1 Monte Carlo Integration (MCI) . . . . . . . . . . . . . . . . . 53 5.2 Two-level Monte Carlo . . . . . . . . . . . . . . . . . . . . . . . . . . 57 5.3 Multilevel Monte Carlo (MLMC) . . . . . . . . . . . . . . . . . . . . 59 5.3.1 The History of MLMC . . . . . . . . . . . . . . . . . . . . . . 61 5.3.2 MLMC Implementation . . . . . . . . . . . . . . . . . . . . . 62 5.3.2.1 MLMC for Stochastic ODEs . . . . . . . . . . . . . . 63 5.3.2.2 MLMC for Stochastic PDEs . . . . . . . . . . . . . . 63 5.3.3 MLMC Algorithm . . . . . . . . . . . . . . . . . . . . . . . . 64 5.4 Applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 5.4.1 Exponential Growth and Decay Equation (Toy example) . . . 67 5.4.2 Advection Equation . . . . . . . . . . . . . . . . . . . . . . . . 71 5.4.3 Buckley-Leverett Equation . . . . . . . . . . . . . . . . . . . . 78 5.4.4 Pressure Equation . . . . . . . . . . . . . . . . . . . . . . . . 81 5.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 86 6 Multilevel Markov Chain Monte Carlo Applied to Uncertainty Quantification 88 6.1 Random Walk Metropolis (RWM) . . . . . . . . . . . . . . . . . . . . 89 6.1.0.1 Effect of Step Size on RWM . . . . . . . . . . . . . 91 6.2 Multilevel Markov Chain Monte Carlo (MLMCMC) . . . . . . . . . . 92 6.2.1 Two-level MCMC . . . . . . . . . . . . . . . . . . . . . . . . . 92 6.3 Output Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95 6.3.1 Autocorrelation . . . . . . . . . . . . . . . . . . . . . . . . . . 95 6.3.2 Effective Samples . . . . . . . . . . . . . . . . . . . . . . . . . 97 6.3.3 Chain Thinning . . . . . . . . . . . . . . . . . . . . . . . . . . 98 6.3.4 Burning-in Period . . . . . . . . . . . . . . . . . . . . . . . . . 98 6.3.5 Convergence Diagnostic . . . . . . . . . . . . . . . . . . . . . 98 6.3.5.1 Gewek Test . . . . . . . . . . . . . . . . . . . . . . . 99 6.3.5.2 Raftery-Lewis Test . . . . . . . . . . . . . . . . . . . 99 iv 6.4 Teal South . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 6.4.1 Sensitivity Analysis . . . . . . . . . . . . . . . . . . . . . . . . 102 6.4.2 Teal South Results . . . . . . . . . . . . . . . . . . . . . . . . 103 6.5 Scapa Field Description . . . . . . . . . . . . . . . . . . . . . . . . . 114 6.5.1 Scapa Results . . . . . . . . . . . . . . . . . . . . . . . . . . . 115 6.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120 7 Multilevel Proxy for Quantifying Uncertainty 122 7.1 Experimental Design . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 7.1.1 Random Sampling . . . . . . . . . . . . . . . . . . . . . . . . 123 7.1.2 Stratified sampling . . . . . . . . . . . . . . . . . . . . . . . . 124 7.1.3 Latin Hypercube Sampling (LHS) . . . . . . . . . . . . . . . . 124 7.1.4 Sobol Sequence . . . . . . . . . . . . . . . . . . . . . . . . . . 126 7.1.5 Radial Basis Function (RBF) . . . . . . . . . . . . . . . . . . 126 7.2 Proxy Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 128 7.3 Error Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 130 7.4 Multilevel Proxy (MLproxy) . . . . . . . . . . . . . . . . . . . . . . . 131 7.4.1 Optimising the Number of Samples . . . . . . . . . . . . . . . 135 7.5 Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136 7.6 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142 8 Multilevel Hamiltonian Monte Carlo for Quantifying Uncertainty in Reservoir Simulation 145 8.1 Hamiltonian Monte Carlo (HMC) Algorithm . . . . . . . . . . . . . . 146 8.1.1 Hamiltonian Dynamics . . . . . . . . . . . . . . . . . . . . . . 147 8.1.2 Leapfrog Method . . . . . . . . . . . . . . . . . . . . . . . . . 148 8.1.3 Relation between the Potential Energy and the Misfit Function152 8.1.4 Advantages and Disadvantages of HMC over MCMC . . . . . 152 8.2 Nadaraya-Watson Kernel Regression . . . . . . . . . . . . . . . . . . 155 8.2.1 Gaussian Kernel . . . . . . . . . . . . . . . . . . . . . . . . . . 156 8.2.1.1 Calculate Gradient of Misfit . . . . . . . . . . . . . 157 v 8.2.2 Polynomial Kernel . . . . . . . . . . . . . . . . . . . . . . . . 157 8.3 Comparison between MLMCMC and HMC . . . . . . . . . . . . . . . 159 8.4 Multilevel Hamiltonian Monte Carlo (MLHMC) . . . . . . . . . . . . 161 8.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 169 9 Conclusion and Future Work 172 9.1 Key Findings . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173 9.2 Future Plan . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 176 A 179 B 182 Bibliography 185 vi List of Tables 5.1 Different choices for the number of sample required to use for solving the SPDE. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64 5.2 Comparison between MCI and MLMC with Lax-Wendroff scheme for solving the advection equation. . . . . . . . . . . . . . . . . . . . . . 77 6.1 Significance of the autocorrelation values . . . . . . . . . . . . . . . . 96 6.2 Parameterisation and prior ranges for Teal South (Hajizadeh et al., 2011) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102 6.3 Parameterisation and prior ranges for Scapa (Farooq, 2011). . . . . . 115 vii