ARTÍCULOS ORIGINALES

Un modelo fractal para estimar la conductividad hídráulica no saturada de rocas fracturadas

 

Leonardo B. Monachesi y Luis Guarracino

CONICET. Facultad de Ciencias Astronómicas y Geofísicas, Universidad Nacional de La Plata, Argentina. Paseo del bosque s/n, 1900 La Plata, Argentina.: luisg@fcaglp.unlp.edu.ar; lmonachesi@fcaglp.unlp.edu.ar

 

---

RESUMEN

Durante las últimas dos décadas el estudio y la modelación del flujo de agua en rocas fracturadas no saturadas ha recibido considerable atención por parte de investigadores de distintas áreas como la geología, geofísica e ingeniería. Una de las razones principales para su estudio radica en la búsqueda de formaciones geológicas de baja permeabilidad para la construcción de repositorios de residuos nucleares. Para modelar el flujo de agua cuando se utiliza la aproximación del continuo es necesario contar con relaciones constitutivas del medio poroso considerado. Las relaciones constitutivas son curvas de saturación (S) y conductividad hidráulica (K) en función de la altura de presión (h). En general, la determinación en laboratorio de la curva de saturación no suele presentar mayores inconvenientes. En cambio, la determinación experimental de K(h) resulta difícil y costosa, por lo que se suele recurrir a modelos teóricos para su estimación. Entre los modelos teóricos más utilizados se destacan los de Burdine (Burdine, 1953) y Mualem (Mualem, 1976), los cuales permiten predecir la conductividad hidráulica K a partir de la curva de saturación S. Estos modelos predictivos han sido desarrollados para medios porosos de tipo sedimentario, no existiendo en la literatura modelos específicos para rocas fracturadas. En este trabajo se presenta un modelo fractal simple para predecir la conductividad hidráulica no saturada de rocas fracturadas. La deducción del modelo se basa en la hipótesis de que la red de fracturas puede ser descrita mediante un objeto fractal clásico denominado carpeta de Sierpinski. La expresión propuesta de K es cerrada y depende únicamente de tres parámetros independientes: la dimensión fractal y las aperturas máxima y mínima de las fracturas. Una de las características que presenta el modelo propuesto, y que permite validarlo en forma teórica, es que la curva de conductividad hidráulica relativa es equivalente a la que se obtiene con el modelo de Burdine cuando se asume un factor de tortuosidad unitario. Finalmente, se realiza un estudio comparativo entre el modelo analítico propuesto y las relaciones constitutivas obtenidas mediante simulación numérica por Liu y Bodvarsson (2001). La comparación muestra que el modelo propuesto puede predecir dichas relaciones constitutivas en un amplio rango de valores de saturación de agua, y constituye por lo tanto una validación adicional. El modelo de K(h) propuesto es único en su tipo ya que ha sido derivado íntegramente a partir de leyes físicas clásicas y propiedades geométricas de la red de fracturas.

Palabras clave: Conductividad hidráulica; Rocas fracturadas; Flujo no saturado

ABSTRACT

Modeling groundwater flow in unsaturated fractured rocks has received considerable attention in the last two decades. One of the main reasons for focusing on the study of water flow in this type of media is the search for potential safe permanent storage facilities for geological disposal of high-level nuclear wastes. Knowledge of constitutive relations is indispensable for the numerical solution of the equations describing water flow in unsaturated porous media. These relations are saturation (S) and hydraulic conductivity (K) curves, both expressed as functions of the pressure head (h). The experimental determination of K(h) is tedious and time-consuming and the measurements are variables, error-prone, and applicable to only a narrow range of pressure head h. An alternative to direct measurement is to use theoretical models which predict K(h) from the saturation curve S(h) that 68 can be easily measured in laboratory. The majority of the models that predicts K(h) from S(h) have been developed for describing unsaturated flow in sedimentary formations (granular porous media). Predictive models specifically designed for fractured hard rocks are virtually nonexistent. In this study, a hydraulic conductivity model for fractured rocks is derived. The proposed K(h) model is based on the assumption that the fracture pattern is self-similar. The fracture pattern is described using the Sierpinski carpet, a classical fractal object that contains a self-similar geometric pattern of pores. The proposed model has a closed form analytical expression with three independent geometric parameters: the fractal dimension of the Sierpinski carpet and the maximum and minimum fracture apertures. One of the main features of the proposed model that allows its validation is that the expression of K(h) is identical to the one obtained by using the Burdine model with a tortuosity factor equal to one. The proposed model can represent the constitutive relations for fractured rocks obtained by Liu and Bodvarsson (2001) using numerical simulation techniques. The proposed K(h) model is the only existing model which has been completely derived from physical concepts and geometric properties of the fracture pattern.

