A sensitivity study of relevant parameters on sand production in unconsolidated sandstone reservoirs

Sand production is a prevalent problem during oil and gas production from weakly consolidated or unconsolidated formations. It can erode downhole equipment and surface facilities, cause pipeline blockage, leakage and damage casing and generate additional need for sand disposal. Decision for appropriate sand control strategy requires engineering analysis of the key parameters affecting sand production to evaluate timing and severity of sanding over the life of field conditions. This paper presents the effects of reservoir and geomechanical parameters including well trajectory, poroelastic stress coefficient, Biot’s factore, maximum horizontal stress, horizontal stress anisotropy ratio, cohesive strength and Uniaxial Compressive Strength (UCS) on sand production from openhole wells. The results indicated that the critical bottomhole flowing pressure increases with increasing the poroelastic stress coefficient, Biot’s factor, maximum horizontal stress and horizontal stress anisotropy ratio, but decreases with increasing the cohesive strength and UCS of rocks. Furthermore, the results show that the wellbore inclination and azimuth have a significant role in sanding potential during production. Also, it is concluded that in normal stress regime the critical bottom hole flowing pressure of a horizontal borehole is greater than the vertical borehole, so the vertical boreholes have less potential for sanding than the horizontal boreholes and almost all the deviated wells. Accurate prediction of the conditions for sand production is critical to the design of cost effective completion.


INTRODUCTION
Production of oil and gas from hydrocarbon-bearing reservoirs can result in a reduction of the reservoir pore pressure (formation pressure) unless pressure support is provided from a gas cap or an aquifer.A reduction in the reservoir pressure, in turn, results in a reduction in the stresses acting within the reservoir.This stress-depletion response is known from simple theoretical calculations and has been observed in field data.Wells drilled into depleting fields will experience this change in the stress acting within the field.This has important implications for well design as an increasing number of wells are being designed as openhole or barefoot completions.Well design is also becoming increasingly customised in order to minimise the number of wells required to drain a field.*Corresponding author.E-mail: m.hayavi2013@gmail.com.Tel: +98 939 183 5971.
Author(s) agree that this article remain permanently open access under the terms of the Creative Commons Attribution License 4.0 International License Consequently, large horizontal well sections are drilled with multilateral branches -several horizontal well sections drilled from one producing well.An important aspect of the design for these customised and multilateral wells is the openhole stability over the lifetime of the field, which is dependent on how the reservoir stresses will change with depletion (Addis, 1997).
In the petroleum industry, the term -sand production‖ refers to the production of solid particles together with the formation fluids from depleted reservoirs.This phenomenon is common mainly in weak sandstone reservoirs and is a possible consequence of the degradation of the mechanical properties of the rock surrounding the wellbore caused by drilling, completion, and production operations particularly, during the production phase, the decrease of pore pressure causes a concentration of stresses around the wellbore and the perforation tips which, in turn, can lead to the failure of the rock.When the right conditions are met, e. g. production rates sufficiently high, then the failed rock can be mobilized and the fluid flow can drag to the surface grains, particles or aggregates of the damaged rock (Volonté et al., 2010).
Sanding can erode downhole equipment and surface facilities, cause pipeline blockage, leakage and damage casing due to formation subsidence, lead to more frequent well intervention and workovers, and generate additional need for sand disposal (Addis, 1997).Decision for appropriate sand control strategy requires engineering analysis of the key parameters affecting sand production to evaluate timing and severity of sanding over the life of field conditions.
Sand production mechanisms can be summarized into the following points (Wang et al., 2005): (1) Shear failure induced by fluid pressure drawdown can lead to the breaking of sand grain bonds and the alteration of the material's mechanical properties; (2) Tensile failure caused by high hydrocarbon production rates can lead to dilation of solid skeleton and the loss of solid particles mechanical interactions through disaggregation; (3) High stresses due to completion cause the formation to fail (in compression) whereas fluid viscous drag forces bring the failed materials from the perforation tunnels into the wellbore.
There are many factors that must be considered to obtain a comprehensive understanding of how and why sand production occurs and such factors include: (1) Geological factors.(2) Rock composition.(3) Mechanical factors.(4) Drilling practices.(5) Production operations (Joseph et al., 2012).Tao et al. (2008) reported that as the flow rate is increased, the wellbore pressure decreases and the effective tangential stress in the elastic zone increases.On the other hand, the effective stresses in the plastic zone decrease with increasing flow rates due to the increase in the pore pressure gradient.
Increasing the flow rate can initiate sanding, therefore the critical production rate at which the cavity becomes unstable should be determined during oil and gas production in weakly consolidated or unconsolidated formations.Morita et al. (1987) by separating the effect of well pressure and local pressure gradient around the cavity, proposed an analytical approach to study the effects of many parameters on sand production.It was pointed out that the following parameters may affect sand production: (1) Boundary loads such as well pressure and in-situ stresses; (2) Fluid flow induced force which is dependent on such factors as flow rate, permeability, viscosities of fluids, relative permeability for multiphase flow and fluid saturation, (3) rock deformation character, (4) rock strength character, (5) perforation cavity geometry and shot density, (6) cyclic loading history.
Sanding onset prediction involves stress calculation at cavity (including wellbore or perforation tunnel) surface.Even though a numerical model, such as finite element model, is more general, analytical or semi-analytical models may be more convenient and easier to use under special conditions.Besides, an analytical model is always useful to verify numerical models.In petroleum engineering, the vertical/horizontal wellbore, perforation tunnel and their adjacent formation are often approximated as thick-walled hollow cylinder.Using this approximation, it can be able to obtain an analytical or semi-analytical solution for the near wellbore/perforation tunnel stress state and use it in sand production prediction (Yi, 2003).
Drilling the wellbore or creating any cavity like perforations, changes the stress pattern in the medium around the cavity (Dai et al., 2014).Increase of drawdown augments effective stresses in an interval around the wellbore.This is attributed to the fact that pore pressure recovers much more slowly compared to total stresses.Therefore, shear stresses increase around the cavities once higher drawdowns are used.In this respect, acts depletion very similarly to drawdown in increasing shear stresses around the cavities.Shear failure induced by fluid pressure drawdown takes place once shear stresses exceed limit shear strength of the intact rock and can lead to the breaking of sand grain bonds and the alteration of the material's mechanical properties.Shear failure mechanism is mainly active around the cavities where two major criteria are fulfilled.First, shear stresses are very high and second, differential deformations are possible (Nouri et al., 2003).
Sanding related studies can go back to the 1920's, and before the 1970's, most studies were associated with sand control, and then the interests shifted to the determination of onset of sand production and a Mohr-Coulomb (M-C) type shear failure was initially postulated (Luo et al., 2012).
The objective of this paper is to present the sensitivity analyses of reservoir and geomechanical parameters on the sand production from open hole wellbores.

