GRAINSCALE SIMULATION OF MULTIPHASE FLUID FLOW IN POROUS

The computational simulation of petroleum production systems gained great importance in the last years. This is due to the interest of the industry in developing new methods of secondary oil recovery, which will allow an increase in the recovery rates at reasonable costs. Within this context, the study of the multiphase flow through porous media, which allows a better comprehension of the physical phenomena inside a complex medium such as a reservoir rock, is of great interest. The present study uses the Moving Particle SemiImplicit (MPS) method to model the multiphase flow through a porous medium. MPS is a fully lagrangean, meshless method that discretizes the simulation domain in particles. It solves the mass and moment conservation equations by using the differential operators derived from a particle interaction model.As a study case, the flow of two immiscible fluids – oil and water – through a sandstone sample was simulated. The geometry of the sandstone sample was obtained by a high resolution CT scan of a real porous rock and has a volume of 1 mm3. Various situations of wettability were simulated to study its influence on oil recovery in the sample. This was accomplished by varying the contact angle between the fluids and the solid matrix of the rock. It was verified that higher recovery rates are obtained when water wets the surface of the grains with more strength than the oil.


INTRODUCTION
The modeling of fluid flow through porous media is of great interest for the petroleum industry, because it allows the comprehension of how the interfacial physics affects the behavior of the fluids and its macroscopic properties.Recently, some studies on numerical approaches to this kind of flow have been reported.The water-oil flow was simulated in a 2D porous media using the Smooth Particle Hydrodynamics (SPH) method in [10] and it was reproduced using a 3D model in [3].Despite the use of a hypothetical porous media with rela-tively high porosity, their results provide many insights about the relevance of surface tension on the effectiveness of oil recovery.
On the other hand, [6] showed the simulation of the single phase flow of water in a real sandstone sample using the Moving Particle Semi-Implicit (MPS) method and calculated the absolute permeability of the sample.The results found were very high for a sandstone sample, but according with the authors this is expected to happen in samples of very small size.
The main goals of this study is to show an implementation of the MPS method to simulate the multiphase fluid flow in porous media and its application to model a system which represents an oil reservoir in order to study the behavior of the multiphase flow during production phase.It will be studied the influence of wetting angles on the productive capacity of the reservoir.

COMPUTATIONAL METHOD
The computational method used in this study is the Moving Particle Semi-Implicit (MPS), which was proposed in [5].This method was initially developed to simulate the dynamic behavior of an incompressible fluid.It is a particle based method and utilizes the Lagrangean description to formulate its governing equations: continuity equation and Navier-Stokes: where ρ is fluid density, u  is the velocity vector, P is pressure, ν is kinematic viscosi- ty, m is particle mass, is surface tension force and f  accounts for other accelerations.All differential operators are derived from a particle interaction model based on the weight function which depends on the relative position of the particles: where r is the distance between two particles and r e is the effective radius, which limits the region where interaction between particles occurs.
For three-dimensional cases analyzed herein, r e was set to 2.1l 0 , where l 0 is the initial distance between particles, to calculate pressure gradient and the particle number density.r e is 4.0l 0 to cases involving the laplacian operator.
To assure the incompressibility condition of the flow, MPS adapts the concept of particle number density (pnd).The pnd of particle i, considering its neighboring particles j, is calculated by the following equation: The incompressibility condition is met when the value of pnd is constant.In other words, pnd must be equal to its initial value (pnd 0 ) at any given moment.
Given a scalar function , the gradient vector and the laplacian of a particle , considering its neighboring particles , are represented by Eqs. ( 5) and (6).The gradient vector used in this study is the formulation proposed in [7].
where is the number of spatial dimensions and is a parameter that represents the growth of variance and may be calculated using the following equation: The MPS method is based on a semi-implicit algorithm where the pressure is calculated implicitly and all other terms such as gravity and viscosity are calculated explicitly.The Poisson's equation can be deducted from the continuity equation and the pressure gradient: where, is the particle number density calculated explicitly.The term of the left can be dicretized using the laplacian model, leading to a system of linear equations.
To calculate the surface tension forces, an approach based on potential force between particles was used.It is an analogy to the inter-molecular forces that produce this effect.We adopted the formulation proposed in [4], in which the potential force is given by: where C is a coefficient that defines the intensity of the potential and can be calculated by the following equation: where σ is the surface tension coefficient and p(r) is a function that defines the shape of the potential and is given by: As a result, the attraction or repulsion force between two particles is given by a very simple expresion: The direction of this force is the same of the vector that connects two particles.The computational calculation of the potential force given by Eq. ( 12) may present some instabilities.To correct this, a smoothing term is added: where C smooth is a weighing constant.Following [4], C smooth = 1.0 was used in the present study.
To calculate surface tension between a solid and a liquid particle, the following equation is derived from Young's equation: where C fs is the coefficient of the potential between fluid-solid and θ is the contact angle in the liquid-solid interface.To calculate the interaction force between a solid and a fluid particle, we substitute C in Eq.( 12) by C fs .
Particles are considered as free surface ones when their pnd is smaller than ( ).According to [5], may vary between 0.80 and 0.99.In the present study, we adopted β=0.97.
Rigid walls are discretized in particles in which there's no movement.Two layers of dummy particles are placed in the side where there's no contact with the fluids.Dummies particles exist only to ensure that the pnd of particles near the walls are correctly calculated.
As inflow boundary condition, fluid particles are injected to the domain by moving wall particles with their dummies at a constant velocity.When they displace a distance equal or larger than the particle distance, these wall particles are replaced by fluid particles.At the same time, the position of wall and their dummies are shifted back one particle distance, restarting the loop.