Keywords: Hydraulic conductivity; Fractured rocks; Unsaturated flow

---

 

INTRODUCCIÓN

La determinación experimental de la conductividad hidráulica en términos de la altura de presión , resulta difícil y costosa. Es por ello que una alternativa a la medición directa es el empleo de modelos teóricos que predicen la conductividad hidráulica a partir de la curva de saturación , cuya determinación en laboratorio es relativamente sencilla. Los modelos predictivos de son ampliamente utilizados y entre ellos se destacan los modelos de Burdine, Mualem y Assouline (Assouline, 2001).
Los modelos para estimar a partir de han sido desarrollados para medios porosos de tipo sedimentario, no existiendo en la literatura modelos específicos para rocas fracturadas. En el trabajo de Liu y Bodvarsson (2001) se ha utilizado el modelo de Burdine para predecir la conductividad hidráulica en rocas fracturadas empleando distintas expresiones para la curva de saturación. En el modelo de Burdine los poros se representan mediante tubos capilares paralelos, mientras que en los modelos de Mualem y Assouline se asumen geometrías más complejas. La elección del modelo de Burdine se basa precisamente en su simplicidad ya que resulta más adecuado para representar el flujo en fracturas.
En este trabajo se deriva un modelo de conductividad hidráulica para rocas fracturadas. El modelo se basa en la hipótesis de que el patrón de fracturas de la roca es de tipo fractal. La autosimilitud es una propiedad característica de los objetos fractales, que ha sido observada en redes de fracturas por varios autores (Turcotte, 1986; Obuko y Aki, 1987; Barton y Zoback, 1992; Berkowitz y Hadad, 1997). Para describir la red de fracturas utilizaremos la carpeta de Sierpinski que es un objeto fractal clásico que posee un patrón de poros autosimilar. La carpeta de Sierpinski ha sido utilizada previamente para describir medios porosos, tanto sedimentarios (Tyler y Wheatcraft, 1990) como fracturados (Guarracino, 2006).
El modelo propuesto se deriva a partir de las propiedades de la carpeta de Sierpinski, la ley de flujo entre 2 placas paralelas y relaciones capilares clásicas. La expresión analítica obtenida para la conductividad hidráulica depende de 3 parámetros geométricos: la dimensión fractal y las aperturas máxima y mínima de las fracturas de la red. La expresión de obtenida es idéntica a la que se obtiene con el modelo de Burdine cuando se asume un factor de tortuosidad unitario. Por otra parte el modelo permite ajustar las relaciones constitutivas simuladas por Liu y Bodvarsson (2001) para un medio fracturado sintético. El modelo de propuesto posee la ventaja de ser el único modelo para roca fracturada derivado íntegramente a partir de leyes físicas y propiedades geométricas de la red de fracturas.

DESCRIPCIÓN DE LA RED DE FRACTURAS

Para derivar el modelo de consideraremos un volumen representativo elemental (REV) de roca fracturada. El REV será conceptualizado como un cubo de lado a formado por una matriz impermeable y una red de fracturas que constituyen los poros del medio (ver Figura 1).

