SIMULACIÓN NUMÉRICA DEL FLUJO BIFÁSICO AGUA-PETRÓLEO EN UN MEDIO POROSO
Johana Lizeth Pinilla Velandia1
1. Escuela de Ingeniería de Petróleos, Universidad Industrial de Santander, UIS, Carrera 27 calle 9, Bucaramanga, Colombia.
Recepción: 9 de Octubre de 2013 Aceptación: 19 de Diciembre de 2013
RESUMEN
Se presenta el modelamiento y simulación del flujo bifásico agua-aceite en un medio poroso homogéneo. Modelamientos de este tipo son necesarios para ayudar a comprender ciertas técnicas como la recuperación secundaria de hidrocarburos. Los dos fluidos considerados son inmiscibles e incompresibles. Se utiliza la ley de Darcy generalizada para el modelamiento del flujo de los fluidos. La resolución numérica se basa en una formulación IMPES del modelo original. La ecuación de la presión se resuelve de manera implícita y la ecuación de la saturación de manera explícita. La discretización espacial es hecha utilizando el método de volúmenes finitos. Las simulaciones numéricas se llevan a cabo en un medio poroso en tres dimensiones, considerando siempre un pozo de inyección y un pozo de producción. Los resultados permiten poner en evidencia los efectos de la presión de inyección, la presión capilar y la difusión.
Palabras clave: Ley de Darcy, flujo bifásico, volúmenes finitos, IMPES, medio poroso, recuperación de petróleo, simulación numérica.
Numerical simulation of oil-water twophase
flow in porous media
ABSTRAC
This work is motivated by the need to better understand the comprehension of two-phase flow water-oil, in the secondary recovery of oil. Two-phase immiscible and incompressible fluids are considered. The generalized Darcy´s is used for modeling fluid flow. For the numerical resolution the separation of the calculation of the pressure and saturation is required. An IMPES formulation of the original model is performed. The pressure equation is solved implicitly and the saturation equation explicitly. The spatial discretization is made by the finite volume method. Numerical simulations are carried out in a three-dimensional porous medium by taking into account injection and production wells. Numerical results highlight the effects of injection pressure, capillary pressure and diffusion.
Keywords: Darcy law, two-phase flow, finite volume, IMPES, porous media, oil recovery, numerical simulation.
1. INTRODUCCIÓN
La comprensión de los mecanismos del flujo bifásico de fluidos en medios porosos es de gran interés en la industria petrolera. Por ejemplo, las técnicas de recuperación secundaria y terciaria implican el desplazamiento del hidrocarburo por otro fluido, como agua en el primer caso o un polímero en el segundo caso. El estudio del flujo de fluidos en medios porosos puede realizarse a tres escalas diferentes [1]: escala del poro, escala de Darcy y gran escala. A la escala del poro, la longitud característica es del orden del diámetro de los poros (1µm a 100 µm); a esta escala, la velocidad del flujo es muy pequeña y la aproximación de Stokes es generalmente válida para modelar el flujo de los fluidos [2]. A la escala de Darcy, el medio poroso es tratado como un medio continuo cuya dimensión característica es del orden de algunos centímetros y el modelo de Darcy es ampliamente utilizado para describir el flujo a esta escala. En el modelamiento a gran escala o escala de cuenca, la longitud característica es del tamaño de las heterogeneidades y el modelamiento se realiza mediante ecuaciones obtenidas a partir de técnicas de cambio de escala, tales como el método de promedio volumétrico [3]. En este trabajo nos hemos interesado en el modelamiento y simulación, a la escala de Darcy, del desplazamiento del petróleo, provocado por la inyección de agua que permite recuperar el hidrocarburo en un pozo de producción. La presentación del artículo está estructurada de la siguiente manera. La sección 2 está dedicada a la presentación de las ecuaciones que intervienen en el modelamiento físico. Luego, en la sección 3 se describe el modelamiento matemático haciendo énfasis en el método IMPES (IMplícito en Presión, Explícito en Saturación). En la sección 4 se presenta de manera detallada la resolución numérica del sistema de ecuaciones y finalmente, en la sección 5 se muestran los resultados de simulaciones numéricas en un medio poroso en tres dimensiones.
2. MODELAMIENTO FÍSICO
En el contexto de la mecánica de fluidos, un medio
poroso puede ser definido como un material en el que
existen cavidades interconectadas por canales, en los
cuales fluyen los fluidos; estas cavidades son llamadas
poros [4]. Consideremos un medio poroso homogéneo e
incompresible, donde se encuentran presentes dos fases:
una fase agua denotada w y otra fase aceite, denotada
o. Sea este medio poroso y
, la frontera
de
. Las dos fases consideradas son no miscibles.
Denotemos por
, la variable
espacio y
, la variable tiempo. Las ecuaciones
que gobiernan el flujo de estos fluidos se describen a
continuación [5], [6], [7].
2.1 Conservación de la cantidad de movimiento
La ecuación de conservación de la cantidad de movimiento se traduce por la ley de Darcy generalizada [8], [9], [10].
Donde es la velocidad de filtración de
la fase
es el tensor de permeabilidad efectiva
de la fase
es la viscosidad de la fase
es la presión de la fase
es la densidad de la fase
y es
constante ya que los fluidos son incompresibles y g es el
vector aceleración de la gravedad. Además se tiene que
la permeabilidades efectiva y relativa se relacionan de la
siguiente manera:
donde K es el tensor de permeabilidad intrínseca del
medio poroso y es la permeabilidad
relativa de la fase y depende de la saturación del agua.
2.2 Ecuación de conservación de la masa
Teniendo en cuenta que el medio poroso es incompresible y el flujo de los fluidos es también incompresible, la ecuación de conservación de la masa para cada fase se escribe:
Donde es constante y corresponde a la porosidad del
medio poroso,
es la saturación de la fase α y qα es un
término relacionado con el flujo del fluido inyectado.
2.3 Ecuación de continuidad
La saturación de un fluido se define como la fracción del volumen de los poros de la roca ocupados por el fluido [7]. Asumiendo que el medio está completamente saturado, tenemos
donde Sw y So son la saturación de agua y aceite, respectivamente.
2.4 Ecuación de la presión capilar
Donde Pc(Sw) es la presión capilar, Pw la presión del agua y Po es la presión del aceite. Empíricamente, se tiene que Pc es una función de la saturación Sw.
2.5 Modelo físico final
Resumiendo, se tiene que el sistema de ecuaciones que gobierna el flujo bifásico incompresible en un medio poroso e incompresible es:
El sistema (5) es un sistema de diez ecuaciones y diez incógnitas. Las incógnitas son: las componentes de las velocidades para cada fase Vw = (uw, vw, ww), Vo = (uo, vo, wo), las saturaciones Sw y So y las presiones Pw y Po.
3. MODELAMIENTO MATEMÁTICO
Notemos que el sistema (5) es un sistema acoplado de ecuaciones diferenciales parciales no lineales, dependientes del tiempo. Para tratar este problema y separar el cálculo de la presión del cálculo de la saturación, utilizamos el método IMPES [11] [12], ampliamente utilizado por su facilidad de implementación y bajo costo computacional respecto a otros métodos, como el de solución simultánea [13]. Una revisión de los métodos existentes puede encontrarse en [14].
3.1 El método IMPES
La idea principal del método IMPES es separar el cálculo de la presión y el de la saturación. La aproximación en tiempo de la ecuación de la presión y de la saturación se realiza mediante un esquema implícito y explícito, respectivamente. Discusiones detalladas acerca de este método pueden encontrarse en [12] y una revisión de mejoramientos al algoritmo original puede encontrarse en [15] y [8]. A continuación se detalla la formulación IMPES del sistema (5). Empecemos por introducir las siguientes notaciones:
Modalidad de fase α:
Modalidad total:
Fracción de flujo de la fase α:
Velocidad Total:
Ahora se escribe (2) para cada fase:
Donde . Se hace la suma de estas dos últimas ecuaciones y luego, reemplazando
(3) en la ecuación obtenida, se tiene que:
Para expresar la velocidad total VT en función de Po y Sw, se escribe (1) para cada fase:
Ahora se hace la suma de estas dos últimas ecuaciones y utilizando (4), se obtiene la ecuación de velocidad total:
Entonces las velocidades para cada fase se pueden escribir de la siguiente manera:
Se reemplaza (7) en (6) y se obtiene la ecuación de la presión:
Escribiendo (2) para α= w y reemplazando Vw por (8), se obtiene la ecuación de la saturación:
Finalmente, las ecuaciones (7), (10) y (11) forman un sistema de ecuaciones equivalentes al sistema (5) llamado formulación IMPES. Ignorando los efectos de la gravedad, dicha formulación se escribe:
Esta nueva formulación permite escribir un sistema de cinco ecuaciones y cinco incógnitas. Las variables desconocidas son: las tres componentes de la velocidad total VT = (u, v, w), la saturación del agua Sw y la presión del aceite Po. A partir de ahora, llamaremos ecuación de la velocidad, de la presión y de la saturación a la primera, segunda y tercera ecuación del sistema (12), respectivamente.
3.2 Condiciones límites
Sea Ω un dominio acotado en . Sea la frontera de
Ω, de normal exterior n y [0, T] el intervalo de tiempo
de estudio del fenómeno físico. La frontera Γ está
particionada así:
Denotamos por Γe la parte de la donde se inyecta el agua, Γs la parte de Γ donde se recupera el aceite, Γi la parte impermeable de Γ y Γp la parte de donde se utilizan condiciones límites de periodicidad. Según la naturaleza de las fronteras consideradas, asociamos las siguientes condiciones límites:
• En el borde de inyección Γe podríamos
utilizar condiciones límites de velocidad donde
corresponde al caudal de inyección. O bien utilizar
condiciones límites sobre la presión,
,
donde Pe es la presión que se impone en el borde
Γe.
• En el borde impermeable Γi:
• En el borde de salida Γs, consideramos que no hay
flujo de presión capilar. Al igual que en borde de
inyección, las condiciones límites de salida pueden
ser sobre la velocidad o sobre la presión. Para la
velocidad, , y corresponde
al flujo total de los fluidos recuperados. O bien para
la presión
.
3.3 Condiciones iniciales
Al instante t=0 la saturación del agua es conocida. Entonces, se tiene:
4. SOLUCIÓN NUMÉRICA
La estrategia de aproximación utilizada se basa en los trabajos presentados en [16] y [17].
4.1 Algoritmo
Sea [0,T], T > 0, el intervalo de tiempo del fenómeno
a estudiar. Sea 0 = t0 < t1 < ... < tN = T una
partición de este intervalo, para N un entero positivo. Al
instante es conocida en todo el dominio
Ω. Entonces, tenemos para n = 0,1...N:
1. Se resuelve implícitamente la ecuación de la presión
del sistema (12) para encontrar la presión .
2. Se calcula con la ecuación de la velocidad del
sistema (12).
3. Con , y conocidos, se calcula
a
partir de la ecuación de la saturación del sistema(12).
4.2 discretización
4.2.1 notaciones y mallado
La discretización espacial es hecha mediante el
método de volúmenes finitos. El dominio Ω es el
paralelepípedo , con Lx, Ly
y Lz constantes positivas (ver figura 1). Se utiliza un
mallado cartesiano rectangular, compuesto de mallas
paralelepipédicas. Sean Nx, Ny y Nz, el número de nodos
en la dirección x, y, z respectivamente. El paso del
mallado espacial en cada una de estas direcciones está
definido por:
Los puntos (xi, yj, zk) son los centros de los nodos Qi,j,k, definidos por (ver figura 2):
El borde Qi,j,k, denotado ∂Qi,j,k,está dado por:
4.2.2 localización de las variables desconocidas
La presión del aceite, Po, y la saturación del agua, Sw, se localizan en el centro de las mallas. Las componentes de la velocidad total VT = (u, v, w) se toman en el centro de las caras de las mallas (ver figura 3).
Para simplificar las notaciones, eliminaremos los subíndices de la presión, la saturación y la velocidad total de las siguiente manera Po= P, Sw= Set VT= V.
4.3 discretización de la ecuación de la presión
Sea . Integrando la ecuación de la
presión del sistema (12) sobre Qi,j,k, tenemos:
Utilizando la fórmula de la divergencia obtenemos
Donde n∂Qi,j,k, es la normal exterior a la superficie ∂Qi,j,k,. Desarrollando para el lado izquierdo de la ecuación (13) obtenemos:
De manera similar, se desarrolla el lado derecho de (13). Las derivadas en cada dirección son discretizadas por diferencias finitas centradas, esto nos permite calcular P en el centro de las mallas. Se obtiene así el siguiente esquema:
Las movilidades, M, se aproximan mediante una media armónica en la interfaz entre las mallas [17]:
La resolución del sistema de ecuaciones obtenido es llevada a cabo mediante método iterativo de gradiente conjungado [18].
4.4 Discretización de la ecuación de la Velocidad Total
Una vez conocida la presión, se calcula cada una de las componentes de la velocidad en el centro de las caras de las mallas:
4.5 Discretización de la ecuación de la saturación
Se desarrolla un algoritmo de dos pasos. Primero se resuelve una ecuación de tipo hiperbólico no lineal y luego una ecuación de difusión. Para Sn dado, buscamos Sn+ ½, solución de:
Conociendo la solución de (14), buscamos Sn+1, solución de:
Para resolver la ecuación (14), procedemos de la siguiente manera. Primero, integramos la ecuación en cada uno de los nodos Qi,j,k:
Utilizando la aproximación de segundo orden:
donde , se obtiene la siguiente
aproximación para el primer término de la integral:
En cuanto al segundo término se
utiliza un esquema numérico de Murman [19], [17]. Este
esquema de primer orden se adapta bien a la resolución de
ecuaciones hiperbólicas no lineales como la ecuación (14).
Ahora vamos a concentrarnos en la segunda etapa del método de pasos fraccionarios, es decir, en la ecuación de difusión (15). Integrando la ecuación (15) sobre Qi,j,k:
Utilizando la fórmula de la divergencia y aproximando las derivadas obtenidas mediante diferencias finitas centradas, se obtiene el siguiente esquema:
El esquema (17) implica la resolución de un sistema de ecuaciones lineales. Para la resolución de este sistema se utiliza el método de gradiente conjugado.
Con el fin de verificar el orden de aproximación del esquema numérico descrito para el cálculo de la saturación, se calcula el error de aproximación E (S) para diferentes mallados. Se utiliza la norma L1 del error:
donde Sh es la solución numérica calculada y S la solución analítica dada por:
El cuadro 2 resume los errores obtenidos para cuatro mallados diferentes. Estas simulaciones permiten constatar que el método utilizado aproxima la solución al orden 1.
5. SIMULACIONES NUMÉRICAS
5.1 Simulación en un medio poroso tridimensional
La primera simulación numérica es llevada a cabo en el paralelepípedo Ω = [0.1 m] x [0.1 m] x [0.1 m] y las condiciones límites están definidas así:
Los parámetros físicos son los siguientes:
Y las condiciones iniciales son:
. Las figuras (4), (5) y (6), muestran respectivamente los resultados de la
evolución de la presión. En la figura (4) vemos como
el agua va llenando poco a poco el medio poroso desde
el sitio de inyección y va desplazando el aceite hacia el
sitio que podríamos llamar un pozo de producción. La
figura (5) nos permite verificar que la dirección del flujo
de los fluidos va desde el pozo de inyección hasta el de
producción y la figura (6) nos deja ver que la presión
dentro del medio poroso, evoluciona en el tiempo.
5.2 Comparación de los efectos de la presión capilar
Con el fin de observar los efectos de la presión capilar, realizamos ahora dos tipos de simulaciones. En la primera simulación, suponemos que el flujo de los fluidos es efectuado en ausencia de efectos capilares y en la segunda, tenemos en cuenta dichos efectos. Los parámetros de las simulaciones son:
Donde fueron definidos en
la sección anterior.
Los parámetros físicos son los siguientes:
Las condiciones iniciales son:
. Y la presión capilar, que es función de la saturación, está dada por Pc (S) = 15000 (S – 1).
La figura (7) muestra el frente de saturación para estas dos simulaciones. Puede observarse que la presión capilar tiende a extender el frente de saturación.
5.3 Comparación de los efectos de la presión de inyección
Para comparar los efectos de la presión, realizamos dos simulaciones variando Pe:
En la figura (8) se muestra el frente de la saturación al mismo instante, en cada uno de los casos. Se observa que los efectos de las difusión capilar son más importantes cuando la presión de entrada es débil. En este caso, la difusión capilar es preponderante sobre los efectos del término de transporte. Estos resultados coinciden con los resultados obtenidos en [16].
6. RECURSOS COMPUTACIONALES
La simulaciones numéricas presentadas en este artículo se llevaron a cabo en un computador con las siguientes características: 2 Quad-core Intel Xeon E5420 de 2.4 GHz y 32 Gb de RAM. El cuadro 2 muestra los tiempos de cálculo para cuatro mallados diferentes:
7. CONCLUSIONES
Hemos presentado el modelamiento y simulación del flujo bifásico agua- aceite en ausencia de gravedad, en un medio poroso homogéneo en tres dimensiones. Se utiliza una formulación IMPES para la resolución numérica del modelo. El sistema a resolver se compone de una ecuación para la velocidad, una ecuación de tipo elíptico para la presión y una ecuación de tipo parabólico para la saturación. Los resultados obtenidos muestran que los efectos capilares ayudan a extender el frente de saturación de la fase inyectada y además, que estos efectos son más importantes cuando la presión de inyección es débil.
AGRADECIMIENTOS
El autor de este artículo agradece al Instituto Nacional de Informática y Automática de Bordeaux, INRIA- Sud Ouest, por el financiamiento de este trabajo; al Instituto de Matemáticas de Bordeaux, IMB, por los medios de cálculo ofrecidos; a Mazen Saad, profesor de la Escuela Central de Nantes, por iniciar el código de Fortran que permitió desarrollar este trabajo y a Charles Henri Bruneau, profesor de la Universidad de Bordeaux, por dirigir este trabajo de maestría.
REFERENCIAS
1. D. E. Tognisso, Écoulements de fluides complexes en milieu poreux : utilisation de micelles géantes pour la Récupération Améliorée du Pétrole. PhD thesis, Bordeaux 1, Nov. 2011.
2. J. L. P. Pinilla, Modélisation et simulation àl´échelle du pore de la récupération assistée des hydrocarbures par injection de polyméres. PhD thesis, Bordeaux 1, Dec. 2012.
3. M. Quintard, H. Bertin, and S. Whitaker, "Two-phase flow in heterogeneous porous media: The method of large-scale averaging applied to laboratory experiments in a stratified system," Society of Petroleum Engineers, Oct. 1989.
4. T. Sochi, "Flow of non-newtonian uids in porous media," Journal of Polymer Science Part B: Polymer Physics, vol. 48, no. 23, p. 2437 - 2767, 2010.
5. G. Chavent and J. Ja_r_e, Mathematical Models and Finite Elements for Reservoir Simulation: Single Phase, Multiphase and Multicomponent Flows through Porous Media. Elsevier,Jan. 1986.
6. J. Jafré, "Formulation mixte d´écoulement diphasiques incompressibles dans un milieu poreux," 1980.
7. D. W. Peaceman, Fundamentals of Numerical Reservoir Simulation. Elsevier, Apr. 2000.
8. Z. Chen, G. Huan, and Y. Ma, Computational Methods for Multiphase Flows in Porous Media. SIAM, Apr. 2006.
9. S. Whitaker, "Flow in porous media i: A theoretical derivation of darcy´s law," Transp Porous Med, vol. 1, pp. 3{25, Mar. 1986.
10. D. Lasseux, M. Quintard, and S. Whitaker, "Determination of permeability tensors for twophase ow in homogeneous porous media: Theory," Transp Porous Med, vol. 24, pp. 107-137, Aug. 1996.
11. J. Sheldon and B. Zondek, "One-dimensional, incompressible, noncapillary, two-phase fluid flow in a porous medium," Society of Petroleum Engineers, vol. 216, pp. 290-296, 1959.
12. H. Stone and A. Garder Jr., "Analysis of gascap or dissolved-gas drive reservoirs," Society of Petroleum Engineers Journal, vol. 1, June 196.
13. C. Chueh, M. Secanell, W. Bangerth, and N. Djilali, "Multi-level adaptive simulation of transient two-phase ow in heterogeneous porous media," Computers & Fluids, vol. 39, pp. 1585-1596, Oct. 2010.
14. J. Kou and Sun, Shuyu, "On iterative IMPES formulation for two-phase ow with capillarity in heterogeneous porous media," International Journal of Numerical Analysis and Modeling, Series B, vol. 1, pp. 20 - 40, July 2010.
15. Z. Chen, G. Huan, and B. Li, "An improved IMPES method for two-phase ow in porous media," Transport in Porous Media, vol. 54, pp. 361- 376, Mar. 2004.
16. M. S. Saad, Propriétés de quelques modéles d´écoulements en milieu poreux/par Mazen Samir Saad. PhD thesis, Bordeaux 1, Jan. 1993.
17. F. Marpeau and M. Saad, "3D simulation of radionuclide transport in porous media," International Journal for Numerical Methods in Fluids, vol. 64, no. 1, p. 44{70, 2010.
18. M. R. Hestenes and E. Stiefel, Methods of conjugate gradients for solving linear systems. National Bureau of Standards, 1952.
19. R. J. LeVeque, Finite volume methods for hyperbolic problems. Cambridge [u.a.: Cambridge Univ. Press, 2003.