Sandstone sample
The geometry that was used in our simulations was obtained by a high resolution CT scan of reservoir sandstone.The dimensions of the sample are 1.0 x 1.0 x 1.0 mm and its absolute porosity is 18.20%.The CT scan generated a model with 300 x 300 x 300 voxels.Figure 1 shows the pores of the sandstone sample.To convert the information given by the CT scan into the input file for our MPS simulation program, a pre-processing program was created.The following subsection describes its algorithm.

Pre-processing
The original sample geometry contains 27 million voxels, which must be converted into a list of particles to be used in MPS simulations.Due to the limitations of memory size and processing time, a resolution reduction was needed.We reduced the original sample resolution, which was 300 pts/mm, to 50 pts/mm.The method employed to do this was to evaluate which kind of particle was the most common in a cube of 6 x 6 x 6 points.The resulting particle distance is 2x10 -5 m.
The original CT scan of the sample identifies only two kinds of region: solid and void.The void points were converted into fluid particles and the solid points were converted into two kinds of solid particles: the ones that have direct contact with fluid are defined as wall particles, while the ones that are in the inner region of the solid are defined as dummies.Moreover, dummies that are placed more than two particle distances from a wall particle are removed in order to improve the efficiency of memory usage.
To assure unidirectional flow, we added walls to the lateral faces of the sample.To simulate fluid injection, an inflow boundary condition was added above the sample, with a 0.4 mm transition region between the inflow boundary and the entrance section of the sample.We extruded the bottom section of the sample to form a 0.2mm transition region at the base of the model to allow the correct measurement of pressure at the exit of the sample.
Figure 2 shows a section of the model used in the simulations.It contains a total of 266157 particles and its absolute porosity is 18.17%,only 0.16% less than the original geometry.
Figure 2. Multiphase flow of water and oil in a porous medium where oil is the wetting phase.

Fluid properties
Simulated fluids were water and oil, and their physical properties are described in Table 1.Due to the huge variety of kinds of oil found around the world, it was decided to use reference values.The chosen oil was the West Texas Intermediate, which is used as benchmark oil in the North-American market.Its properties were extracted from [11], at 15ºC and 0% weathering.According to [1], the values of contact angles of an oil/water/rock system vary greatly along the reservoir.In this study, we considered three distinct situations of wettability, in which the water is the wetting phase, the oil is the wetting phase, and a third situation, where both of the fluids wet the rock with the same strength.Table 2 shows the contact angles of water and oil used in each situation.

Simulation of fluid flow
Two different situations of fluid flow were simulated to verify the influence of the contact angles on the capacity of a fluid to displace the other: water injection in an oil saturated sample and oil injection in a water saturated sample.Altogether, we carried out simulations of six different situations, described in Table 3.The simulations were carried out with a constant fluid inflow of 0.0125 m/s and a time step of 10 -6 s.The single phase flow simulations were used to calculate absolute permeability of the porous media.Considering that the steady flow is reached very fast in these cases, 0.05s of simulation were enough to obtain results, while the cases of multiphase flow were simulated for 0.25s.

Pressure drop in multiphase flow
During the multiphase flow simulations, pressure drop was measured over time for both situations of injection (water or oil), each one with the three wettability types.To be sure that the displaced fluid reached its irreducible saturation, we measured the flow for each fluid at the exit of the sample.Figures 3 and 4 show pressure drop and flow measurements for the situations with injection of water in an oil saturated sample and injection of oil in a water saturated sample, respectively.Table 4 shows the measured pressure drop for each situation.Both figures show that the shape of pressure drop is very similar when the same fluid is injected.Although a slight increase in pressure drop can be seen when oil is the wetting phase, the influence of the contact angles on the pressure drop is not very significant.