Figura 1: Volumen representativo elemental de roca fracturada.

Las fracturas se asumen verticales y su patrón de distribución se describe mediante una carpeta de Sierpinski. Este modelo de roca fracturada es idéntico al utilizado por Guarracino (2006) para derivar una curva de contenido de agua. La carpeta de Sierpinski que describe la distribución espacial de las fracturas se genera mediante un algoritmo recursivo (Turcotte, 1997; Tyler y Wheatcraft, 1990). Este algoritmo se aplica a un cuadrado de lado a correspondiente a una sección transversal del REV. El primer nivel de recursión se obtiene al subdividir el cuadrado original en cuadrados elementales de tamaño . El patrón de fracturas se genera removiendo una cantidad de estos cuadrados de manera que formen una red de fracturas de apertura 1 x como s e muestra en la Figura 2.

Figura 2: Dos niveles de recursión de una carpeta de Sierpinski de D=1.84

El patrón de fractura elegido debe repetirse en cada nivel de recursión para que el patrón resultante sea autosimilar. Para un nivel de recursión arbitrario , el cuadrado original se subdivide entonces en cuadrados elementales de tamaño , siendo . Puede demostrarse que en cada instancia de recursión la carpeta de Sierpinski verifica la siguiente relación (Mandelbrot, 1983):

donde es el número de cuadrados de tamaño necesarios para cubrir el área total ocupada por las fracturas de apertura , y es la dimensión fractal de la carpeta. A modo de ejemplo, en la Figura 2 se ilustran dos niveles de recursión para una carpeta de Sierpinski de dimensión fractal
.

Para un valor arbitrario de l , el área acumulativa de las fracturas de apertura está dada por:

donde representa el área de un cuadrado elemental. Reemplazando (1) en (2) se obtiene la siguiente expresión del área acumulativa en término de los parámetros de la carpeta de Sierpinski:

expresión válida para .
El área ocupada por las fracturas cuya apertura se encuentra en el rango se obtiene diferenciando (3) con respecto a :

El signo negativo de la ecuación anterior implica que el área de las fracturas disminuye con el aumento de la apertura. Resulta importante destacar que la carpeta de Sierpinski es un objeto fractal y por lo tanto admite infinitos niveles de recursión. En consecuencia, el rango de variación de la apertura x es infinito. Luego, para que este objeto fractal pueda ser utilizado para representar un medio fracturado real es necesario considerar valores de corte superior e inferior para la apertura. En este trabajo llamaremos y a las aperturas máxima y mínima de las fracturas dentro del REV.

MODELO DE CONDUCTIVIDAD HIDRÁULICA

En esta sección derivaremos una expresión analítica para la conductividad hidráulica de una roca fracturada. Para ello supondremos que el REV definido en la sección anterior se encuentra a una presión efectiva constante . Para este valor de presión, la teoría capilar clásica nos permite definir una apertura de fractura efectiva (Bear, 1988):

donde es la tensión superficial del agua, el ángulo de contacto, la densidad del agua y la aceleración de la gravedad.
El modelo propuesto asume que todas las fracturas de apertura se encuentran completamente saturadas de agua y se consideran fracturas activas (Liu et al, 1998). Las fracturas de apertura se encuentran parcialmente saturadas y se asume que no contribuyen al flujo (Guarracino, 2006). Dado que el rango de aperturas de fracturas se encuentra acotado por y , la ecuación (5) resulta válida para el rango de altura de presión , siendo y . Para todas las fracturas son activas, el medio se encuentra completamente saturado y el caudal es máximo. Contrariamente, si ninguna fractura se considera activa y por lo tanto la saturación y el flujo resultan nulos.
Dado que el REV se encuentra a una altura de presión constante , el gradiente efectivo de la altura de presión es nulo y por lo tanto el agua sólo fluirá por gravedad. Bajo esta hipótesis, el flujo de agua en una fractura puede expresarse mediante la conocida ecuación de Hagen-Poisseuille (Bear, 1988):