STRESS CONCENTRATION AROUND A WELLBORE AT PRODUCTION CONDITION
The stress concentration around a well drilled in an isotropic, elastic medium under anisotropic in-situ stress condition (Maximum and minimum horizontal stresses are different) is described by the Kirsch equations.The general expressions for the stresses at the wellbore wall for a deviated well in the production situation are Fjaer et al. (2008): Where is the radial stress, is the tangential (hoop) stress, is the axial stress induced around the wellbore, P wf is the bottomhole flowing pressure, θ is measured from the azimuth of maximum horizontal stress (Degree) , is the Poisson's ratio and is the farfield pore pressure.The shear stresses at the wellbore wall are denoted , and , while the in-situ stresses in (x, y, z) coordinate system, denoted , , , , and , and they are defined as (Al-Shaaibi et al., 2013): (2) Where and are the maximum and minimum horizontal stresses, respectively, is the vertical stress, i is wellbore inclination and is the azimuth angle due to the maximum horizontal stress ( ) direction (Degree).Figure 1 shows the stress transformation system in a deviated borehole where is the rotation angle around the z′-axis (measured from the x′-axis) and i is the rotation angle around the y′-axis (measured from the z′axis) and is the poroelastic stress coefficient defined as (Al-Shaaibi et al., 2013): (3) Where is Biot's factor.The shear failure known as sanding is expected to happen at the wellbore wall and at the point of maximum tangential stress (θ=90 °) where the rock is under maximum compression.For a vertical borehole, the inclination angle (i) is set to zero and the x-axis is oriented so that it coincides with the major horizontal principal stress axis (that is, =0 °).However, for a vertical well the maximum stress values will always be at θ=90° for any values of the in-situ stresses and Equation (1) becomes: Where is the effective radial stress, is the effective tangential (hoop) stress, is the effective axial stress induced around the wellbore.For a laterally large reservoir compared to its thickness, the change in vertical stress is considered negligible and therefore it is usually kept constant (Rahman et al., 2008).The maximum and minimum horizontal stresses are updated as follows, respectively: (5) (6) Where ( 7) and and are the maximum and minimum horizontal stresses at current production condition, respectively.and are the initial and current reservoir pressures, respectively.

SANDING ONSET PREDICTION MODEL
To predict the well pressure at which sanding will occur the stress equations at the borehole wall must be compared against the failure criterion.The stress equations based on linear poroelastic constitutive laws for deviated and vertical wellbores were introduced earlier.These equations will be inserted in to the Mohr-Coulomb criterion to obtain the critical wellbore pressure for sanding onset (Pw) (Al-Shaaibi et al., 2013).
The Mohr-Coulomb criterion is the simplest and most commonly used strength criterion for geomaterials.The Mohr-Coulomb criterion can be expressed in terms of shear stress and effective normal stress on the shear plane as follows (Zhang et al., 2010): Where is the shear stress, the effective normal stress, and and , respectively, the cohesive strength and the internal friction angle of the rock.The Mohr-Coulomb criterion can also be expressed in terms of the major and minor effective principal stresses, that is, the tangential stress and the radial stress: Where the uniaxial compressive strength of the rock and q is a parameter related to internal friction angle ( ): (10) Introducing Equations 4 into 9 (major and minor effective principal stresses are tangential and radial stresses, respectively), the critical bottom hole flowing pressure for prevention of sand production from a open hole wellbore is: (11)

