- Home
- Conferences
- Conference Proceedings
- Conferences
ECMOR XII - 12th European Conference on the Mathematics of Oil Recovery
- Conference date: 06 Sep 2010 - 09 Sep 2010
- Location: Oxford, UK
- ISBN: 978-90-73781-89-4
- Published: 06 September 2010
1 - 50 of 117 results
-
-
Upscaling of Fractured Oil Reservoirs Using Homogenization Including Non-equilibrium Capillary Pressure
Authors H. Salimi and J. BruiningRecovery in incompletely water-wet fractured reservoirs can be extremely low. In the laboratory, these systems are often mistaken for oil-wet reservoirs, because imbibition only starts after removal of the oil layer, which originally covers the grains. The long time required to remove the oil film will be referred to as delay time. There are two theories that describe the delay time necessary for removal of an oil film, leading to a capillary pressure that depends on time. One of the theories is developed by Barenblatt et al. and it modifies both the capillary pressure and the relative permeabilities. The other theory is developed by Hassanziadeh et al. and it only deals with the non-equilibrium effect for the capillary pressure. No attempt has yet been made to model non-equilibrium effects in fractured reservoirs for a field-scale problem and this is an innovative aspect of this paper. To examine whether the non-equilibrium effect has any effect on larger-scale problems, we apply homogenization to derive an upscaled model for fractured reservoirs in which the non-equilibrium effects are included. We formulate a fully implicit three-dimensional upscaled numerical model. Furthermore, we develop a computationally efficient numerical approach to solve the upscaled model. We use the simulation to determine the range of delay times for which discernable effects occur in terms of oil recovery. It is shown that at low Peclet numbers, i.e., when the residence time of the fluids in the fracture is long with respect to the imbibition time, incorporation of delay times of the order of few months have no significant effect on the oil recovery. However, when the Peclet number is large, the delay times reduce the rate of oil recovery. We will discuss for which values of the delay time (Barenblatt) and capillary-damping coefficient (Hassanizadeh), equivalent results are obtained. This is the first time that such a comparison is made for a field scale project and it shows that both approached show the importance of taking to account delay effects in the capillary pressure behavior.
-
-
-
Incorporation of Global Effects in Two-phase Upscaling for Modeling Flow and Transport with Full-tensor Anisotropy
More LessUpscaling is often applied to coarsen highly-detailed geological descriptions. The upscaling of multiphase flow functions is challenging due to their strong dependency on global flow effects. In this work, we present two approaches to generate upscaled two-phase flow functions with global flow effects incorporated. In global two-phase upscaling, the upscaled two-phase functions are directly computed from global fine-scale two-phase solutions. In a new adaptive local-global two-phase upscaling approach, local boundary conditions (for both pressure and saturation) are determined from time-dependent global coarse-scale solutions. This avoids solving global fine-scale two-phase flow. In both approaches, the upscaled two-phase flow functions are adapted to a specific flow scenario, thus providing more accurate coarse models than the use of generic flows. The methods are applied to heterogeneous permeability fields with full-tensor anisotropy. We demonstrate the impact of full-tensor effects on the global flow dependency in two-phase upscaling. It is also shown that the full-tensor effects in two-phase flow are appropriately accounted for by the upscaled two-phase functions with the global flow effects incorporated. The use of those upscaled functions in conjunction with a two-point flux approximation can effectively capture the two-phase crossflow due to the full-tensor anisotropy on the coarse scale.
-
-
-
Meshless Upscaling Method and Its Application to a Fluid Flow in Porous Media
By A. LukyanovIt is known that the two-point flux approximation, a numerical scheme used in reservoir simulators, has O(1) error when grids are not K-orthogonal. The multi-point flux approximations have found significant interest in the research community. However, non-physical oscillations can appear in the developed multi-point flux approximations when the anisotropy is really strong. In this paper, the meshless multi-point flux approximation (MMPFA) for general fluid flow in porous media is proposed. The MMPFA is based on a gradient approximation commonly used in meshless method and can be extended to include higher-order terms in the appropriate Taylor series. The MMPFA is combined with the mixed corrections which ensure linear completeness. The mixed correction utilizes Shepard Functions in combination with a correction to derivative approximations. Incompleteness of the kernel support combined with the lack of consistency of the kernel interpolation in conventional meshless method results in fuzzy boundaries. In corrected meshless method, the domain boundaries and field variables at the boundaries are approximated with the default accuracy of the method. The resulting normalized and corrected MMPFA scheme not only ensures first order consistency O(h) but also alleviates the particle deficiency (kernel support incompleteness) problem.
-
-
-
The Self-flattening Nature of Trailing Shocks in Augmented Waterflooding – Segregation-in-flow Reestablish Self-sharpness
Authors A.M. AlSofi and M.J. BluntOne of the challenges for reservoir simulation is numerical dispersion. For waterflooding applications the effect is controlled due to the self-sharpening nature of a Buckley-Leverett shock. In this paper, we show that for augmented waterflooding – due to the coupling of compositional dispersion with fractional flow – the trailing shock is no longer self-sharpening. Thus, the simulation of such processes suffers from even higher numerical dispersion effects compared to pure waterflooding. Rather than implementing a higher-order discretization method, we propose a simple scheme based on segregation within a grid block. Compared to current mixing schemes, it differs in that segregation not only affects fluid properties but the transport, too. The scheme is shown to reestablish self-sharpness across the trailing shock. Simulations are performed that demonstrate the effectiveness of the scheme. The simulations also illustrate the possible effects of numerical dispersion on the predictions of augmented waterflooding. Finally, we discuss the extension of this technique to compositional simulation through a coupled limited-flash/segregation-in-flow assumption. Preliminary results demonstrate the potential of the approach as a heuristic method to control numerical dispersion for the simulation of miscible and near-miscible gas injection processes.
-
-
-
A Multiscale Framework for Simulating Multi-phase Flow on Flexible Grids
Authors A. Sandvin, J.O. Skogestad, E. Keilegavlen and J.M. NordbottenThe geological structure of porous rocks include irregular geometries and heterogeneities on a multiple of scales. Reservoir simulations (e.g. oil recovery, CO2 storage, ground water flow) often involve large spatial scales, where local fine-scale heterogeneities might have an important impact on the global flow. Discretisations of the governing equations of multi-phase flow in general render large systems of non-linear equations. The multiscale nature of the geology is inherited by the discrete mathematical problem, which can be ill-conditioned and thereby hard to solve. In this work, we apply domain decomposition (DD) as a solution strategy for the equations. The splitting of the global system into local subproblems can be done either prior to, or after linearisation of the non-linear system. Comparison between DD as a non-linear and linear solver strategy indicates that splitting into subproblems should be considered for use as preconditioners, also for non-linear systems of equations. When applied after linearisation, DD is best suited as a preconditioner in an iterative solution procedure. In this work, we consider a two-level DD framework for linear systems. The framework is flexible, and can be used both as an upscaling procedure, and as a preconditioner.
-
-
-
On the Stokes-Brinkman Equations for Modeling Flow in Carbonate Reservoirs
Authors I. Ligaarden, M. Krotkiewski, K.A. Lie, M. Pal and D. SchmidCavities and fractures can significantly affect the flow paths of carbonate reservoirs and should be accurately accounted for during flow simulation. Herein, our goal is to compute the effective permeability of rock samples based on high-resolution 3D CT-scans containing millions of voxels. Hence, we need a flow model that properly accounts for the effects of Darcy flow in the porous material and Stokes flow in the void volumes on relevant scales. The presence of different length scales and large contrasts in the petrophysical parameters lead to highly ill-conditioned linear systems that make such a flow model very difficult to solve, even on large-scale parallel computers. To identify simplifications that render the problem computationally tractable, we analyze the relative importance of the Stokes and Darcy terms for a wide variety of parameter ranges on an idealized 2D model. We find that a system with a through-going free flow region surrounded by a low permeable matrix can be accurately modeled by ignoring the Darcy matrix and simulating only the Stokes flow. Using this insight, we are able to compute the effective permeability of a specific model from a CT-scan that contains more than eight million voxels.
-
-
-
Effective Local and Extended Local Single-phase Upscaling
Authors S. Du, S.A. Hosseini, Y. Zhou and M.J. KingEffective strategies for upscaling of flow in porous media rely upon (1) the characterization of flow in the reservoir (intercell transmissibility), (2) flow near wells (effective permeability) and (3) an a priori analysis of error. Transmissibility upscaling represents approximately constant gradient Darcy flow and, compared to permeability upscaling, maximizes the spatial resolution of intercell fluid flow available within a simulator. Near well permeability upscaling replies upon the line source approximation and the calculation of well productivity for radial flow near wells. It provides a characterization of the local reservoir quality. It can be performed without knowledge of the location or rates of physical wells. The error analysis has three components. It relies upon a connectivity analysis, an analysis of variance of velocity or slowness, and an estimate of off-diagonal flow terms in the permeability tensor. All three of these error components can be estimated using local flow information. This a priori error analysis allows a practicing engineer or a spatially adaptive algorithm to adjust the resolution of the upscaled model in order to minimize the most egregious upscaling errors. All three of these calculations can be performed using local or extended-local boundary conditions. This is in contrast to global upscaling approaches which have been previously used in the literature but which rely upon knowledge of a global flow field. Although global upscaling approaches can be quite accurate, their utility is limited as specific well locations and flow rates are not available at the point in the life of a project when the geologic model is typically upscaled in preparation for reservoir flow simulation. The utility of this combined upscaling and error analysis approach is demonstrated using the SPE10 test model, which is simple structurally but which has a fairly complex and heterogeneous permeability field.
-
-
-
Upscaling of Two-phase Immiscible Flows in Communicating Stratified Reservoirs
Authors X.Z. Zhang, A.S. Shapiro and E.H.S. StenbyA new analytical upscaling method is described. It is applicable for stratified reservoirs or with vertically distributed properties. The inter-layer communication is assumed to be perfect, corresponding to viscous or gravity dominant regimes. We apply asymptotic analysis to 2D flow equations of two incompressible immiscible fluids and reduce 2D Buckley-Leverrett system to 1D flow problem. Flow velocities in the layers depend on saturations of other layers, which reflects interaction between layers; the exchange terms between the layers may be expressed explicitly. The resulting quasilinear hyperbolic equations allow a self-similar solution, which makes it possible to express the pseudo-fractional flow function as a function of average saturation. The system is solved by the finite difference method. Both cases of discontinuous and continuous (lognormal) distributed layer permeabilities are studied. As comparison, the complete 2D waterflooding problems with very good vertical communication are solved in COMSOL. Generally, saturation profiles of the pseudo 1D model are only slightly different from the 2D simulation result by Comsol. The difference between recovery curves is marginal. This proves validity of the proposed method. Our method is much more advantageous over the classical Hearn-Kurbanov procedure. It is better at predicting displacement fronts and oil recovery.
-
-
-
Analytical Solutions for Co- and Countercurrent Imbibition of Sorbing – Dispersive Solutes in Immiscible Two-phase Flow
Authors K.S. Schmid, S. Geiger and K.S. SorbieWe derive a set of analytical solutions for the transport of adsorbing solutes in an immiscible, incompressible two-phase system. Our analytical solutions are new in two ways: First, we fully account for the effects of both capillary and viscous forces on the transport for arbitrary capillary-hydraulic properties. Second, we fully take hydrodynamic dispersion for the variable two-phase flow field and adsorption effects of the solutes into account. All previously obtained results for component transport in immiscible two-phase systems account for changes in the flow field due to the components’ presence but they ignore dispersive effects. This is surprising given the enormous practical importance of diffusive effects both in hydrological settings and for situations that arise in the context of oil production for chemical flooding, polymer flooding and the transport of wettability altering agents. In both cases, one often faces the situation of fractured media or thin-layer structures where the diffussive transport processes of components dominate over viscous forces. We consider a situation where the components do not affect the flow field and focus on dispersive effects. For the purely advective transport we combine a known exact solution for the description of flow with the method of characteristics for the advective transport equations to obtain solutions that describe both co- and countercurrent flow and advective transport in one dimension. We show that for both cases the solute front can be located graphically by a modified Welge tangent and that the mathematically obtained solutions correspond to the physical notion that the solute concentrations are functions of saturation only. For the advection-dispersion case, we derive approximate analytical solutions by the method of singular perturbation expansion. We give some illustrative examples and compare the analytical solutions with numerical results.
-
-
-
Classic and Hybrid MCMC Methods to Approximate Pore Size Distributions of Carbonate Reservoirs Using Pore-network Models
Authors J.E Juri, M.I.J. van Dijke, K.S. Sorbie and S.K. MasalmehCarbonate reservoirs contain a large fraction of the remaining oil reserves, but their highly heterogeneous structures, at the pore scale, make them difficult to characterize and this is one of the main reasons for their low recovery factors. Pore space imaging and reconstruction techniques are often used to characterize the pore space topology and geometry. However, the wide range of pore sizes encountered in carbonates, including submicron microporosity, usually renders these techniques unsuitable. A classic technique to approximate the pore size distribution (Ritter and Drake) is based on inversion of the mercury intrusion capillary pressure (MICP) curve. However, this technique is inaccurate due to the lack of accessibility of the large pores before the percolation threshold is reached. In this paper, we use a simple lattice network, characterized by an average co-ordination number, a pore volume-vs-radius power law and a pore (radius) size distribution (PSD) with an arbitrarily large range of radii, to match the MICP curve. Efficient optimization methods have been implemented to estimate these principal parameters, as outlined below. The overall aim of this work is to predict rock flow functions for different floods or saturations paths that are difficult and expensive to obtain from laboratory experiments, based on the PSD estimations. In the present approach the PSD has a given number of non-constant radius classes (bins), for which the probabilities need to be determined. As an initial guess, we take the Ritter and Drake PSD, based on which the non-fixed bin widths are determined, using principal theorems derived in information theory and data compression. Then, the pore-network parameters including the PSD class probabilities are optimized to match the experimental MICP, using two methods that are based on the link between the Boltzmann - Maxwell equation describing the molecular state probabilities in statistical physics and the Bayes method for posterior distribution. The first method is the classic Metropolis Markov Chain Monte Carlo method and the second, hybrid, method adds Hamiltonian Dynamics to avoid the random nature of the classic Metropolis algorithm, thus enhancing the convergence of the classical approach. The methods have successfully been applied to two synthetic cases and three carbonate samples from a Middle East reservoir presenting three different carbonate facies. The hybrid method improves the Markov Chain Monte Carlo efficiency in terms of the accepted system parameters states by a factor of 4 to 5.
-
-
-
Pore-scale Modelling of Chemically Induced Effects on Two-phase Flow
Authors Y. Zaretskiy, S. Geiger, K. Sorbie and M. FoersterWe present a finite element - finite volume simulation method for modelling fluid flow and solute transport accompanied by chemical reactions in experimentally obtained 3D pore geometries. We solve the stationary Stokes equation on the computational domain with the FE method using the same set of nodes and the same order of basis functions for both velocity and pressure. The resulting linear system is solved by employing the algebraic multigrid library SAMG. To simulate large 3D samples we partition them into subdomains and treat each separately on a different computing node. This approach allows us to use meshes with millions of elements as input geometries without facing limitations in computer resources. We apply this method in a proof-of-concept study of a digitized Fontainebleau sandstone sample. We use the calculated velocity profile with the finite volume procedure to simulate pore-scale solute transport and diffusion. This allows us to demonstrate the correct emerging behaviour of sample’s hydrodynamic dispersivity. Finally, we model the transport of an adsorbing solute and the surface coverage dynamics is demonstrated. This information can be used to estimate the local change of a sample wettability state and the ensuing changes of the two-phase flow characteristics.
-
-
-
Determination of Connectivity in Vuggy Carbonate Rock Using Image Segmentation Techniques
Authors T.P.G. Gurholt, B.F.V. Vik, I.A. Aavatsmark and S.I.A. AanonsenIt is commonly understood that carbonate reservoirs hold much of the worlds' oil and gas reserves. In this study we examine a heterogeneous vuggy pore class. The vuggy pore space is imaged using x-ray computed tomography (CT). Micro-CT ($mu$CT) is used to better image parts of the interparticle (matrix) porosity. Image segmentation is applied to obtain a three-dimensional (3D) map of the vuggy pore space, which is used to determine its connectivity. We evaluate the performance of selected image segmentation techniques applied for segmenting CT-scans of vuggy carbonate rocks at varying resolutions. Computations are performed for both two-dimensional (2D) and three-dimensional (3D) image segmentation. The computations show that the methods perform well on high resolution $mu$CT scans. For a CT scan with low and non-uniform resolution they wrongly classify a large portion of the voxels. This makes it challenging to determine the connectivity. The computed porosity is, due to image resolution, underestimated when compared to measured porosity from laboratory experiments.
-
-
-
Roles of Transient and Local Equilibrium Foam Behavior in Porous Media – Traveling Wave
Authors E. Ashoori, D. Marchesin and W.R. RossenIn foam EOR, complex dynamics of bubble creation and destruction controls foam properties. We assume that local equilibrium applies throughout a foam displacement on the field scale, with the exception of an entrance region and at shock fronts, where saturations and bubble size change abruptly. We find a range of conditions in which the local-equilibrium condition applies even within the shock front. In a waterflood the width of a shock transition zone is determined by capillary-pressure gradients. For foam, this equation is joined by one for evolving foam texture. If there is no gas ahead of the foam, we prove that foam texture is everywhere at local equilibrium within the shock, regardless of the foam model. If there is gas initially in the formation, slow foam generation and coalescence processes can narrow the shock from that assuming local equilibrium. In other cases, the dynamics of the traveling wave leads to oscillations near the shock; these are not numerical artifacts, but reflections of the models. Multiple steady states seen in experiment for some injection rates can be predicted by certain foam models. The approach of solving for the traveling wave can give rule out some of these states.
-
-
-
Extended K-value Method for Multi-contact Miscible Displacements
Authors G. Rannou, D.V. Voskov and H.A. TchelepiAccurate predictions of gas-injection processes usually require a compositional formulation based on an Equation of State (EoS). Because the thermodynamic behavior of multi-component multiphase systems is highly nonlinear and coupling to the flow equations is quite complex, EoS-based simulations can become computationally prohibitive. For immiscible displacements an efficient approximation of the general compositional problem can be used if the phase-equilibrium ratios (K-values) are assumed to be weak functions of composition; however, this assumption is questionable for near-miscible displacements. The standard K-value approach suffers from significant difficulties, both in terms of robustness and accuracy, in the critical region. Here, we describe an extended K-values method that takes advantage of the solution-path invariance in the compositional space with respect to the hydrodynamic properties. Specifically, we propose an additional degree-of-freedom, which captures the composition dependence of the phase behavior, for use in the tabulation and interpolation of the K-values. An important aspect of the method is the use of the so-called Minimal Critical Pressure criteria (MCP), which indicates when a given composition becomes super-critical. We compare results obtained with EoS- and K-values based simulations for several (isothermal) compositional problems, and we demonstrate the efficiency and accuracy of the proposed method.
-
-
-
A New Variable-set for Nonlinear Flow Simulation Based on Compositional Space Parameterization
By D.V. VoskovWe consider numerical modeling of compositional two-phase flow in porous media, and we propose a nonlinear formulation that employs a variable-set based on compositional space parameterization. In the formulation, the phase-fraction and saturation change 'continuously' in the immiscible region of the compositional space. Inside the two-phase region, these variables are identical to the saturation and phase-fraction of the standard approach. In the single-phase regions, however, these variables can become negative, or larger that unity. We demonstrate that when this variable set is used, the EOS computations are resolved completely within the general Newton loop. That is, no separate phase-stability or flash computations are necessary. The number of general Newton iterations grows only slightly, and overall savings lead to more efficient simulations. We discuss using this variable set, which can be thought of as an extension of the natural variable, in two ways. The first scheme honors the nonlinear dependence of the overall density on phase-fractions and saturation, and the second employs a linearized relation for the overall density. Both schemes are compared with the standard natural-variables formulation using several challenging compositional problems. We also describe how the proposed approach can be used for modeling multi-contact miscible displacement processes.
-
-
-
Continuation of the Tie-simplex Space for Three-phase Thermal-compositional Systems
Authors A. Iranshahr, D. Voskov and H. TchelepiModeling the phase behavior associated with compositional flow simulation of systems that form more than two phases (e.g., steam, or CO2 injection) is a challenging problem. In addition to the coupling of thermodynamics with the nonlinear equations of flow and transport, accurate phase-state identification for mixtures that form three (or more) phases, as function of pressure, temperature and composition, is the subject of active research. We describe a general negative-flash method for multi-component thermal systems that can form three, or more, fluid phases, including a convergence proof. This extended negative-flash approach integrates quite nicely with flow simulation based on adaptive tie-simplex compositional space parameterization. We prove that tie-simplexes change continuously with pressure, temperature, or along a continuous trajectory in compositional space. Continuation of the tie-simplex space provides theoretical justification for the tie-simplex based flow modeling, in which interpolation in pressure and temperature for a limited number of tie-simplexes is performed during a simulation run. We also show how a tie-simplex associated with a given composition can degenerate as a function of pressure and temperature. We focus on the complex behaviors of the tie-triangles and tie-lines associated with three-phase, steam-injection problems. The algorithms that capture the degeneration of the tie-triangles into tie-lines are described in detail. Interpolation in the parameterized compositional space is used to identify the phase-state and proves to be as reliable as a 'conventional' phase-stability test for three-phase mixtures. This extended negative-flash algorithm has been integrated into our tie-simplex based flow simulation framework. We demonstrate the effectiveness of the framework using several compositional steam-injection problems with complex behaviors.
-
-
-
New Models for Heater Wells in Reservoir Simulation
Authors G. Aouizerate, L.J. Durlofsky and P. SamierDownhole electrical heating can be used to achieve the high temperatures required for in-situ upgrading of oil shale. Heater-well models are needed if this process is to be simulated accurately. The traditional Peaceman approach used for fluid injection and production wells may not be applicable because it does not capture transient effects, which can be important in downhole heater models. Here we develop two new models for representing heater wells in reservoir simulators. The first model is applicable for homogeneous systems with properties that are not temperature dependent. For such cases we develop a semi-analytical procedure based on Green’s functions to construct time-dependent heater-well indices and heater-block thermal transmissibilities. For the general case, which can include both fine-scale heterogeneity and nonlinearity due to the temperature dependence of rock properties, we present a numerical procedure for constructing the heater-well model. This technique is essentially a near-well upscaling method and requires a local fine-scale solution in the well region. Local boundary conditions are determined using a local-global treatment. The accuracy of the new heater-well models is demonstrated through comparison to reference solutions for two example problems.
-
-
-
Towards Temperature Transient Analysis in Horizontal Wells
Authors K.M. Muradov and D.R. DaviesTemperature modeling and analysis of horizontal wells has become increasingly important due to the extensive deployment of downhole temperature monitoring devices. Such sensors allow continuous monitoring of zonal production with a sufficiently high resolution. Published descriptions of the necessary analysis tools have concentrated on well temperature models and convective heat transfer in porous media. They rarely describe transient, reservoir temperature changes during conventional production. Here, fluid expansion effects are superimposed on heat conduction to the surrounding layers coupled to a transient inflow regime. This paper will discuss why these already available solutions to related problems cannot be employed directly. It will propose a comprehensive model of the reservoir temperature around a horizontal, liquid production well. The above task is treated as a multidimensional, heat transfer problem. Both an asymptotic, early time, solution in real space and an accurate, integral solution can be derived. The asymptotic, analytical solution has been tested against a scientific, simulation softwares. The applicability and practical importance of the solution will be discussed. It forms the basis of a Temperature Transient Analysis tool for horizontal, oil producing wells where both the flow properties of the reservoir fluid and the inflow distribution are to be evaluated.
-
-
-
Analysis of Coupled and Fully Integrated Models for Low-frequency Electrical Heating Assisted Heavy Oil Recovery
Authors J.A. Torres, I.I. Bogdanov, V. Dabir and A.M. KampWe have studied low frequency heating for heavy oil recovery accompanied by salt water recirculation around electrode wells. The multiphysics nature of this problem naturally calls for a code coupling solution. The main objective of this work is to develop an efficient methodology for such coupling. We address the coupling paradigm which includes questions about the degree of coupling (with respect to synthesis), grid and solver type used in the simulation of each separate problem, problem stability and accuracy, interpolation strategy, parallelization, etc. The implemented solution was validated and applied. We coupled a finite volume (FV) reservoir simulator (Stars) to a general purpose finite element (FE) simulator (COMSOL Multiphysics) used to compute an approximation to Maxwell’s equations. These two simulators are coupled by an in-house coupler written in Matlab, being its critical step the bi-direction mapping of interpolation/integration of data between the FE and FV mesh. Selection of solver parameters has been confirmed to be critical in order to limit the effect of numerical dispersion. Coupling solutions, such as the one described in this paper, allow developing simulation tools that reuse as much as possible specialized and existing tools. Their application to practical problems proves to be useful.
-
-
-
Low-temperature In-situ Combustion of Light Oil
Authors A.A. Mailybaev, D. Marchesin and J. BruiningWe consider flows possessing a combustion front when a gaseous oxidizer (air) is injected into the porous medium, a rock cylinder thermally insulated preventing lateral heat losses and filled with light or medium viscosity oil. When oxygen reacts with hydrocarbons at low temperatures, a series of reactions occur that will convert a part of hydrocarbons to oxygenated hydrocarbons and gaseous product (water, carbon dioxide etc.). The oxygenated hydrocarbons are compounds like ketones, alcohols, aldehydes and acids. This process is termed low temperature oxydation (LTO). Indeed, LTO only involves some 25 % of the possible sites that can react with oxygen. Therefore also the reaction heat per volume of fuel can at most be 25 % of full hydrocarbon combustion and consequently the temperature in the LTO reaction zone is very low. Upstream of the LTO reaction zone evaporization occurs. We formulate conservation laws for liquid oil, gaseous oil, oxygen and inert gas components (combustion products and nitrogen) that includes the reaction terms. Moreover we give the energy conservation equation. We give an analytical solution to the equations. It turns out that the solution consists, from upstream to downstream, of a thermal wave, an LTO wave that combines oxidation and evaporation and a Buckley-Leverett saturation wave. It is shown that in the solution a major role is played by the existence of a resonance line at which the derivatives of the oxygen and evaporated oil flux versus the oil saturation vanish. It means that the derivative of the oil flux versus saturation is positive upstream of the resonance line and negative downstream of the resonance line. The complete solution is described for typical parameters of LTO oil combustion. Computations show that only a small part of the oil vaporizes or reacts. The oil velocity is close to the LTO speed. Thus the LTO wave represents a mechanism of almost complete oil displacement by means of the temperature increase in the LTO wave, which leads to a decrease of oil viscosity and increase of the gas flux in the wave.
-
-
-
Analysis of Grid Resolution for Simulations of CO2 Storage in Deep Saline Aquifers
Authors G.E. Pickup, M. Kiatsakulphan and J.R. MillsThe simulation of CO2 injection into deep saline aquifers provides challenges to reservoir modellers. Estimates of storage and trapping of CO2 must be made for large-scale aquifers for planning purposes, and these simulations are generally undertaken on coarse models with grid cells of several 100m in length. However, some of the physical processes which arise during CO2 storage can only be represented on finer grids. In this work, we used a variety of 2D models to study the effect of grid resolution on a range of processes which take place when CO2 is injected into a saline aquifer. A reasonably fine grid is required to model the buoyant rise of the CO2 plume and the migration under the caprock. Cell sizes of several meters in the horizontal and approximately 1m in the vertical are adequate. However, an even finer grid is required to limit the effects of numerical dispersion and correctly simulate the dissolution of CO2 in brine. When CO2 dissolves in brine, the brine becomes more dense and convection may take place, enhancing the dissolution. Estimates of time to onset of convection and critical wavelength have been determined by a number of people with variable results, but for a typical aquifer, a grid resolution of less than 1m may be required. If the cell size is too large, the amount of convection will be underestimated, and therefore the amount of dissolution will be underestimated. However, our results showed that the error in the amount of CO2 dissolved was greater at the end of injection than at the end of the simulation (100 years after injection ceased), due to numerical dispersion. On the other hand, our simulations showed that the grid resolution had little effect on the build-up of pressure, which indicates that coarse grids may be sufficient for initial assessments of storage potential.
-
-
-
Accurate Discretization of Vertically-averaged Models of CO2 Plume Migration
Authors K.A. Lie, I.S. Ligaarden and H.M. NilsenWhen CO2 is injected into a deep formation, it will migrate as a plume that moves progressively higher in the formation, displacing the resident brine. The invasion front is driven by gravity, and the upward movement of the plume is limited by a low-permeable caprock. Several authors have recently proposed to make a sharp-interface assumption and only describe the plume migration in a vertically-averaged sense. For inhomogeneous permeability, the plume migration is then described by a system of conservation laws with spatially discontinuous flux. If one disregards dissolution and residual trapping, the system reduces to a scalar conservation law with a spatially dependent flux function, which may exhibit different solutions depending on the entropy condition that is enforced to pick a unique solution. We propose a certain set of assumptions that lead to the so-called minimum-jump condition and derive the corresponding solutions to the Riemann problem. Solutions to this problem are fundamental when developing accurate Godunov-type schemes. Here, we take a slightly different approach and present an unconditionally stable front-tracking method, which is optimal for this type of problem. Moreover, we verify the well-known observation that a standard upstream mobility discretization can give wrong solutions in certain cases.
-
-
-
Effects of Non-Darcy Flow in CO2 Injection into Saline Aquifers
Authors A. Mijic and T.C. LaForceOne of the most important aspects of flow in porous media is the flow velocity. It defines the convective part of the transport and usually is assumed to be given by Darcy law. However, for the analysis of near-wellbore phenomena, such as salt precipitation during CO2 injection, it is necessary to include nonlinear flow behaviour in the zone of interest. This paper presents the modification of existing solutions for the immiscible non-Darcy displacement in order to account for mutual solubility of fluids and salt precipitation. The flow is governed by the Forchheimer type equation for the two phase flow. The analytical solution is implemented in CO2 injection modelling to investigate the effects of a nonlinear flow regime on phase saturation profiles. The results enable the determination of a drying front and hence the zone where the salt precipitation will occur. By implementing the approximate solution for the solid salt saturation it is shown that non Darcy displacement has a significant influence on obtained saturation values and the skin factor than previously estimated. Finally, the non-Darcy displacement dependency on injection rates further contributes to the modification of salt precipitation pattern.
-
-
-
Multicomponent Multiphase Transport Solutions for Application to CO2 Storage
Authors A.L. Goater and T. LaForceReliable prediction of subsurface flow will be vital as Carbon Dioxide (CO2) storage grows into a global-scale industry. Therefore, it is essential to achieve fundamental objectives such as identifying and numerically modelling the physical solutions to three-phase transport problems. In this paper we present a one-dimensional (1D) conservation law for a three-phase, three-component transport problem of interest for CO2 storage in an oil reservoir. An example analytic solution for the purely advective problem is proposed. The componentwise Essentially Non-Oscillatory (ENO) method is used to minimise numerical error as compared with single-point upstream weighting (SPU). The simulated solutions both converge to the proposed dispersion-free analytic solution. However, the ENO approach exhibits only first order convergence, which is the same as SPU. Fine-grid ENO simulations of the advection/dispersion equation in three-phase fluid flow with non-negligible capillary pressures are also considered. It is demonstrated that the physically-realistic dispersion caused by capillary pressure alters the simulated solution in a fundamentally different way than numerical effects, dispersing some shocks while having little impact on others.
-
-
-
Multidimensional and Higher Order Edge Based Upwind Schemes for Hyperbolic Systems Including Gravity in Porous Media on Structured and Unstructured Grids
Authors S. Lamine and M.G. EdwardsStandard reservoir simulation schemes employ single-point upstream weighting for approximation of the convective fluxes when multiple phases or components are present. These schemes introduce both coordinate-line numerical diffusion and crosswind diffusion into the solution that is grid and geometry dependent. Families of locally conservative higher-order multidimensional upwind schemes are presented for convective flow approximation in porous media that reduce both directional and crosswind diffusion. The schemes are coupled with full-tensor Darcy flux approximations and handle general flow conditions including counter current gravity driven flows and systems of hyperbolic equations. Characteristic vector upwind approximations are proposed and compared with the simulation standard single-point upstream weighting schemes. When dealing with systems of hyperbolic equations, upwind characteristic wave decomposition is used in combination with different limiting strategies involving conservative, primitive and characteristic variables. Alternate wave vector tracing approximations are proposed based on phase velocities and characteristic velocities and comparisons are presented. The higher order multidimensional schemes are designed for general grids. Conditions for a local discrete maximum principle are presented that ensure solutions are free of spurious oscillations. Benefits of the resulting schemes are demonstrated for gravity segregated flow and polymer flood systems. The test cases involve a range of unstructured grid types with variations in orientation and permeability that lead to flow fields that are poorly resolved by standard simulation methods. The higher order multidimensional formulations are shown to effectively reduce numerical diffusion, leading to improved resolution of concentration and saturation fronts.
-
-
-
Monotone Multi-dimensional Upstream Weighting on General Grids
Authors E. Keilegavlen, J. Kozdon and B. MallisonThe governing equations for multi-phase flow in porous media often have a mixed elliptic and (nearly) hyperbolic character. The total flux for each phase consists of two parts: a geometry and rock dependent term that resembles a single-phase flux and a mobility term representing fluid properties and rock-fluid interactions. The geometric term is commonly discretized by two or multi point flux approximations (TPFA and MPFA, respectively). The mobility is usually treated with single point upstream weighting (SPU), also known as dimensional or donor cell upstream weighting. It is well known that when simulating processes with adverse mobility ratios, e.g. gas injection, SPU yields grid orientation effects. For these physical processes, the governing equations are unstable on the scale at which they are discretized, rendering a challenging numerical problem. These challenges must be addressed in order to avoid systematic biasing of simulation results and to improve the overall performance prediction of enhanced oil recovery processes. In this work, we present a framework for multi-dimensional upstream weighting for multi-phase flow with gravity on general two-dimensional grids. The methodology is based on a dual grid, and the resulting transport methods are provably monotone. The multi-dimensional transport methods are coupled with MPFA methods to solve the pressure equation. Both explicit and fully implicit approaches are considered for treatment of the transport equations. The results show considerable reduction of grid orientation effects compared to SPU, and the explicit multi-dimensional approach allows larger time steps. For the implicit method, the total number of non-linear iterations is also reduced when multi-dimensional upstream weighting is used.
-
-
-
On the Upstream Mobility Finite Difference Scheme
Authors A. Adimurthi, J. Jaffré, S. Mishra and G.D. Veerappa GowdaIncompressible two-phase flow in porous media when capillarity effects are neglected can be modeled as a nonlinear scalar conservation with unknown one of the phase saturations. The standard scheme used in reservoir simulation or in hydrogeology uses the upstream mobility flux. The solution calculated with such a scheme was proven to be converging to the entropy solution. This scheme is used in the homogeneous case as well as in the case of a change in rock types when the mobilities functions are changing from one rock type to the other. However in this latter case, while the scheme is giving the right solution most of the time, we will exhibit some examples of mobility functions for which the scheme does not give the right solution. An alternative scheme is to use an extension of the Godunov flux for which we give a simple formula which can be used when a change of rock types occurs. Thus the solution given by this scheme can be proven to be correct even in cases where the upstream mobility scheme fails.
-
-
-
Finite Volume Schemes for Multiphase Flow Simulation on Near Well Grids
Authors C. Guichard, J. Brac, R. Eymard and R. MassonIn reservoir engineering a proper well modeling requires to simulate accurately multiphase flow taking into account the singular pressure distribution in the well vicinity and the large difference of scales between the wellbore radius and the reservoir dimension. This paper investigates the numerical simulation of a 3D near-well model using a radial mesh exponentially refined down to the well boundary. The radial and the reservoir CPG grids are matched using either hexahedra or both tetrahedra and pyramids. Various multipoint finite volume methods are compared such as the O and L schemes, and the newer G and GradCell schemes which are briefly described. Our first objective is to introduce a new scheme which includes the advantages of compact stencil, symmetry and convergence. Our second objective is to study the behavior of these schemes on the above two types of meshes. First, their convergence behavior are compared for the deviated well single-phase flow analytical solution described in the literature with an homogeneous anisotropic permeability tensor. Second, the influence of the type of mesh is studied on a CO2 injection test case in an aquifer with solubility of CO2 in the water phase and anisotropic permeability fields.
-
-
-
Horizontal Simulation Grids as Alternative to Structure-based Grids for Thin Oil-zone Problems
By O. PettersenAs a general rule, the layering in reservoir simulation grids is based on the geology, e.g. structure tops. In this paper we investigate the alternative of using horizontal layers, where the link to the geology model is by the representation of the petrophysics alone. The obvious drawback is the failure to honor the structure in the grid geometry. On the other hand a horizontal grid will honor the initial fluid contacts perfectly, and horizontal wells can also be accurately represented. Both these issues are vital in thin oil-zone problems, where horizontal grids may hence be a viable alternative. To investigate this question, a number of equivalent simulation models were built for a segment of the Troll Field, both geology-based and horizontal, and various combinations of these. In the paper it is demonstrated that the horizontal grid is able to capture the essentials of fluid flow with the same degree of accuracy as the geology-based grid, and near-well flow is considerably more accurate. For grids of comparable resolution, more reliable results were obtained by a horizontal grid than a geo-grid. A geo-grid with local grid refinement and a horizontal grid produced almost identical results, but the ratio of computing times was more than 20 in favor of the horizontal grid. In the one-phase regions of the reservoir, relatively coarse cells can be used without significant loss of accuracy.
-
-
-
A Model for Conductive Faults in Heterogeneous Media for Non-matching Grids
Authors X. Tunc, I. Faille, T. Gallouët, M.C. Cacas and P. HavéIn the standard Corner Point Grid approach, a fault zone is supposed to be sufficiently thin to be represented as a surface with non-matching interfaces, due to the throw, across the fault. Following this point of view, we study an approach where fluid flow along the fault is modelized by a lower dimensionnal model, the fault thickness becoming a model parameter.A fault zone is represented by two sets of faces, corresponding to the fault surface, each set being conforme with its neighbouring matrix grid-blocks. Fault zone properties are assigned to each face. Finite Volume Discretization of fluid flow equations leads to an additional mass balance equations in each fault face, with three types of fluxes. The first one represents mass exchange between fault faces along each set. The second one is between neighbouring non-conforme fault faces across the fault and finally, the last one, is between fault faces and their neighbouring matrix grid-blocs. In this paper, we mainly focus on how to handle the non-matching grids and therefore consider only one-phase fluid flow. However, the extension of this approach to complex multiphase flow is straightforward as it enters into the Finite Volume framework.
-
-
-
Modification of a Reservoir Grid to Achieve Monotonicity of the Control Volume Method
Authors R.A. Yorgova and I. AavatsmarkSpurious oscillations may occur for nonmonotone discrete methods. Therefore, monotonicity is a desired property for control-volume methods. Control-volume MPFA methods are only conditionally monotone. Sufficient criteria for monotonicity of the MPFA O-method on quadrilateral grids in 2D for the pressure equation has previously been developed. Tests on rough quadrilateral grids show that the monotonicity conditions are often violated. In this work we present modification of rough quadrilateral grids on homogenous and on heterogeneous media obtaining grids where the MPFA O-method is monotone. The reconstruction is based on short moves of the grid nodes. The moves are defined by weighted gradients on functions based on the monotonicity criteria. The main steps of the reconstruction algorithm are described. Examples of rough quadrilateral grids and their modified grids are included in order to show the efficiency of the presented algorithm.
-
-
-
Practical Gridding Algorithms for Discrete Fracture Modeling Workflows
Authors B.T. Mallison, M.H. Hui and W. NarrThis paper describes grid generation algorithms for simulating flow and transport in fractured porous media. The methods are designed to only capture details of the fracture network geometry larger than the specified grid resolution. The final grids honor fractures approximately while maintaining good cell quality. Improved representations of the input geometry can be obtained by generating a new grid with finer resolution. Several numerical examples are presented in two and three dimensions that demonstrate that the algorithms are robust and practical for industrial applications.
-
-
-
Accurate Locally Conservative Discretizations for Modeling Multiphase Flow in Porous Media on General Hexahedra Grids
Authors M.F. Wheeler and G. XueFor many years there have been formulations considered for modeling single phase ow on general hexahedra grids. These include the extended mixed nite element method, and families of mimetic nite difference methods. In most of these schemes either no rate of convergence of the algorithm has been demonstrated both theoret- ically and computationally or a more complicated saddle point system needs to be solved for an accurate solution. Here we describe a multipoint ux mixed nite element (MFMFE) method [5, 2, 3]. This method is motivated from the multipoint ux approximation (MPFA) method [1]. The MFMFE method is locally conservative with continuous ux approximations and is a cell-centered scheme for the pressure. Compared to the MPFA method, the MFMFE has a variational formulation, since it can be viewed as a mixed nite element with special approximating spaces and quadrature rules. The framework allows han- dling of hexahedral grids with non-planar faces by applying trilinear mappings from physical elements to reference cubic elements. In addition, there are several multi- scale and multiphysics extensions such as the mortar mixed nite element method that allows the treatment of non-matching grids [4]. Extensions to the two-phase oil-water ow are considered. We reformulate the two- phase model in terms of total velocity, capillary velocity, water pressure, and water saturation. We choose water pressure and water saturation as primary variables. The total velocity is driven by the gradient of the water pressure and total mobility. Iterative coupling scheme is employed for the coupled system. This scheme allows treatments of different time scales for the water pressure and water saturation. In each time step, we first solve the pressure equation using the MFMFE method; we then Center for Subsurface Modeling, The University of Texas at Austin, Austin, TX 78712; [email protected]. yCenter for Subsurface Modeling, The University of Texas at Austin, Austin, TX 78712; [email protected]. 1 solve the saturation using discontinuous Galerkin (DG) method by taking multiple small time steps within the large time step. In addition, the MFMFE method allows effcient computations of the total and capillary velocity since the method gives the local velocity approximation in terms of surrounding pressure degrees of freedom. Both theoretical and computational results are discussed and presented. Exten- sions to advection-diffusion equations and non-Newtonian polymer ooding [6] are also discussed.
-
-
-
Quasi-positive Continuous Darcy-flux Approximation on Cell-centred Triangular Grids
Authors M.G. Edwards and H.A. FriisA new family of cell-centered finite-volume schemes is presented for solving the general full-tensor pressure equation on arbitrary unstructured triangulations. The new schemes are flux continuous and have full pressure support (FPS) over each subcell with continuous pressure imposed across each control-volume sub-interface, in contrast to earlier formulations. The earlier methods are point-wise continuous in pressure and flux with triangle-pressure-support (TPS) leading to a limited quadrature range. An M-matrix analysis identifies bounding limits for the schemes to posses a local discrete maximum principle. Positive definite conditions are also given. Comparisons show that the earlier pointwise TPS methods can induce strong spurious oscillations in solutions for problems involving strong full-tensor anisotropy where the M-matrix conditions are violated, leading to unstable solutions in such cases. The earlier unstructured cell-centered TPS schemes are investigated in terms of stability and are shown to be decoupled at high anisotropy ratio. In contrast to TPS, the new FPS formulation leads to well resolved stable solutions that are essentially free of spurious oscillations. A substantial degree of improved convergence behavior, for both pressure and velocity, is observed in all convergence tests.
-
-
-
Convergence of the MPFA L-method – Strengths and Difficulties
More LessThe convergence of the Multi Point Flux Approximation (MPFA) O-method on general grids has recently been proved by rewriting the method in a weak form using quadratures. We will examine in detail the quadrature of the MPFA L-method in 2D. To obtain the correct quadrature one must introduce auxiliary variables that impose the continuity of the pressure over an interface – a characteristic of the L-method. The auxiliary variables can successively be eliminated from the quadrature. Compared with the O-method, the L-method obtains a stronger diagonalization of the quadrature matrix. This is important as the convergence proof of the method relies on the assumption that the symmetric part of the quadrature matrix has positive eigenvalues. On the other hand, the diagonalization also means that the L-method has a stronger inter-element coupling, as pressure continuity is imposed along element interfaces and not within the partial elements themselves. We will examine and compare the assumption for convergence for the O-method and the L-method in 2D aiming for a better understanding of the merits of the two methods and how anisotropy and grid deformation influence their behaviour.
-
-
-
On Rough Grids – Convergence and Reproduction of Uniform Flow
Authors R.A. Klausen and A.F. StephansenIn reservoir simulation the flow driven by pressure differences and gravity is usually approximated by means of the empirically derived Darcy's law. The elliptic pressure equation that results is preferably solved with a mass conservative method, for which several candidates are available. We distinguish between full field methods and raw field methods. On the one hand we have continuous full field methods which approximate the velocity with a field that is defined inside the element itself. On the other hand we have discrete methods that approximate the flux, and for which no velocity field is defined inside the elements. We denote this raw field. We will examine some mass-conservative methods on the basis of two different characteristics: uniform flow reproduction and convergence on rough grids. Analysis apart, we propose an interpolation on polygons that can reproduce uniform flow by using edge basis functions constructed by means of barycentric functions. This permits us to obtain an H(div) field starting from discrete fluxes evaluated at the interfaces of the elements. We also show that a similar interpolation can be constructed to obtain an H(curl) field on polygons in 2D.
-
-
-
Locally Adaptive Timestepping in Reservoir Simulators
Authors H. Mc Namara, G. Bowen and P.J. DellarAn algorithm for locally adapting the step-size for large scale finite volume simulations of multi- phase flow in petroleum reservoirs is suggested which allows for an “all-in-one” implicit calculation of behaviour over a very large time scale. Some numerical results for simple two-phase flow in one space dimension illustrate the promise of the algorithm, which has also been applied to very simple 3D cases. A description of the algorithm is presented here along with early results. Further development of the technique is hoped to facilitate useful scaling properties.
-
-
-
Improved Gravity Splitting for Streamline and Reordering Methods
Authors H.M. Nilsen and J.R. NatvigFor heterogeneous reservoirs, fast, stable and accurate methods are hard to obtain. Large changes in the velocity field leads to severe time step restrictions in explicit schemes or expensive time steps in implicit schemes. In the absence of gravity, the exact velocity field will be loop free in the sense that there are no closed integral curves. For streamline methods, the absence of closed integral curves ensures that all streamlines have endpoints in wells. Likewise, this property implies that the Jacobian matrix of an implicit scheme with an upwind flux approximation can be reduced to a triangular matrix by a permutation, or reordering, of the unknowns. For the above methods the effect of gravity is usually handled by operator splitting in the transport equation. Gravity will introduce rotation in the total velocity, which may yield closed integral curves. Even when the effect of gravity is small, it can limit the efficiency and the robustness of streamline and reorder methods. To overcome this problem, we split the total velocity field in a loop-free part without the effect of gravity, and a part driven by gravity only. This is achieved by solving the pressure equation two times with different righthand sides.
-
-
-
A Robust Streamline Tracing Method for Systems with Non–neighbor Connections
By D. KachumaThe use of complex corner point systems in field flow reservoir simulators allows the modelling of a wider variety of cell geometries and hence a better representation of faulted reservoirs. If faults are involved, then a cell face may be next to two or more cell faces. This phenomenon can also be observed in situations where the cells in a particular region have been refined to improve accuracy. This creates complex connections between cell faces which require special attention and representation within the simulator. Conventionally, this is handled by the use of the so-called non-neighbour connections (NNC) whereby the simulator includes, on top of the regular neighbour fluxes, the fluxes between the extra cell connections. We present a robust technique which is applicable regardless of the flow configuration. Our technique involves subdividing each cell with non-neighbour fluxes into a local regular subgrid. Using this grid the local streamfunction is constructed numerically. By reducing the size of the local subgrid, the accuracy of the traced streamlines improves. We demonstrate this technique first on a detailed 2D synthetic example that is designed to show the robustness of the technique. The practical utility of our algorithm is then demonstrated in a structurally complex and heavily faulted full model of a North Sea field which includes cells with several non-neighbour configurations in different faces. This treatment is contrasted with the usual approach and also other techniques designed to trace streamlines in heavily faulted systems.
-
-
-
Streamline Tracing on Irregular Geometries
More LessThe accurate and efficient tracing of streamlines is fundamental to any streamline-based simulation method. For grids with irregular cell geometries, such as corner-point grids with faults or Voronoi-diagram (pebi) grids, most efforts to trace streamlines have been focused on subdividing irregular cells into sets of simpler subcells, typically hexahedra or simplices. Then one proceeds by reconstructing the velocity field on, and tracing through, the sets by a more basic algorithm. One such basic algorithm applies to incompressible flow on simplices. In that case there is a cell-wise constant velocity that is consistent with given face fluxes, as long as those face fluxes sum to zero for the cell. We give an efficient and simple formulation of this algorithm using barycentric coordinates. Another approach to irregular cell tracing is computing the streamline directly on the complex cell geometry. We give a new method based on generalized barycentric coordinates for direct tracing on arbitrary convex polygons, which generalizes the corner velocity interpolation method of Hægland et al (2007). The method generalizes to convex polytopes in 3D, with a restriction on the polytope topology near corners that is shown to be satisfied by several popular grid types.
-
-
-
A Stochastic Sharpening Method for the Propagation of Phase Boundaries in Multiphase Lattice Boltzmann Simulations
Authors T. Reis and P.J. DellarExisting lattice Boltzmann models that have been designed to recover a macroscopic description of immiscible liquids are only able to make predictions that are quantitatively correct when the interface that exists between the fluids is smeared over several nodal points. Attempts to minimise the thickness of this interface generally leads to a phenomenon known as lattice pinning, the precise cause of which is not well understood. This spurious behaviour is remarkably similar to that associated with the numerical simulation of hyperbolic partial differential equations coupled with a stiff source term. Inspired by the seminal work in this field, we derive a lattice Boltzmann implementation of a model equation used to investigate such peculiarities. This implementation is extended to different spacial discretisations in one and two dimensions. We shown that the inclusion of a quasi-random threshold dramatically delays the onset of pinning and facetting.
-
-
-
Parallel Sparsified Solvers for Reservoir Simulation
By H.M. KlieThe premise of the present work lies on the fact that a large number of system coefficients may be disregarded without compromising the robustness of the overall solution process. That is, given a linear system it is possible to construct a preconditioner by dropping a large number of off-diagonal nonzeros and use it as a suitable proxy to approximate the original matrix. This proxy system can be in turn factored to generate a new class of block ILU preconditioners, approximated inverses and algebraic multigrid implementations. We propose two parallel algorithms to sparsify a given linear system: (a) random sampling sparsification (RSS), and (b) percolation-based sparsification (PBS). The former one relies on the idea that coefficients are included into the sparsified system with a probability proportional to its effective transmissibility. The latter relies on capturing highly connected flow paths described by whole set of transmissibility coefficients. Depending on the case, the RSS and PS algorithms have the potential to reduce in orders of magnitude the number of floating point operations associated with the preconditioner. Results confirming the benefits of sparsified solvers are illustrated on a wide set of field cases arising in black-oil and compositional simulations.
-
-
-
High Performance Manycore Solvers for Reservoir Simulation
More LessThe forthcoming generation of many-core architectures compels a paradigm shift in algorithmic design to effectively unlock its full potential for maximum performance. In this paper, we discuss a novel approach for solving large sparse linear systems arising in realistic black oil and compositional flow simulations. A flexible variant of GMRES (FGMRES) is implemented using the CUDA programming model on the GPU platform using the Single Instruction Multiple Threads (SIMT) paradigm by taking advantage of thousands of threads simultaneously executing instructions. The implementation on the GPU is optimized to reduce memory overhead per floating point operations, given the sparsity of the linear system. FGMRES relies on a suite of different preconditioners such as BILU, BILUT and multicoloring SSOR. Additionally, the solver strategy relies on reordering/partitioning strategies algorithms to exploit further performance. Computational experiments on a wide range of realistic reservoir cases show a competitive edge when compared to conventional CPU implementations. The encouraging results demonstrate the potential that many-core solvers have to offer in improving the performance of near future reservoir simulations.
-
-
-
Multiple Kernel Learning Approach for Reservoir Modelling
Authors V. Demyanov, L. Foresti, M. Kanevski and M. ChristieThe paper proposes a novel data driven approach for modelling petrophysical properties in oil reservoir. We aim to improve realism of reservoir models with a more intelligent way of integrating the raw data and geological knowledge. Multiple Kernel Learning (MKL) provides enhanced interpretability of the model by using separate kernels for input variables performs kernel/feature selection to solve a regression problem in a high-dimensional feature space. MKL has an advantage of rigorous control over the model complexity to achieve the right balance between data fit and prediction accuracy. The MKL reservoir model was designed to integrate data and prior knowledge, which describe geological structure at multiple scales. Geological structures can be detected by applying convolution filters on noisy seismic data to capture changes in gradients, orientations, sizes of meandering channels. Such "geo-features" are added as input variables into the MKL model, which optimises the weighted combination of kernels to fit to the available data. MKL application to a synthetic meandering channel reservoir has demonstrated capacity of selecting the relevant input information for detecting the channel structure. Experiments with noisy seismic inputs highlighted feature selection skills of MKL which was able to filter them out.
-
-
-
MP Simulations Without Computing MP Statistics
Authors G. Mariethoz, P. Renard, J. Straubhaar and J. CaersIn recent years, multiple-point simulation has become an invaluable tool to integrate geological concepts in subsurface models. However, due to the high CPU and RAM demand, its use is restricted to relatively small problems with limited structural complexity. Moreover, it only allows for the simulation of univariate fields. We present an alternative method that produces conditional realizations honoring the high-order statistics of uni- or multivariate training images. It is based on a sampling method introduced by Shannon (1948), strictly equivalent to the original method of Guardiano and Srivastava (1993), but that does not need to compute conditional probabilities and to store them. In the sampling process, we use a distance between data configurations that allows simulating both discrete and continuous variables. As a result, the simulation algorithm is drastically simplified and has more possibilities. Since nothing is stored, neighborhoods can have virtually any size. Moreover, the neighborhoods are not restricted to a template, making multiple-grids unnecessary. Multivariate data configurations can be considered, allowing to generate realizations presenting a given multivariate multiple-point dependence. In addition to having virtually no RAM requirement, the method is straightforward to parallelize. Hence it can produce very large and complex realizations.
-
-
-
Application of Stochastic Partial Differential Equations to Reservoir Property Modelling
Authors R. Potsepaev and C.L. FarmerExisting algorithms of geostatistics for stochastic modelling of reservoir parameters require a mapping (the 'uvt-transform') into the parametric space and reconstruction of a stratigraphic co-ordinate system. The parametric space can be considered to represent a pre-deformed and pre-faulted depositional environment. Existing approximations of this mapping in many cases cause significant distortions to the correlation distances. In this work we propose a coordinate free approach for modelling stochastic textures through the application of stochastic partial differential equations. By avoiding the construction of a uvt-transform and stratigraphic coordinates, one can generate realizations directly in the physical space in the presence of deformations and faults. In particular the solution of the modified Helmholtz equation driven by Gaussian white noise is a zero mean Gaussian stationary random field with exponential correlation function (in 3-D). This equation can be used to generate realizations in parametric space. In order to sample in physical space we introduce a stochastic elliptic PDE with tensor coefficients, where the tensor is related to correlation anisotropy and its variation is physical space.
-
-
-
Spartan Random Fields and Applications in Spatial Interpolation and Conditional Simulation
More LessGeostatistics plays an important role in the estimation and simulation of reservoir parameters, and in forecasting reservoir production. This presentation focuses on Spartan spatial random fields (SSRFs), which provide a new framework for geostatistical applications. The idea motivating SSRF development is that spatial correlations can be represented by local interactions. A brief overview of SSRF mathematical properties will be given. SSRF variogram models, which in the isotropic case contain 3-4 parameters and thus offer more flexibility than standard models, will be presented. It will be argued that empirical variogram estimation is not necessary to infer SSRF parameters. SSRF spatial interpolation of scattered data will be discussed. It will be shown that approximate but closed-form, efficiently calculable expressions become possible. Hence, numerical solution of kriging linear systems is avoided, leading to improved numerical complexity. Connections will be drawn between to Markov random fields, statistical physics, minimum curvature estimation, and local random fields. An SSRF application to 2D conditional simulation will be presented. Finally, directions for the future development of SSRFs will be discussed, including extensions to non-Gaussian data distributions, using discretized random field models with “spin” interactions.
-
-
-
Isometric Unfolding of Stratigraphic Grid Units for Accurate Property Populating – Mathematical Concepts
Authors S. Horna, C. Bennis, H. Borouchaki, C. Delage and J.F. RainaudIn traditional methods used to populate stratigraphic units, the distortions can be very important and affect the setting up of the static and dynamic parameters necessary to the reservoir simulation, hence the simulation results. These distortions result from the mapping between the original curvilinear stratigraphic grid and the intermediate cartesian grid in which the property populating is processed. To minimize the deformation and improve the populating process, we propose a new original isometric unfolding process based on the minimization of the elastic tensor deformation. This method could be applied for every type of deposit: horizontal, parallel to top, parallel to bottom, proportional. Starting from a structural model defined into a coordinate line grid, the user chooses a reference iso-chronological level represented by a triangulated surface. The contacts between this surface and fault surfaces are explicitly extracted as coincident edges. These coincident edges are used to constraint an unfolding process, respecting the above constraints, join the fault lips opened by geological tectonic events. In this paper we will focus on mathematical concepts of the algorithms used in the whole unfolding process and present some results throw actual case studies.
-
-
-
History Matching of a Stochastic Multifractal Subesismic Fault Model
Authors M. Verscheure, A. Fourno and J.P. ChilèsMany geostatistical methods have been developed to generate realistic models of subseismic fault networks. The problem is that these models are difficult to history match. In this work, we present an original methodology in which history matching is performed through a modification of fault positions. We first propose a multifractal methodology to generate the faults. The method has been specifically developed to allow history matching. The model parameters are derived from the analysis of the seismic faults. A stochastic algorithm is used to generate 3D subseismic fault realizations that are constrained to the seismic faults. The fracture network is then discretized on a dual porosity simulation grid. Equivalent flow parameters are computed using an analytical method. Last, full field simulations are performed using a multiphase flow simulator. Then, we introduce a method to gradually change the locations of faults while preserving multifractal properties. Changes in locations are driven from a reduced number of parameters. These parameters are gradually modified to optimize the geometry of the realization and compel the initial fault model to reproduce the hydrodynamic behaviour observed on the field. The potential of the methodology is demonstrated on a 3D case study.
-
-
-
Fault Displacement Modelling Using 3D Vector Fields
Authors P. Røe, F. Georgsen, A.R. Syversveen and O. LiaIn history matching and sensitivity analysis it is useful to study the effect of changing the position of horizons near faults and the resulting facies juxtaposition due to change in fault displacements. A new algorithm for calculating a 3D displacement field is developed. This is applicable to a wide range of faults due to a flexible representation and gives the possibility to apply this field to change the displacement and thereby moving horizons. The fault is represented as a surface in a local coordinate system defined by a centre point, dip and strike angles. The displacement of points associated with the fault outside the fault surface is described by a 3D vector field. The displacement on the fault surface can be found from the intersection lines between horizons and the fault surface (fault lines). Away from the fault surface the displacement field is defined by a monotonic decreasing function. The displacement for the fault can be changed by a user-defined factor. Then the whole displacement field is changed and points on horizons around the fault are moved by applying the modified displacement field. The interaction between faults influencing the same points are taken care of by truncation rules.
-