donde z es la coordenada vertical y la viscosidad.
El caudal volumétrico total Q a través de una sección horizontal del REV puede calcularse sumando la contribución de todas las fracturas activas:

Reemplazando (4) y (6) en (7) se obtiene la siguiente expresión para el caudal volumétrico total:

Teniendo en cuenta la ecuación (5), la (8) puede expresarse en términos de la altura efectiva :

Por otro lado el flujo volumétrico en el REV también puede calcularse a partir de la ley de Buckingham-Darcy (Buckingham, 1907):

donde es la conductividad hidráulica y representa el área de la sección transversal a la dirección del flujo. Combinando (9) y (10) obtenemos finalmente la siguiente expresión para la conductividad hidráulica:

La conductividad hidráulica saturada se obtiene considerando en (11):

Finalmente, la conductividad hidráulica relativa vendrá dada por el cociente entre las ecuaciones (11) y (12):

La expresión anterior es válida para . Para valores de menores a . Inversamente, para valores de mayores que .
La conductividad hidráulica relativa suele expresarse en términos de la saturación S . Para el tipo de medio considerado la expresión de fue calculada por Guarracino (2006):

Finalmente, reemplazando (14) en (13) se obtiene:

La ecuación (13) (o su análoga (15)), junto con la ecuación (14) conforman un modelo constitutivo para roca fracturada. Este modelo ha sido derivado enteramente a partir leyes físicas básicas y de propiedades geométricas de la red de fracturas y sólo depende de tres parámetros independientes: .

Validación del modelo
En esta sección la conductividad hidráulica relativa (15) será comparada con la expresión que se obtiene a partir del modelo predictivo de Burdine. Por otra parte se realizará un ajuste del modelo propuesto con las relaciones constitutivas obtenidas en forma numérica por Liu y Bodvarsson (2001). Las expresiones de la conductividad hidráulica se obtienen generalmente a partir de la curva de retención utilizando los modelos predictivos de Burdine ó Mualem. El modelo de Burdine asume que el flujo se desarrolla a través de un grupo de tubos capilares paralelos de diámetro constante, mientras que en el modelo de Mualem los tubos capilares son tortuosos y de diámetro variable. Resulta evidente que el flujo a través de una fractura puede representarse en forma adecuada mediante una serie de tubos capilares de diámetro constante y es por ello que el modelo de Burdine ha sido utilizado por distintos autores para obtener expresiones de para rocas fracturadas (Liu y Bodvarsson, 2001; Guarracino, 2006; Guarracino y Quintana, 2008).
La expresión del modelo de Burdine en términos de es la siguiente (Burdine, 1953):

donde es el factor de tortuosidad, cuyos valores oscilan entre 0 y 1. Cuando T =1 la trayectoria del flujo de agua entre dos puntos del medio poroso es recta mientras que si T = 0 la trayectoria entre dichos puntos es infinita. Luego, una expresión de para una roca fracturada puede obtenerse utilizando la curva de saturación (14) en (16):

En el modelo fractal propuesto el flujo de agua se produce a lo largo de trayectorias rectas definidas por las fracturas activas, por lo cual T =1 para todo el rango de saturación. Para este valor de tortuosidad se obtiene exactamente la expresión de la conductividad hidráulica relativa (13) correspondiente al modelo fractal. Resulta importante destacar que el modelo de Burdine fue derivado para medios sedimentarios, donde las trayectorias seguidas por las partículas de agua no son rectilíneas y dependen fuertemente del grado de saturación. Sin embargo la comparación del modelo fractal con el modelo de Burdine nos permite extender la validez de este último para el caso de rocas fracturadas a través de la modificación del factor de tortuosidad.
La modificación del factor de tortuosidad del modelo de Burdine fue propuesta originalmente por Liu y Bodvarsson (2001) para estimar la conductividad hidráulica en un medio fracturado con baja conectividad entre las fracturas. En la Figura 3 se ilustra el factor de tortuosidad propuesto por Liu y Bodvarsson (2001) junto con los factores del modelo de Burdine y del modelo propuesto. El modelo fractal representa el caso ideal en que la trayectoria del flujo de agua es rectilínea y no es modificada por la geometría del medio ni por el grado de saturación. El factor de tortuosidad de Liu y Bodvarsson (2001) representa un caso intermedio entre un medio fracturado sin tortuosidad (T =1) y un medio sedimentario homogéneo .