CONSTRUCTION OF MECHANICAL EARTH MODEL
A Mechanical Earth Model (MEM) is an explicit description of the mechanical properties of the reservoir and overburden formations, including rock strength and elastic properties, the state of in-situ stresses and their direction as well as pore pressure.It forms the basis for any geomechanical analysis, such as wellbore stability analysis, sanding prediction evaluation, hydraulic fracture design, mechanical characterization of fractured formation, fault seal evaluation, reservoir compaction and subsidence evaluation, etc.Therefore, the reliability of sanding prediction analysis largely relies on robustness of the MEM (Mohiuddin et al., 2009).The field and laboratory data required for this study are provided from a number of production wells in one of fields located in the south west of Iran.

Calculation of rock mechanical properties
The mechanical properties of formations and dynamic elastic constants of subsurface rocks can eb derived from the measurement of elastic wave velocities and density of the rock.Sonic logging and waveform analysis provide the means for obtaining continuous measurements of compressional and shear velocities.These data, in conjunction with a bulk density measurement, permit the in-situ measurement and calculation of the mechanical properties of the rock.The elastic moduli relationships, in terms of elastic wave velocities (or transit times) and bulk density can be calculated from following equations (Alipour and Mirzaahmadian, 2010): Where is the dynamic Poisson's ratio , E d is the dynamic Young's modulus (psi), Δt s is shear wave travel time (ft/s), Δt c is compressional wave travel time (ft/s), k B is dynamic bulk modulus (psi), K R is the rock modulus (psi), ρ b is the bulk density (gr/cm 3 ), ρ gr is the grain density (gr/cm 3 ), and -b‖ is the constant coefficient which is equal to 1.34*10 10 .For the Asmari formation of mentioned oilfield, an equation developed for estimation of shear wave travel time by Nabaei et al. (2009) was used: (13) Dynamic data cannot directly be utilized to develop mechanical models.So, they should be first converted into static data through some calculation changes made and then used in geomechanical model (Wang, 2000).Poisson's ratio and static Young's modulus are both calculated via the following relations in south west of Iran.The results show good conformity with laboratorial data (Abdideh and Fathabadi, 2013): Where the static Poisson's ratio and Es is the static

In-situ stresses and pore pressure
In-situ stress magnitudes and orientations play a very important role in geomechanical analysis, and they are the most basic parameter inputs in analysis of sand production.Vertical stress is induced by the weight of the overlying formations.The vertical stress can be calculated by integration of rock densities from the surface to the depth of interest based on Equation 16.In fact, density log can be used to calculate overburden stress: Where σ v is vertical stress (psi), z is depth of interest (ft), ρ(z) is the density as a function of depth (gr/cm 3 ), g is gravitational acceleration (ft/s 2 ) and ̅ is the mean overburden density of rocks (gr/cm 3 ).Rocks of Asmari formation have an average density of 2.5 g/cm 3 .By considering horizontal strain and deformation effect, Hooke's law can be applied to derive the horizontal stresses and strains relationships (Perchikolaee et al., 2010).The following equations are obtained, and are used to calculate the minimum and maximum horizontal stress with tectonic strain effects (Al-Qahtani and Rahim, 2001).
Where ε 1 and ε 2 are strains due to tectonic forces in maximum and minimum directions and considered 1 and 1.5, respectively.Based on drilling information pore pressure gradient is estimated 0.365 psi/ft.Figure 2 shows the stresses and pressure profiles in the study area.The study area is associated with a normal faulting stress regime, where the maximum principal stress is vertical stress ( h < H < v ).

Uniaxial compressive strength
The uniaxial compressive strength (UCS) is defined as following correlation based on rock lithology and core measurements for Asmari formation, which is described by an exponential relationship between uniaxial compressive strength (psi) and sonic wave travel time (ft/s) (Khaksar et al., 2009):

Effect of cohesive strength
The behavior of sand free envelope plot based on the Critical Bottom Hole Flowing Pressure (CBHFP) of wellbore for series of cohesive strength from 290 to 1370 psi at the depth of 8700 ft are shown in Figure 3. Table 1 shows the geomechanical and reservoir properties of target depth.It can be seen that, as the cohesive strength of target sand is increased, the amount of sand free drawdown and depletion becomes larger and if the value of cohesive strength reaches 1370 psi, no sanding risk considered for this well at mentioned depth and was termed Idealized cohesive strength (ICS).Figure 4 show the values of initial and promoted cohesive strength (S oi and S op ) for prevention of sand production at current production condition as function of depth for Ahwaz sandstone member of Asmari formation.