Irreducible saturation
To obtain the irreducible saturation of a fluid, we considered a sample saturated with a fluid and injected another fluid in the sample until the first fluid stopped being produced.The quantity of each fluid remaining inside the pores of the sample was measured.

Oil irreducible saturation
Water was injected in an oil saturated sample until irreducible saturation was achieved.After that, the quantity of particles of each fluid present in the pores of the sample is calculated, and the results are shown on Table 5.In Figure 5 we can see the distribution of the particles in the pores of the sample.It is noted that wettability have a significant influence over the capacity of water displacing oil.In Figure 5, oil, which completely fills the sample in the beginning of the simulation, is represented by grey particles.The distribution of waterrepresented by blue particlesgives a good insight on how the flow advances through the complex network of pores.Even after steady flow is reached, it is possible to observe that some regions still contain the fluid that filled them at the beginning of the simulation, meaning that they are either a pore without communication or a dead-end.
Regions nearby the edges of the sample also contain portions of the original fluid entrapped in them, clearly showing the blocking effects of the walls that were added to assure unidirectional flow.To minimize this effect, it is advised to use a larger sample.
In the situation where oil wets the rock it is possible to observe that water flows through oil, and a substantial portion of oil remains adhered to the grains of the rock.When water is the wetting phase, it displaces a larger portion of oil.It is also noted that a portion of oil is found in isolated pores, without communication.As the sandstone sample is quite small and considering that this oil cannot be removed from the sample, the results might overestimate the measurement of the irreducible saturation.

Water irreducible saturation
To measure water irreducible saturation, we injected oil in a water saturated sample.When water ceased to flow, the particles of each fluid present inside the sample were counted, and the results are listed in Table 6 and illustrated in Figure 6.Comparing these results with the ones obtained in the previous subsection, it is concluded that oil displaces water with much more ease than the opposite situation.This is probably due to the fact that oil has a much higher viscosity than water, giving it a much lower mobility.This way, it is easier for the oil to displace the water instead of flowing through it.

CONCLUSIONS
In the present study, we developed an implementation of the Moving Particle Semi-Implicit (MPS) method capable of simulating the multiphase fluid flow through porous media.As case study, it was used a sandstone sample in which both oil and water were injected.Moreover, irreducible saturations for each fluid and the absolute permeability of the rock were calculated, and the multiphase flow was studied over time.
In multiphase flow simulations, it is verified that oil displaces water with much more ease than water displaces oil.This phenomenon is attributed to the difference of viscosity between the fluids, which makes the oil much less mobile.According to the numerical results, the contact angles of the oil/water/rock system did not influence significantly pressure drop through the sample.Considering that we employed a sample with relatively coarse resolution, it may be possible that the particle distance was too large to effectively represent the surface tension effect.
As future work, we suggest the simulation of simultaneous injection of water and oil in the sample.Considering that the MPS method allows the measurement of flow for each fluid separately, the relative permeabilities can be easily calculated.Moreover, in order to allow a more accurate representation of the fluid flow, it is suggested to utilize samples with higher resolutions than the one that we employed.Also, we suggest the employment of a larger porous rock sample in order to compare our results with the ones obtained by measurements using real samples.

ACKNOWLEDGMENTS
The author would like to express their gratitude to PETROBRAS for the financial support to the research.

Figure 1 .
Figure 1.Visualization of the pores of the sandstone sample.

Figure 3 .
Figure 3. Pressure drop and flow measurements over time for the situation in which an oil saturated sample was injected with water.

Figure 4 .
Figure 4. Pressure drop and flow measurements over time for the situation in which an water saturated sample was injected with oil.

Figure 5 .
Figure 5. Fluid particle distribution inside the sample after oil irreducible saturation was attained.In situation (a) oil is the wetting phase, in (b) there is no wetting phase and in (c) water is the wetting phase.

Figure 6 .
Figure 6.Fluid particle distribution inside the sample after water irreducible saturation was attained.In situation (a) oil is the wetting phase, in (b) there is no wetting phase and in (c) water is the wetting phase.

Table 1 .
Fluid properties used in the simulations.

Table 2 .
Wettability situations addressed in this study.

Table 3 .
multiphase fluid flow situations simulated in this study.

Table 4 .
Pressure drop measured in each situation of multiphase fluid flow.

Table 5 .
Oil irreducible saturation in the oil saturated sample after injection of water for each situation.

Table 6 .
Water irreducible saturation in the sample for each situation.