Figura 3. Factores de tortuosidad de los modelos de Burdine, Liu y Bodvarsson y el propuesto.

La determinación experimental de las relaciones constitutivas para rocas fracturadas presenta numerosas dificultades técnicas, por lo que los datos disponibles para la validación del modelo son prácticamente inexistentes (Liu y Bodvarsson , 2001; Tuller y Or, 2002). En este trabajo validaremos los resultados mediante su comparación con las relaciones constitutivas simuladas por Liu y Bodvarsson (2001) para una roca fracturada. Dichos autores diseñaron un experimento computacional similar al procedimiento de laboratorio realizado para determinar relaciones constitutivas. Para ello consideran un REV impermeable con fracturas aleatorias conceptualizadas como medios porosos bidimensionales. Para un valor de presión capilar en los bordes del REV resuelven la ecuación de Richards y calculan la saturación efectiva y la conductividad hidráulica relativa . Mediante este procedimiento, Liu y Bodvarsson (2001) obtienen relaciones constitutivas para dos redes de fracturas sintéticas, las cuales serán utilizadas para validar las estimaciones de ( ) r K S del modelo propuesto.
Los parámetros del modelo fractal ( min max D, h , h ) se obtienen a partir del ajuste de la curva de saturación (14) a los datos simulados. En la Tabla 1 se listan los valores de los parámetros ajustados para las dos redes de fracturas. Utilizando la ecuación (15) con los parámetros ajustados se estimó la conductividad hidráulica relativa en función de la saturación efectiva. En las Figuras 4 y 5 se muestran las estimaciones del modelo propuesto junto con los resultados numéricos obtenidos por Liu y Bodvarsson (2001) para las dos redes de fracturas consideradas.

Tabla 1. Parámetros ajustados del modelo fractal

Figura 4. Comparación del modelo propuesto con los valores simulados por Liu y Bodvarsson para la red de fracturas 1.

-
Figura 5. Comparación del modelo propuesto con los valores simulados por Liu y Bodvarsson para la red de fracturas 2.

En ambos casos el modelo fractal predice muy bien los datos de para saturaciones altas e intermedias. Para valores bajos de saturación (menores a 0.3 aproximadamente) el modelo propuesto no logra predecir los valores calculados por Liu y Bodvarsson (2001). Esta discrepancia se debe a que las curvas del experimento numérico poseen una tortuosidad dependiente del grado de saturación mientras que el modelo propuesto no contempla en forma explícita este fenómeno. Para saturaciones bajas el flujo se desarrolla a través de las fracturas de menor apertura que presentan menor conectividad entre ellas. Esta situación corresponde a un valor del factor de tortuosidad bajo que se aparta del modelo propuesto (T =1).
Por otra parte, para saturaciones altas e intermedias la mayor contribución al flujo corresponde a las fracturas de mayor apertura, que son las que poseen mayor conectividad y un factor de tortuosidad cercano a 1. Es de destacar que, a pesar de su simplicidad, el modelo propuesto logra predecir la conductividad hidráulica en más de tres órdenes de magnitud, en el rango de saturaciones altas e intermedias para las cuales el flujo de agua resulta significativo. Resulta evidente que un mejor ajuste a los datos simulados podría obtenerse introduciendo en el modelo un factor de tortuosidad dependiente de la saturación. Sin embargo el análisis de este fenómeno (que puede ser significativo en el rango de bajas saturaciones) se encuentra más allá de los objetivos del presente trabajo.