Effect of uniaxial compressive strength
The behavior of sand free envelope plot for series of uniaxial compressive strength (UCS) from 1122 to 8000 psi at the depth of 8700 ft is shown in Figure 5.It can be seen that, As the UCS of target sand is increased, the amount of sand free drawdown and depletion becomes  strength (IUCS).Figure 6 show the values of initial and promoted uniaxial compressive strength (UCS i and UCS p ) for prevention of sand production at current production condition as function of depth for Ahwaz sandstone member.

Effect of well trajectory
Figure 7 shows the two dimensional plot of critical bottom hole flowing pressure of the wells with different inclination and azimuth in Ahvaz sandstone member.The vertical axis is the wellbore inclination in degree, horizontal axis is the azimuth in direction of maximum horizontal in degree and the colored regions indicate the CBHFP.It is concluded that in normal stress regime the CBHFP of a horizontal borehole is greater than the vertical borehole, so the vertical boreholes have less potential for sanding than the horizontal boreholes and almost all the deviated wells.In this case, the best drilling trajectory is a well with inclination of 30° and azimuth of 90° (Zare-Reisabadi et al., 2012).It is also obvious that, drilling parallel to the minimum horizontal stress direction is the best trajectory in this case.In addition, it shows that the CBHFP or sanding potential is highly sensitive to the inclination in all direction or azimuth.

Effect of horizontal stress anisotropy
Figure 8 show the effects of maximum horizontal stress and horizontal stress anisotropy ratio (HSAR) on CBHFP.
The horizontal stress anisotropy ratio is defined as: (20) As Figure 8 depicts, an increase in maximum horizontal stress leads to increase in the CBHFP in an invariable HR.It can be seen that increasing HSAR tends to increase in the CBHFP.

Effect of Biot's factor and poroelastic stress coefficient
Figure 9 indicates the influences of the Biot's factor and poroelastic stress coefficient on CBHFP.It can be   concluded that the CBHFP increases by increasing of the poroelastic stress coefficient in particular Biot's factor.Also in a fixed poroelastic stress coefficient an increase in Biot's factor increases the CBHFP.

Conclusion
This paper presents the effects of reservoir and geomechanical parameters including well trajectory, poroelastic stress coefficient, Biot's factor, maximum horizontal stress, horizontal stress anisotropy ratio, cohesive strength and uniaxial compressive strength (UCS) on sand production from openhole wells.The results indicated that: 1.In normal stress regime the critical bottom hole flowing pressure of a horizontal borehole is greater than the vertical borehole, so the vertical boreholes have less potential for sanding than the horizontal boreholes and almost all the deviated wells.In this case, the best drilling trajectory associated with inclination of 30° and azimuth of 90°.
2. The critical bottomhole flowing pressure increases with increasing the poroelastic stress coefficient, Biot's factor, maximum horizontal stress, horizontal stress anisotropy ratio, but decreases with increasing the cohesive strength and UCS of rocks.

INTRODUCTION
Water-flooding is a popular secondary recovery method which is responsible for most oil produced beyond the primary recovery mechanisms (Willhite, 1986).It is widely used because it is highly effective and water is readily available in huge quantities and is inexpensive.For this reason, water-flooding has commanded considerable attention from the petroleum industry in attempts to improvise the process of oil extraction from the subsurface.Water flooding is used in oil reservoirs to increase the rate of oil production and oil recovery.Modelling this process, especially in heterogeneous reservoirs such as the ones experience vertical stratification, is very tedious and cannot be done analytically since all governing forces cannot be included.In order to anticipate the performance of such process, we need a very good understanding of the subsurface in addition to accurately capture fluid behaviour under the reservoir conditions and rock fluid interaction properties.
With these has been said, huge cost of subsurface sampling (fluid and rock), well testing, 4D seismic acquisition, permeant downhole gauges for real time reservoir management, sophisticated data assimilation algorithms and powerful processors needed to develop good picture of reservoir and construct updated representative model which utilized to make decision of implementing such process as perfect secondary recovery techniques for system and manage it throughout the project.
The aim of this review is to predicate water flooding performance in stratified reservoirs using cheap qualitative correlation derived from data driven proxy model.While there are other predication techniques, data-driven proxy model is relatively modern and incorporates numerous parameters governing such techniques and leading to achievement of close and precise estimations of ultimate recovery factor.It also implicitly includes relevant theories and their applicability.

THEORETICAL FRAMEWORK
Since the discovery of water flooding technology in Pennsylvania in the early days of oil mining, several theories have been developed to explain the phenomenon.Each of the theories has their particular assumptions, limitations and drawbacks which need to be taken into account when used to predict the performance.In 1941, Leverett had presented in his paper the concept of fractional flow which became very important concept that has been used widely to estimate the performance of water-flooding.Buckley and Leverett (1942) were the first to develop a theory that was published 1942.Basically, Buckley and Leverett's theory is considered to be one dimension displacement model.To cope with the reservoir heterogeneity, Stiles (1949) published the first technique to predict the performance of the water flooding in stratified reservoirs.His method is solely valid for reservoirs with unity mobility ratio or very close to unity.The main assumption in his method is that the absolute permeability for each layer is the only parameter control for the velocity (the velocity is solely proportional to layers' absolute permeability) and the system is noncommunicating.Dykstra and Parson (1950) extended Stiles' work in 1951 in attempts at estimating the flood sweep efficiency for non-communicating stratified reservoirs for different mobility ratios.
This was the first work that published to explain the vertical stratification reservoirs with mobility ratio larger or Omar et al. 61 lower than unity.Additionally, they provided basic equations which were used to predict the performance and a method to handle stratification.Besides, they introduced their coefficient (VDP) to characterise the degree of permeability variation in the vertical direction.
In reality, it is rare to find such stratified reservoirs consisting of discrete layers that are completely isolated.Even though in most circumstances there is a thin impermeable layer which lies between the outer layers due to depositional process, it is usually discontinuous and there is spatial communication between layers.Such reservoirs are called communicating stratified reservoirs.Ignoring cross flow between layers is also assumed in the all models presented above.Hiatt (1958) was the author of the first paper predicting of water-flooding performance in the stratified reservoirs with complete cross flow.El-Khatib (1985) presented an extension to Hiatt's work by publishing a model which takes into account the variation in the other rock parameters other than the variation in permeability.He investigated the difference between the communicating and the non-communicating system and compared their impacts on the water-flooding performance.All the above mentioned theories deal with the horizontal reservoirs.Naturally, due to tectonic movements and depositional environments the oil traps, reservoirs usually exhibit an inclination from the horizontal.The reservoirs which dip from the horizontal with an angle cause the water to be injected in an up or down dip direction.El-Khatib (2010) presented a mathematical model to account for the effects of the gravitational forces in terms of dipping angle and the density difference between the driving and the driven fluids.El-Khatib's model (the model is demonstrated in the Figure 1) includes equations to predict the performance of the inclined reservoirs and their permeability variation and relies on log-normal distribution.He accounted for the inclination by including the effect of the density difference and the dipping in form of dimensionless number.All aforementioned models and other mathematical models failed to incorporate all the forces that are responsible for immiscibly of displaced hydrocarbons by water in production wells.Conventionally, the geological properties are not uniformly spatially distributed.For this reason, the models ignore either one or more forces to simplify the mathematical approach towards solving this mechanism.The only available tool that can model all the main forces is the reservoir simulation tool.The construction of the full model to simulate behaviour of the reservoir can be done by coupling laws of mass conservation, momentum and energy.Fundamentally, the two phase continuity flow equations for the oil and water, Equations 2 and 3 can be *Corresponding author.E-mail: Ahmed.omar@ucalgary.ca.
Author(s) agree that this article remain permanently open access under the terms of the Creative Commons Attribution License 4.0 International License solved to evaluate pressure distribution in a typical reservoir as function of time.The results can be coupled with Darcy laws to calculate phase flows and time relationship.Since there are four unknowns in the equations (So, Sw, Po and Pw), the simulator needs two more equations to couple them with flow equations for oil and water and solve those equations simultaneously for predicting the performance. (2) (3) Where kx, ky and kz are the absolute permeability in I, j and k respectively.are oil phase and water phase potentials respectively.Differential oil potentials are defined as presented in the equations below: From the above equations, oil pressure, water pressure and the saturations for water and oil need to be involved in the equations.Capillary pressure, (p c =p o -p w ) and total saturation (s o +s w =1) are also included in solving two phase differential equations (Hebb, 1949).
The simulator is used daily by reservoir engineers since it yields good estimations that simulate the dynamic behaviour of the reservoir.However, despite its prominence and importance the reservoir simulator still contains considerable degree of uncertainty due to geological heterogeneity and the parameters are not accurately measured.Moreover, reservoir parameters are not static and keep changing throughout life of the reservoir.As a result, this leads to ill model that has to be updated with respect to time by using data assimilation for tuning the model to reproduce the past (history matching process).Due to the large number of parameters in model, a huge number of combinations are included to tune the model which may accurately match the past but at same time forecast totally wrong future.Clearly the reservoir simulation model parameters are definitely not well defined or inexact as measured at the initial conditions.This parameterization issue is continuous with time as the physical properties of the system change with time not to mention the complex illposed inverse problem of the reservoir simulation (El-Khatib, 2010).
There are newer techniques applied to solve these problems like downscales or full data-driven black box models.Lately, they yield very promising result.Nowadays, neural network technique is increasingly implemented to handle reservoir behaviour complexity and nonlinearity.Data-driven proxy model is a neurocomputing technique of regression that was introduced by McCulloch and Pitts (1943).They replicated the way of processing information by human brain to perform some arithmetic operations.They presented the looselyconnected processors ability of sufficiently performing any logic operation with minimal percentage of error.The implementation of the modern neuro-computing can be traced back to 1960s.The neural network consists of the following (Braunschweig and Day, 1995): 1) Multiple processors of either independent or connected computers each possessing local memory or being implemented in software to process data independently 2) Unidirectional connections to connect the processors for transferring information either by physical wire or simulated by software.From the above equations, oil pressure, water pressure and the saturations for 3) Each processor has many inputs.As such, information from other processors is summed up to evaluate the mathematical equation that is used as a transformational function.The result of the transformation function is the only single output from the processor.4) Lastly, the learning rule to train the network either by adjusting the local memory of physical processor or modifying the connections strength to optimise it and find the best solution for the problem.
This model is also known as the linear threshold gate.As such, it consists of a set of inputs of a given output, say y.
The inputs have the following order I 1 , I 2 , I 3 , I 4 ,…,I m..This model classifies the set of inputs into two different categories and in so doing render the output y binary.It can be illustrated in Figure 2 as follows: The artificial neural network (ANN)-based models are applied widely in industry to model nonlinearity such as the identification of injector-producer relationship identification and flow units characterization between injection-production wells (Lee et al., 2008) .In the light of production optimization field, the ANN has been widely implementing in anticipating and optimizing the well performances in different kind of reservoirs with wide range of heterogeneity degrees (Bansal and Ertekin 2013).It is also used in correlating the reservoir fluid properties without physically relating to the independent variables to lower the cost of field sampling and expensive laboratory experimental procedures.Naseri et al. (2012) published ANN-based correlation to predict the viscosity of crude oil at different operation conditions and isothermal subsurface temperatures and it gave a promising result in estimating the viscosity of the large set of Iranian crudes.Nada et al. (2012) conducted a survey by collecting 104 real data sets to validate propagation of neural network to derive correlation that predicts bubble point pressure for Iraqi fields with average error percentage about 6.5%.This data processing techniques has recently been used as a hybrid with other artificial intelligent techniques such as fuzzy logic (FL) and genetic algorithm (GA) to enhance neuro-computing technique.Abdelrigeeb et al. (2014) have presented novel comprehensive study derived from hybrid system of a genetically optimized neural network (GA-ANN) and Neuro-Fuzzy (NF) modelling techniques to estimate bubble point pressure parameters and eventually they introduced an expert system.
AI techniques are used widely in improving oil recovery (IOR) as well.Ghoodjani and Bolouri (2011) correlated the final oil recovery factor as function of porosity, pore throat size, viscosity ratio for stratified reservoir at constant injection rate for pore scale water-flooding process.EOR processes have drawn too much attention from the industry lately due to costly and time consuming experiments to investigate fluid and rock behaviour under wide range of operating conditions.As a result, the minimum miscible pressure (MMP) is considered as the main parameter that governs the efficiency of miscible Omar et al. 63 gas projects, there are many predictive MMP correlations presented for miscible CO 2 flooding based on regression of laboratory data.Amin et al. (2013) used modern neurocomputing techniques called Least-Squares Support Vector Machine (LSSVM) for impure and pure CO 2 systems and their MMP result showed near perfect match with slim tube displacements and rising bubble apparatus (RBA) experiments result.
For the reservoir management and surveillance topic, a complex simulation model with acceptable level of representation is the cornerstone.The uncertainty associated with developing geological models, fluid behaviour and reservoir conditions affect the output significantly.In addition to data dimensionality and model complexity, the manner of processing the data leads to huge computational costs.The artificial intelligent techniques coupled with data mining are lately utilized to reduce the reservoir simulation's problems without causing under-modelling issues.They are collectively known as surrogates for full reservoir simulation model either on well level by reproducing only production and injection data or mimic the model dynamic parameters by allocation the pressure and saturation in all spatial discretization grids for every time step (Mohaghegh et al., 2012(Mohaghegh et al., , 2014)).The novelty is in the implementation of such techniques to systematically model the flooding mechanisms by correlating the properties affecting the main governing forces coupled with reservoir physical characterizations for rapidly and accurately predicting its feasibility.