CONCLUSIONES

En este trabajo se ha derivado un modelo analítico para predecir la conductividad hidráulica en rocas fracturadas parcialmente saturadas. El modelo se basa en una descripción fractal de la red de fracturas utilizando una carpeta de Sierpinski. La expresión de la conductividad hidráulica propuesta, junto con la correspondiente curva de retención, conforman el único modelo constitutivo para rocas fracturadas derivado íntegramente a partir de propiedades geométricas y leyes físicas. El modelo posee tres parámetros independientes que definen la geometría de la red de fracturas: la dimensión fractal y las aperturas máxima y mínima. La expresión de la conductividad hidráulica relativa coincide con la curva que predice el modelo clásico de Burdine cuando se considera un factor de tortuosidad unitario. Por otra parte el modelo fue validado mediante la comparación con relaciones constitutivas obtenidas en forma numérica, obteniendo un muy buen ajuste.

REFERENCIAS

1. Assouline S., 2001. A model for soil relative hydraulic conductivity based on the water retention characteristic curve. Water Resour. Res. 37 (2), 265-271.         [ Links ]

2. Barton, C.A., Zoback, M.D., 1992. Self-similar distribution and properties of macroscopic fractures at depth in crystalline rock in the Cajon pass scientific drill hole. J. Geophys. Res. 97, 5181-5200.         [ Links ]

3. Bear, J., 1988. Dynamics of Fluids in Porous Media, Elsevier, New York.         [ Links ]

4. Berkowitz, B., Hadad, A., 1997. Fractal and multifractal measure of natural and synthetic fracture networks. J. Geophys. Res. 102 (B6), 205-218.         [ Links ]

5. Buckingham, E., 1907. Studies on the movement of soil moisture. Washington DC, U.S. Department of Agriculture. Bureau of Soil Bulletin 38, 29-61.         [ Links ]

6. Burdine, N.T., 1953. Relative permeability from pore size distribution data. Trans. Am. Inst. Min. Metall. Pet. Eng. 198, 71-78.         [ Links ]

7. Guarracino, L., 2006. A fractal constitutive model for unsaturated flow in fractured hard rocks. J. Hydrol. 324, 154-162.         [ Links ]

8. Guarracino, L. Quintana, F., 2008. A constitutive model for flow in unsaturated fractured rocks. Hydrological Processes. In press.         [ Links ]

9. Liu, H. H., Bodvarsson, G.S., 2001. Constitutive relations for unsaturated flow in fracture network. J. Hydrol. 252, 116-125.         [ Links ]

10. Liu, H. H., Doughty, C., Bodvarsson, G.S., 1998. An active fracture model for unsaturated flow and transport in fractured rocks. Water Resour. Res. 12 (3), 513-522.         [ Links ]

11. Mualem, Y., 1976. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 12 (3), 513-522.         [ Links ]

12. Mandelbrot B. B., 1983. The Fractal Geometry of Nature. W. H. Freeman, Ney York, 488 pp.         [ Links ]

13. Obuko, P.G., Aki, K., 1987. Fractal geometry in the San Andreas fault system. J. Geophys. Res. 92, 345-355.         [ Links ]

14. Tuller M., Or D., 2002. Unsaturated hydraulic conductivity of structured porous media: a review of liquid configuration-based models. Vadose Zone Journal 1:14-37.         [ Links ]

15. Turcotte D.L., 1986. Fractals and fragmentations. J. Geophys. Res.91, 1921-1926.         [ Links ]

16. Turcotte D.L., 1997. Fractals and Chaos in Geology and Geophysics. Cambridge University Press, New York, 398 pp.         [ Links ]

17. Tyler, S.W., Wheatcraft, S.W., 1990. Fractal process in soil water retention. Water Resour. Res. 26 (5), 1047-1054.         [ Links ]

Recibido: 26 de agosto de 2009
Aceptado: 9 de junio de 2010