Prologue
In this study, objectives were accomplished using CMG (WINPROP, IMEX and CMOST) which is a black oil three-dimension reservoir simulator.Firstly all massive simulation studies were run and presented here were designed to investigate the impact of the vertical heterogeneity, mobility ratio, dipping angle, layers permeability ordering, injection rate and density difference between the displacing and displaced fluids.The cases are aimed at capturing the effects of these parameters on the main three forces which govern the performance of the water-flooding in stratified reservoirs.After all parameters allocated, they are used as input to neural network model and optimize it to accurately duplicate the final recovery factor anticipated.

Simulation model description
The geological model (Figure 4) of reservoir is described by a model with dimensions 1500*1500*300 in the X*Y*Z dimensions respectively.The model is divided in terms of grid cells (reservoir gridding) to generate 15 cells in x direction, 15 cells in Y direction and 10 simulation cells in Z direction in order to get rid of all numerical dispersion problems.The reservoir is constructed from discrete layers.Each layer has different horizontal and vertical permeability.The permeability variation in vertical direction is characterized by Log-Normal distribution with VDP=0.4 demonstrated in Figure 3 and the geometric permeability average is  equal to 475 md.The reservoir model contains one well to inject the water and another one to produce the oil as shown in Figure 4.The injection well (I) was located in the cell (1,1) to inject the water up dip direction in the cases where the reservoir is inclined and the production well (P) was located at cell (15,15).The top of reservoir is located at 4300 ft in term of reservoir depth and the initial pressure at 4450 ft is set to be 4000 psi.The shallowest depth with 100% water saturation being at 5200 ft.The production was started simultaneously with the injection operation at (10 Sep 2015).The reservoir was depleted at maximum fluid production rate equal to 3500 bbl/day and the water was injected at 3900 bbl/day of water.These operation constrains were applied to keep depleting the reservoir above its saturation pressure.The model was running for 15 years.The reservoir porosity is set to be constant 0.3 throughout the reservoir.PVT properties for oil and the gas which used in the model are created using the winpop software and the hydrocarbon composition used are presented in Table 1.In term of interaction fluid-rock properties, the model of well sorted consolidated sandstone exponents=3 is used.

Sensitivity cases description
Seven simulation studies were conducted to examine the effects of vertical heterogeneity, the factors which control the gravity dimensionless number, cross flow between the adjacent layers, the layers permeability order and the mobility ratio on the performance of the water-flooding is such reservoirs which exhibited the stratification: 1) The first three cases were designed to examine the impact of reservoir dipping angles with including the effect of mobility ratio.All the rest of properties are set to be replica to the base case model.
In addition to horizontal model, the water-flooding in 30, 45, -30 and -45 (Figure 5 demonstrate the inclination cases) inclined models were simulated for mobility ratios 5, 1 and 0.5.
2) The case 4 was solely designed to investigate the impact of density different between displacing and displaced fluids.The simulation study was conducted on horizontal communicating models with the density of the crude oils are (40, 55 and 60 lb/ft 3 ).The mobility ratio set to be 0.5 for all studies.
3) The purpose of study number 5 is to examine the gravity cross flow which depend on the actual permeability ordering of the layers in the reservoir.The evaluation of the bouncy effect can be achieved by randomly assigning the layers permeability without any order.In addition to the base case and fining downward case models (all the cases are shown in Table 2).
4) The injection rate considered to be one of the main control factors and have major effect on the gravity forces.However increasing of the injection rate normally leads oil recovery factor to get higher.Virtually injection rate should not exceed the critical value.The values of injection rates used in this study are shown in Table 3 and all the cases were carried in horizontal model with mobility ratio equal to 0.5.5) The horizontal permeability variation in the vertical direction has significant effect on the performance of the water-flooding.To investigate the effect of the degree of the reservoir heterogeneity and since the permeability is represented by log-normal distribution.
The layers permeability was generated using log-normal distribution with VDP equal to 0.2 and 0.6 in addition to the base case.The horizontal base case model was used to examine the effect of the permeability variation degree.The permeability distribution functions are demonstrated in Figure 6.

ANALYSIS OF SIMULATION RESULT
There was a positive correlation between the recovery factor and the reservoir dipping angle.It was observed that an increase in recovery factor was positively related to the dipping angle of reservoirs.Figure 7 is a plot of gas recovery factor against time with different reservoir inclinations and a mobility ratio of 0.5.The increased recovery factor is attributed to the gravity forces particularly in this permeability depth configuration by withdraw the water from the top layers to adjacent ones that led to reduction of the water velocity in those layers by (ρg sinα).Furthermore, reshaping the flood front by gravity cross flow eventually sharpened the flood front leading to delayed water breakthrough and displacement of more crude from the rocks to wells.Oil recovery has an indirect proportional relationship with mobility ratio as can be noticed in Figure 8.As such, an increase in mobility ratio leads to a decrease in oil recovery and vice versa.Gravity decreases the water velocity since the up-dipping injection is more pronounced in the case with higher mobility ratio (M=5) than with a favourable mobility ratio (M=0.5). Figure 9 illustrates this argument.With higher mobility ratio, the gravity force sharpens the front and delays the water breakthrough.On the other hand, with lower mobility ratio, the gravity has minor impact since the flood front is already sharp and recovered oil is high.Typically, the oil recovery increases as the mobility ratio decreases as clearly demonstrated in Figure 10.
The increase in density difference between the displaced and displacing fluid enhances the oil recovery, Figure 11.The term which represents the density difference appears in the fractional flow equation for inclined reservoirs, Equation 7. As this difference increases, it lowers the water cut (Figure 12) and increases the oil recovery. (7) Basically, the force of gravity increases water saturation in the bottom layers through withdraws of water from the top layers towards the bottom.The domination of the gravity force is dependent on permeability-depth configuration and Kv/Kh ratio.In the case where the increases the oil recovery.
Basically, the force of gravity        permeability decreases with depth, the front advances slower in the top layers due to downward drainage.This effect sharpens the front and enhances the performance as illustrated in Figure 13.Additionally, the cross flow between the layers due to the bouncy effect is the greatest in the fining downward case.The bouncy drives the oil to flow upward while the water segregates downward which delays the front advance.The magnitude of the bouncy effect on the performance significantly relies on the area which opens to the gravity segregation.The extension of this area depends on the permeability contrast between the adjacent layers.For the opposite permeability configuration, the performance gets worse since the gravity effect improves the water velocity in bottom layers which accelerates the water breakthrough.In principle, the gravity force improves the recovery as permeability decreases with depth and vies versa.In the case where the permeability is randomly distributed with depth (that is, permeability alternately increase and decrease with depth), the incremental oil recovery as permeability decreases with depth tends to be compensated with the reduction in the recovery brought about by the opposite permeability configuration.
As the result of this compensation, the random permeability distribution tends to have minor effect on the performance.The fining upward case is the worst case since it accelerates the water breakthrough and increase the water cut (Figure 14).
The performance of the flood worsens as the Dykstra-Parson coefficient increases because the reservoir became more heterogeneous as the permeability variation increases (Figure 15).
Lastly the, injection rate sensitivity, the recovery factors improve as the water injection rate increases and it is summarized in Table 3 and demonstrated in Figure 16.The optimum recovery factor is obtained at injection rate 10000 b/d.Principally, similar displacement processes depended on the viscous and gravity forces, the higher injection rates such 12000 and above, tend to accelerate the water velocity in layers with high absolute permeability.This rapid movement leads to reduction of the gravitational cross flow ending up with earlier breakthrough and poor performance.The critical value of surface rate for this system is (3200 b/d).

NEURAL NETWORK DESCRIPTION
Based on the above sensitivity study, the physical properties that affect the cross flow due to viscous and gravity forces were randomly changed by CMOST software with Latin Hypercube method used as the sampling method to generate more than 3000 solved examples.Those solved cases were supplied to the network for training and being used the neural network input layer (Table 4).Regarding training algorithms, back   The statistical summary for polynomial regression is shown in Table 5.Even for the deployment data set, The anticipated values by artificial neural network agreed with the simulations results (Figure 21), on the other hand the polynomial regression result fluctuate with good modelling at some ranges to very bad matching at others with high variation domain (Figure 22).For this reason, it gives the proxy model a basis for predicting the oil recovery factor for reservoirs implementing water flooding as enhanced oil recovery method just by plugging in the reservoir model description and injection operations parameters.
Although this data driven model gives a good first estimation of the recovery factory it under-models all the reality in term of the physics of the reservoir and not take into account any parameter outside the training data set or beyond its quality which is solely dependent on quality of the training.In the course of the accuracy, it cannot accurately estimate any physical phenomena especially complex unforeseen parameters such as water break through time.All the models' parameters are meaningless from the physics point of view.

Conclusion
In summary, predicating water flooding performance using the data-driven proxy model provides more reliable results than the traditional models that were proposed earlier.The main disadvantage of the correlation is their failure to predict correctly spatial and temporal pressure    and fluid saturations as well as other reservoir phenomena's.
In this survey, artificial neural networks model was used and the main advantage is that it incorporates most of the affecting parameters that can be easily measured at the sites.In fact, regression by ANN has extremely low error margins as compared to statistical regression as shown in figure 18 and 19.From Figure 21, it can be seen that simulation results significantly agreed with those obtained via ANN.However, it does not take into account parameters not included in the training set which is a major shortcoming.Additionally, ANN cannot estimate particularly unanticipated reservoir parameters such as water breakthrough.In the light of the high degree of agreement between simulation results and ANN, it can be concluded that the precision, accuracy and capacity of ANN in predicating water flooding performance is reliably dependable.

Future studies
1. Including the effect of the cross flow between the strata due to capillary pressure forces which dominated in specific kind of the reservoirs where there is high capillary pressure contrast in the adjacent layers.In such reservoirs, the system tries to establish the capillary equilibrium by drawing the water from the loose layers to the tighten ones.2. Enhancing the performance of the network by extending the training parameters to include others like (irregularity of the reservoir area, continuous heterogeneity in X,Y and Z directions, number of injections and producing wells and flood pattern).

Figure 2 .
Figure 2. Pore pressure and in-situ stresses profiles of Asmari formation.

Figure 3 .
Figure 3. Sand free operating envelope plots for different values of cohesive strength.

Figure 6 .
Figure 6.Prediction of UCS for Ahvaz sandstone member.

Figure 7 .
Figure 7. CBHFP for various wellbore trajectories in normal stress regime.

Figure 8 .
Figure 8. Effects of maximum horizontal stress and horizontal stress anisotropy on CBHFP.

Figure 14 .
Figure 14.Water cut vs Time for model 0° inclination with different permeability ordering cases.

Figure 21 .
Figure 21.Simulation and ANN regression values.

Figure 22 .
Figure 22.Simulation and normal regression values.

Table 1 .
Geomechanical and reservoir properties of target depth.