Hybrid percolation and finite element model for the micromecanical analysis of soils

María de la Luz Pérez–Rea^1,*, Jaime Horta–Rangel¹, Carlos Santiago López–Cajún¹, Teresa López–Lara¹, Juan Bosco Hernández–Zaragoza¹ y Víctor Manuel Castaño²

¹ Universidad Autónoma de Querétaro, División de Estudios de Posgrado, Facultad de Ingeniería, Cerro de las Campanas s/n, Col. Niños Héroes, 76010 Querétaro, Qro., México. * perea@uaq.mx

² Centro de Física Aplicada y Tecnología Avanzada, Universidad Nacional Autónoma de México, Boulevard Juriquilla 3001, 76230 Querétaro, Qro., México.

Manuscrito recibido: Marzo 1, 2008 ]]> Manuscrito corregido recibido: Enero 22, 2009
Manuscrito aceptado: Enero 30, 2009

RESUMEN

Los modelos de percolación han sido utilizados en muchas áreas de la ciencia para modelar diversas propiedades de los materiales. En este trabajo se presenta una revisión de dichos modelos y se hace una descripción técnica de la teoría de percolación. También se describen aspectos micromecánicos de suelos, definidos por medio de diversas técnicas experimentales. Las ventajas de combinar la teoría de percolación con la simulación por elemento finito se ejemplifican por medio del análisis del comportamiento esfuerzo–deformación de modelos de redes regulares, que representan la microestructura de suelos arcillosos, para estimar el módulo elástico que tendrá macroscópicamente el suelo. Se muestra que los resultados obtenidos son muy cercanos a los determinados experimentalmente. Finalmente se discuten los alcances de esta nueva metodología.

Palabras clave: suelo, percolación, elemento finito, micromecánica, simulación computacional, modelado.

ABSTRACT

Percolation models have been used in many areas of science for the modeling of diverse properties of materials. In this paper, a revision of such models, as well as a technical description of the percolation theory is presented. Micromechanical aspects of soils, defined by means of diverse experimental techniques, are also discussed. The advantages of combining the percolation theory with the finite element method is exemplified through the analysis of the stress–strain behaviour in regular network models, which represent the microstructure of clayey soils, to estimate the elastic module of the macroscopic soil. The obtained results closely resemble experimentally determined values. The scope of this new methodology is also discussed.

Key words: soil, percolation, finite element, micromechanics, computational simulation, modeling.

INTRODUCCIÓN

En estudios recientes se ha intentado relacionar la conductividad hidráulica de suelos finos compactados con cambios en la microestructura del suelo en tres diferentes niveles: microestructural, miniestructural, y macroestructural. El primero considera los pequeños poros entre las partículas individuales; el segundo, considera los poros entre los grumos y, el tercero, considera los macroporos tales como grietas, fisuras, vacíos productos de las raíces, etc. (Collis y McGown, 1974; Mitchell et al. 1973). Desde los inicios de la Mecánica de Suelos como ciencia, se ha destacado que la microestructura de los suelos juega un papel importante en el comportamiento macroscópico de los mismos, y aunque la mayoría de los tratados de esta ciencia introducen grandes capítulos sobre mineralogía y micromorfología de suelos, no hacen un vínculo con el comportamiento a gran escala. Estructuralmente, sólo las grietas y fisuras visibles en suelos han sido objeto de investigación para ser tomadas en cuenta sin llegar a conclusiones claras. La complejidad en la composición mineralógica del suelo y en las configuraciones microestructurales que éste adquiere, ha dificultado su modelado micromecánico ya que se le observa como un medio desordenado. Sin embargo, con el desarrollo de poderosos equipos de cómputo y mayores capacidades de cálculo, la ingeniería de suelos comienza a enfocarse al modelado de los fenómenos cuyo comportamiento es muy importante como el flujo de agua y el comportamiento mecánico. Este modelado ha sido tomado desde el punto de vista de la Mecánica de Medios Continuos, sin tomar en cuenta la microestructura. Sin embargo, algunas teorías modernas de modelado probabilístico pueden ser usadas para encontrar representaciones adecuadas de microestructura de materiales. En las últimas dos décadas, el desarrollo de una familia de métodos teóricos poderosos ha permitido interpretar observaciones experimentales y predecir muchas propiedades de los sistemas desordenados (Sahimi, 1994). La teoría de percolación, en particular, ha jugado uno de los papeles más importantes para el entendimiento de sistemas desordenados y sus propiedades. Los problemas relativos a medios porosos han recibido particular atención. Por ejemplo, los problemas de flujo bifásico han sido ya, en buena medida, entendidos, y, al menos bajo ciertas condiciones, pueden ser interpretados como un fenómeno de percolación. La física estadística de sistemas desordenados permite desarrollar métodos que deriven las propiedades macroscópicas de tales sistemas a partir de las leyes que gobiernan el mundo microscópico, o bien, deducir propiedades microscópicas a partir de la información macroscópica que puede observarse con técnicas experimentales. El rápido desarrollo de estas teorías se debe a que se ha logrado apreciar el importante papel de la interconectividad física y química de los elementos microscópicos y sus efectos en las propiedades macroscópicas.

Si bien la mayoría de los trabajos reportados en la literatura tratan problemas de conductividad hidráulica y conductividad eléctrica, existen otras propiedades que también podrían ser modeladas, tales como las propiedades elásticas. Reportes de análisis de materiales constituidos por celdas octaédricas repetidas muestran que el método del elemento finito puede aplicarse para la predicción de propiedades elásticas efectivas en materiales porosos como las espumas metálicas (Deshpande et al. 2001).

ANTECEDENTES

Teoría de percolación

La teoría de percolación ha sido inspirada en el proceso sencillo de elaboración de café en una cafetera percoladora, el cual en principio es el flujo del agua a través de un medio poroso. En términos más generales podría decirse que es la transmisión de una propiedad a través de un medio, donde se toman en cuenta las propiedades del medio. Según reporta Sahimi (1994), el uso moderno de los procesos de percolación fue desarrollado por Flory (1941) y Stockmayer (1943) para describir la reacción de moléculas ramificantes que formaban macromoléculas, trabajo que a Flory le valió posteriormente el Premio Nobel. Este proceso de polimerización puede conducir a la formación de una red muy grande de moléculas conectadas por enlaces químicos. Después, Broadbent y Hammersley (1957) trataron con el concepto de la dispersión de partículas de fluido hipotéticas a través de un medio aleatorio, utilizando, precisamente, la teoría de percolación. Los términos fluido y medio se propusieron de manera totalmente general: fluido puede ser un líquido, vapor, flujo de calor, corriente eléctrica, infección, sistema solar o cualquier fluido o propiedad que pueda moverse a través de un medio. El medio, donde se lleva al fluido, puede ser el espacio poroso de una roca, un suelo, un arreglo de árboles o el mismo universo. La dispersión de un fluido a través de un medio desordenado involucra algunos elementos aleatorios, pero el mecanismo puede ser uno de dos diferentes tipos. En el primer tipo, la aleatoriedad se atribuye al fluido en un proceso de difusión y en el segundo, se atribuye al medio en un proceso de percolación.

La teoría de percolación también se utiliza para describir el modelo simple de la penetración de mercurio en un medio poroso y sus propiedades durante el proceso (Chatzis y Dullien, 1977). La penetración con mercurio es el método más común para medir la porosidad del suelo. Schmidt y McDonald (1979), mientras estudiaban las areniscas, analizaron el proceso diagenético de las mismas haciendo uso de la teoría de percolación. Barton y Hsieh (1989) estudiaron patrones de agrietamiento de formaciones montañosas usando esta teoría y Sahimi et al. (1990) reportaron numerosas aplicaciones de la teoría de percolación al estudio del flujo a través de medios porosos.

Los modelos de percolación se inician con la construcción de una red que representa al medio o a la propiedad que se quiera modelar. Estas redes pueden representar estructuras de los materiales o bien, las formas en que las propiedades se distribuyen en el medio. Las redes se constituyen de celdas unitarias repetidas en 2D o 3D (redes regulares) o de elementos como bastones o esferas dispuestos aleatoriamente (redes irregulares), lo que origina dos entidades principales: sitios (nodos) y lazos (Sahimi, 1994).

Análogamente, en el problema de percolación de sitio, los sitios de la red están ocupados con probabilidad p y las vacantes con probabilidad 1–p. Los dos sitios vecinos más cercanos se llaman "conectados" si están ambos ocupados, y los racimos conectados en la red se definen, de nuevo, de la manera obvia. Existe un umbral de percolación de sitio p_cs, por encima del cual un racimo infinito de sitios ocupados atraviesa la red. Se debe apuntar que, dependiendo de una aplicación específica, se pueden desarrollar muchas variantes de la percolación enlazada o de sitio, o combinada.

La deducción de los valores exactos de p_cb y p_cs ha sido posible sólo para ciertas redes relacionadas a la llamada "red de Bethe" y para algunas cuantas redes bidimensionales. El umbral de percolación de una red tridimensional ha sido calculado numéricamente por simulaciones Monte Carlo o por otras técnicas.

La Figura 1 muestra algunas redes bidimensionales regulares. En las Tablas 1 y 2 (Sahimi, 1994) se resumen estimaciones típicas de p_cb y p_cs para redes bidimensionales y tridimensionales, respectivamente. En general, p_cb < p_cs.

El valor de la probabilidad crítica decrece conforme aumenta la dimensión; por ejemplo, un agregado no percolante en una red en 2D puede ser una sección transversal de un agregado percolante en 3D. Si se supone una red cuadrada, aún cuando todos los agregados fueran finitos (p < p_cr = 0.5), podrían ser parte de un agregado infinito en 3D (p < p_cr = 0.2488). Otra diferencia fundamental entre la percolación en 2D y 3D es que en 2D la percolación ocurre en una sola fase a la vez y nunca podrá ocurrir simultáneamente en las dos fases (la fase ocupada y la fase vacía). En cambio, en 3D, la percolación puede ocurrir simultáneamente en dos o más fases. Y aunque el comportamiento de los fenómenos en 3D es muy diferente al comportamiento en 2D, frecuentemente se hacen simplificaciones para reducir los problemas tridimensionales a dos dimensiones debido a las limitaciones de procesamiento de los equipos de cómputo existentes.

Además del umbral de percolación, las propiedades topológicas de las redes de percolación se caracterizan por varias cantidades importantes:

Probabilidad de percolación P(p). Esto es, la probabilidad de que cuando una fracción de ligas ocupadas es p, un sitio dado pertenezca al racimo (percolante) infinito de ligas ocupadas.

Fracción accesible X^A(p). Esto es, la fracción de ligas ocupadas pertenecientes al racimo infinito.

Fracción de soporte X^B(p) Esto es, la fracción de ligas ocupadas en el racimo infinito, el cual lleva el flujo o la corriente. El soporte de un sistema de percolación juega un papel fundamental en sus propiedades de transporte, porque la tortuosidad de las trayectorias del transporte se controla por la estructura del soporte.

Longitud de correlación ξ_p(p). Esto es, el radio típico de los racimos conectados para p < p_c, y la escala de longitud sobre la que la red aleatoria es macroscópicamente homogénea (e.g., la escala de la longitud sobre la cual las propiedades del sistema son independientes de su tamaño lineal L) para p > p_c. Entonces, en cualquier simulación Monte Carlo de percolación se debería tener L >> ξ_p para que el resultado sea independiente de L.

Conductividad eléctrica efectiva g_e. Esto es, la conductividad eléctrica de una red de resistencia aleatoria en la cual una fracción p de ligas es conductora y el resto son aislantes. Análogamente, si una red representa el espacio del poro en un medio poroso en el cual una fracción p de los poros está abierta al flujo o la fusión, una difusividad efectiva Dr y una permeabilidad hidrodinámica k se pueden definir.

Módulo elástico efectivo E. Esto es, el módulo elástico de la red en el cual una fracción p de las ligas son elementos elásticos (e.g., resortes), mientras que el resto no tienen rigidez alguna (están cortados).

Las propiedades topológicas tales como la fracción accesible o de apoyo se calculan normalmente por simulaciones Monte Carlo. Por ejemplo, Stauffer et al. (1982) proporcionan un programa de computadora para calcular X^A(p), mientras Liem y Jan (1988) discuten un método para calcular X^B(p). Las propiedades de transporte, tales como la conductividad efectiva o el módulo elástico, se pueden calcular por simulaciones Monte Carlo o por aproximaciones analíticas.

Donde los exponentes de las ecuaciones son los exponentes críticos y caracterizan a las propiedades de las redes. Isichenko (1992) presenta un resumen de algunos de los exponentes críticos de percolación. Estos exponentes dependen solamente de la dimensión en el espacio y no del tipo de red o de la clase de problema de percolación; por esta razón se les llama también "universales". Los exponentes críticos se obtienen de las leyes de potencia que definen el comportamiento de las propiedades de las redes. En la mayoría de los casos y utilizando las ecuaciones anteriores, sólo dos de ellos se calculan a partir de los principios básicos. Se considera al exponente de longitud de correlación v y al exponente de densidad del agregado infinito β como los exponentes de percolación básicos. La elección de estos exponentes es relativa, ya que se consideran básicos porque fueron los primeros en ser determinados analíticamente en 2D. Los demás exponentes se obtienen usando argumentos de escalamiento o de autosimilitud.

En dos dimensiones (D=2), estos exponentes se conocen analíticamente, v = 4/3 y β = 5/36. Para D=3, sólo se tienen disponibles aproximaciones numéricas: v0.90 y β0.40. A estos dos exponentes críticos se les llama también exponentes "estáticos" porque caracterizan sólo la geometría y distribución de agregados. El resto de los exponentes pueden expresarse como función de los exponentes básicos y se les denomina "dinámicos". Éstos expresan las propiedades físicas correspondientes a redes de percolación que representan medios desordenados tales como la conductividad, la elasticidad, o la permeabilidad. Los agregados de percolación se utilizan usualmente para modelar objetos físicos, como sólidos amorfos, materiales compuestos, rocas porosas y polímeros, entre otros. El concepto de universalidad de los exponentes críticos reduce el enorme número de modelos diferentes a unos pocos, cuyos conjuntos de exponentes críticos son calculados y tabulados. El cálculo de los exponentes no es fácil. En cierto sentido, es una parte importante que estudia la física de fenómenos críticos. En estos términos, el problema de percolación pertenece a una de las clases de universalidad más simples, la cual es usualmente llamada la "clase de universalidad de percolación aleatoria, o no–correlacionada". Los exponentes críticos universales v y β no dependen del tipo de red (cuadrada, triangular, o cualquiera de las descritas) sino únicamente de la dimensión D. Esta afirmación ha sido confirmada por una vasta variedad de experimentos numéricos reportados en la literatura.

En este sentido, una de las sugerencias originales del presente trabajo es el que los modelos de elemento finito pueden utilizarse para analizar el comportamiento mecánico de redes percolantes.

Propiedades elásticas de sistemas percolantes

El problema de elasticidad de redes aleatorias percolantes se considera análogo al problema de la conductividad eléctrica (de Gennes, 1976), basándose en el modelo de Born para elasticidad microscópica de una celda. Para celdas cúbicas simples en cualquier dimensión, el umbral de rigidez ocurre en p=1. Los estudios de este problema se han limitado a ciertos tipos de celdas. Para tales celdas (e.g., una celda triangular en 2D) el valor del exponente de elasticidad es numéricamente mayor que el exponente de conductividad. Está claro, sin embargo, que el problema de percolación asociado con este tipo de elasticidad de celda difiere para una percolación regular de sitio o de lazo. Kantor y Webman (1984) propusieron un modelo para la elasticidad de una red percolante. En su modelo, el umbral de rigidez era idéntico al geométrico. La estructura del racimo percolante infinito ha sido analizada junto con un modelo de cadenas elásticas y se obtuvo el exponente crítico . Los resultados mostraron que el comportamiento elástico de la red percolante pertenece a una clase de universalidad diferente a la conductividad y que el exponente crítico es considerablemente mayor que el exponente de conductividad, t. Guéguen et al. (1997) estudiaron las propiedades físicas de rocas realizando, entre otros, modelos de percolación que tomaban en cuenta la microestructura de las mismas. Para los modelos de percolación utilizaron redes homogéneas regulares de lazos y sitios y encontraron que entre los modelos que analizaron, la teoría de percolación presentaba el marco más adecuado para tomar en cuenta algo más que simplemente la porosidad de las rocas en la predicción de las propiedades mecánicas y la localización de las fracturas. Taya (1999), reportó modelos micromecánicos de compósitos fibrorreforzados inteligentes SMA (shape memory alloy), los cuales fueron modelados mediante redes de percolación aleatorias, las cuales mostraron ser muy cercanas a resultados experimentales limitados. Por su parte, Pabst et al. (2007) efectuaron durante más de cinco años estudios en materiales cerámicos porosos utilizando diferentes modelos de predicción, incluyendo la teoría de percolación, y demostraron que la mayoría de los modelos predictivos para las propiedades efectivas de suspensiones y materiales porosos son leyes de potencia, o relaciones exponenciales en función del umbral del percolación. Estudios similares para observar el comportamiento de propiedades elásticas en compósitos de FRP (Fiber–reinforced polymers) realizados por Rastellini et al. (2008), dejan entrever que existe una relación intrínseca entre el volumen de la fracción de fibrorrefuerzo y el módulo elástico efectivo del compósito.

En opinión de los autores, el modelo de celdas proporciona una descripción correcta del comportamiento elástico de materiales compuestos no homogéneos, macroscópicamiente construidos con regiones localmente rígidas (sólidos) y regiones localmente muy suaves (vacíos). Cerca de p_c la rigidez macroscópica del material estará determinada por la elasticidad de largos y tortuosos canales de material rígido, el cual está contenido en el esqueleto principal del racimo percolante.

El material a modelar en este trabajo es el suelo, pero la metodología puede aplicarse a otros medios porosos. El comportamiento del suelo ha sido tradicionalmente modelado con la teoría de elasticidad con resultados aceptables (Poulos y Davis, 1974; Davis y Selvadurai, 1996). Las propiedades básicas del suelo pueden ser determinadas utilizando técnicas de caracterización tales como la microscopía electrónica, la difracción de rayos X y el análisis térmico, las que, debidamente interpretadas, proporcionan información que permite modelar más cercanamente las condiciones reales del suelo. Una vez que la porosidad del suelo ha sido modelada con la teoría de percolación, el método del elemento finito (FEM, por las siglas en inglés) puede aplicarse para el análisis del comportamiento mecánico. El análisis de micromecánica facilita la comparación entre los resultados obtenidos con el modelo de elemento finito y los reales.

Para este trabajo se consideraron dos redes regulares en 2D: la red cuadrada y la de estructura hexagonal o panaloide, en sistemas de percolación de lazos. La red cuadrada por ser la más utilizada para diferentes tipos de modelo y por la facilidad de encontrar referencias y la red hexagonal por ser la que más se asemeja a la microestructura que presenta el suelo, según las teorías basadas en los reportes existentes de microscopía electrónica (Juárez y Rico, 1984). Cada lazo representa partículas individuales de arcilla, las cuales tienen dos dimensiones pequeñas y una grande, de 1 y 10 micrómetros de tamaño medio, respectivamente. Las dimensiones fueron determinadas mediante microscopía electrónica de barrido. La suposición de una sección rectangular para estas microvigas obedece a que las partículas de arcilla tienen forma de placas y los poros tienen paredes planas, lo cual se dedujo indirectamente de la forma de la isoterma de adsorción en la caracterización de la porosimetría con adsorción de nitrógeno líquido (Pérez–Rea et al. 2003). Ambas redes proporcionan una descripción correcta del comportamiento elástico de un material compuesto macroscópicamente no homogéneo, constituido de regiones rígidas y regiones que son localmente muy suaves (estas regiones suaves, representan los vacíos del suelo). Cada una de las partículas de arcilla (montmorilonita, en este caso), tiene sus propiedades individuales, como módulo elástico, E=150 GPa, y relación de Poisson, v=0.35, que han sido medidas y reportadas en la literatura (Lee, 2006). Considerando que los granos de arcilla son prácticamente indestructibles, debido a los enlaces iónicos con que se encuentran ligadas las láminas de arcilla (Keedwell, 1984), cada uno de estos granos puede comportarse individualmente como una viga elástica, lo que no necesariamente se reflejará en el comportamiento macroscópico del suelo. La hipótesis principal supone que la transmisión de los esfuerzos dentro de la masa del suelo se lleva a cabo, a nivel microscópico, a través de los puntos de contacto entre los granos, los cuales podrían ser considerados como uniones de tipo rígido (Keedwell,1984) en suelos con bajos contenidos de humedad. A partir de esta hipótesis, demostrada macroscópicamente en suelos granulares (Juárez y Rico, 1984) y bajo la suposición de que todos los poros se encuentran interconectados entre sí, se tiene un compósito de dos fases: la fase sólida, a la cual constituyen los granos de arcilla dispuestos en diferentes estructuras, y la fase vacía, a la cual constituye el aire dentro de los poros. Esto es válido partiendo del hecho de que la mayoría de los esfuerzos que transmite un edificio al suelo no saturado son absorbidos casi inmediatamente por la parte sólida del suelo y los efectos de deformación asociados a la dificultad para desalojar el agua de los poros son mínimos, comparados con las deformaciones iniciales.

El método del elemento finito

El método del elemento finito es un método numérico que ha llegado a ser una herramienta poderosa en la solución de problemas de ingeniería para conocer el estado de esfuerzo y de deformación en un medio continuo. El medio puede ser heterogéneo y anisótropo; sólido, líquido o gas. Las aplicaciones van desde el análisis esfuerzo–deformación de sólidos como automóviles, aeronaves y edificios, hasta problemas complejos de flujo que pueden modelarse con relativa facilidad. En este método de análisis, una región compleja que define un continuo se discretiza en formas geométricas simples llamadas elementos finitos. Las propiedades del material y las relaciones gobernantes, son consideradas sobre esos elementos y expresadas en términos de valores desconocidos en las fronteras del elemento. Considerando debidamente las cargas y restricciones, se obtienen un conjunto de ecuaciones cuya solución proporciona el comportamiento aproximado del continuo.

El método del elemento finito se originó con la necesidad de realizar análisis estructurales en aeronaves (Hrenikoff, 1941) y posteriormente otros autores hicieron modificaciones al método de Hrenikoff para modelar problemas de torsión, armaduras y otros elementos estructurales (Courant, 1943, Turner et al. 1956). No fue sino hasta 1960 que se empleó el término elemento finito para describir el método (Clough, 1960). Las bases matemáticas se fijaron en la década de 1970 y el desarrollo actual de las computadoras y los procesadores han puesto el método al alcance de la mano, con aplicaciones en micromecánica para el análisis del comportamiento de materiales estructurados (Bolander y Saitot, 1997; Deshpande et al. 2001; Murad y Moyne, 2002, Shan y Gokhale , 2002).

El método del elemento finito para el análisis de estructuras (cálculo de los desplazamientos) consiste en dividir el continuo en líneas imaginarias que constituye cada uno de los elementos finitos (mallado). Los elementos están conectados entre sí mediante un número discreto de nodos situados en sus fronteras. Los desplazamientos de estos nodos son las incógnitas fundamentales del problema. Se toma un conjunto de funciones que definen de manera única el campo de desplazamientos dentro de cada elemento finito en función de los desplazamientos nodales (funciones de forma). Estas funciones de desplazamiento definen el estado de deformación del elemento en función de los desplazamientos nodales. Estas deformaciones, en conjunto con las propiedades constitutivas del material, definen el estado de esfuerzo en todo el elemento y consecuentemente, en sus fronteras. Se determina un sistema de fuerzas concentradas en los nodos, tal que equilibre los esfuerzos en las fronteras y cualesquiera cargas repartidas, lo cual se expresa con la relación F=ku ; donde F es el campo de fuerzas, k es la rigidez y u es el campo de desplazamientos. La solución de la ecuación que integra a todo el sistema permite obtener los desplazamientos nodales.

El elemento finito más sencillo, es quizás el elemento finito viga, el cual puede ser usado en análisis 1D, 2D y 3D. Es un elemento lineal con dos nodos, que puede tener varios grados de libertad en los nodos. Es el elemento más utilizado para el análisis estructural. También se pueden utilizar elementos finitos en el plano en forma de triángulos o cuadrados. Para análisis 3D son utilizadas formas volumétricas como tetraedros o cubos. La combinación de estas formas geométricas, así como la cantidad de nodos y las características de cada elemento finito permite modelar cualquier morfología de continuo, haciendo las mallas más refinadas donde se requiera conocer el comportamiento simulado con mayor detalle.

MATERIALES Y MÉTODOS

El material utilizado para la calibración del modelo corresponde al suelo obtenido de un pozo a cielo abierto en el Fraccionamiento Jurica de la Ciudad de Querétaro, (Figura 2). Se realizaron sondeos hasta 1.50 m de profundidad y se obtuvieron muestras inalteradas de suelo a las cuales se les midieron sus propiedades físicas. El módulo elástico y la relación de Poisson fueron medidas mediante pruebas de compresión triaxial realizadas en especímenes remoldeados, en los que se reprodujo la humedad natural, variando la compacidad. Se dejaron en reposo durante cinco días para que, de acuerdo con las propiedades trixotrópicas de la arcilla, se recuperara la estructura. La presión confinante en todas las pruebas fue la misma y equivalente a la presión horizontal in situ.

Modelo de elemento finito aplicado a redes

Para la aplicación del método del elemento finito (MEF) a las redes se recurrió al software comercial ANSYS^MR, versión 7.1, el cual es una herramienta computacional diseñada para proporcionar a los ingenieros la facilitad de realizar tareas, entre otras, como la construcción de modelos de computadora para estructuras, productos, componentes o sistemas; aplicar cargas de operación u otras condiciones de diseño y estudiar las respuestas físicas, tales como esfuerzos, distribución de temperaturas o campos electromagnéticos. Un análisis típico de ANSYS consta de tres pasos fundamentales: construcción del modelo, aplicación de cargas y obtención de la solución y, finalmente, revisión de los resultados. La construcción del modelo es la parte del análisis que requiere la mayor parte del tiempo. En esta etapa deben definirse el tipo de elemento o elementos que van a utilizarse, las constantes reales de los elementos, las propiedades materiales y la geometría del modelo. El sistema de ANSYS no asume un sistema de unidades para el análisis. Se pueden escoger las unidades que se deseen pero debe tenerse cuidado de que dichas unidades sean consistentes en todos los datos de entrada. La librería de ANSYS contiene más de 150 tipos de elemento finito. Cada elemento tiene un número único y un prefijo que identifica la categoría del elemento: e.g., BEAM4, PLANE77, SOLID96. Los elementos vigas son usados para crear idealizaciones unidimensionales de estructuras 3–D. ANSYS ofrece soluciones computacionales eficientes para este tipo de elementos comparados con elementos sólidos y placas. Cuando se selecciona el tipo de elemento, automáticamente quedan determinadas algunas propiedades como el conjunto de grados de libertad, lo cual también queda implícito al elegir la modalidad de análisis: estructural, térmico, magnético u eléctrico; y si el elemento elegido se encuentra en el espacio bidimensional o tridimensional. El modelo de comportamiento del material también puede elegirse desde un menú. Los modelos disponibles son estructural, térmico, electromagnético, acústico, de fluidos, piezoeléctricos.

Una vez que se ha definido el tipo y las opciones de análisis, el siguiente paso es aplicar las cargas. La palabra "carga" como la utiliza ANSYS incluye las condiciones de frontera (restricciones, apoyos o especificaciones de campos en la frontera), así como las fuerzas aplicadas externas o internas. Cuando la solución ha sido calculada, se pueden utilizar los post–procesadores de ANSYS para revisar los resultados. Puede obtener desplegados de contornos, configuraciones deformadas y listas tabuladas para revisar e interpretar los resultados del análisis. Incluye la estimación de error, combinaciones de cargas y cálculos entre resultados. El análisis estructural es probablemente la aplicación más común del método del elemento finito. El término estructural implica no solamente estructuras de ingeniería civil como los puentes y los edificios, ya que también abarca estructuras navales, aeronáuticas y mecánicas, tales como embarcaciones, cuerpos de aviones y casas de máquinas, así como los componentes mecánicos individuales como pistones, partes para maquinarias y herramientas. ANSYS puede llevar a cabo siete tipos de análisis estructural. Las primeras incógnitas calculadas (grados nodales de libertad) en un análisis estructural son los desplazamientos. Otras cantidades como las deformaciones, esfuerzos, y reacciones, se derivan de los desplazamientos nodales.

ANSYS utiliza el lenguaje paramétrico APDL que se utiliza para efectuar tareas automáticas o bien, construir el modelo en términos de parámetros (variables). El lenguaje APDL también tiene otras características que permiten evitar la repetición de comandos, hacer macros, utilizar cláusulas if–then–else, hacer ciclos y efectuar operaciones vectoriales y matriciales. Un conocimiento básico de APDL facilita el manejo del programa y permite hacer modelos más adecuados al fenómeno que se esté estudiando.

Los parámetros son las variables de APDL y son muy similares a las variables o parámetros utilizados en FORTRAN. ANSYS cuenta además con una guía de referencias teóricas cuyo propósito es dar a conocer al usuario de ANSYS las bases teóricas con que se construyó el programa. Entendiendo la teoría, el programa puede usarse más inteligentemente conociendo sus alcances y sus limitaciones. Por razones de espacio no se abunda más en este sentido, pero el lector podrá acceder a la ayuda disponible en ANSYS para consultar las bases teóricas o bien, consultar las referencias ahí especificadas.

Para este modelo en particular se utilizaron elementos viga de propiedades elásticas para el modelado (BEAM3). Este tipo de viga es un elemento con capacidades de tensión, compresión y flexión. El elemento tiene tres grados de libertad en cada nodo: desplazamientos nodales en las direcciones x e y, además de rotación respecto al eje z. Se consideró una temperatura constante en todo el cuerpo del elemento finito durante todo el proceso de carga y deformación.

Las funciones de forma para este elemento finito son:

Donde i e j identifican los nodos; u y v son los desplazamientos en las direcciones x e y, respectivamente, del sistema de coordenadas local para cada elemento. L es la longitud del elemento, s es la pendiente y θ_z es la rotación en la dirección z.

donde A es el área de la sección transversal; E es el módulo elástico; L es la longitud del elemento; e I es el momento de inercia

El factor de forma es , que toma en cuenta el efecto de la rigidez por cortante para vigas prismáticas de sección rectangular. G es el módulo de cortante; A^s es el área de cortante (A^s = A / F^s) y Fs es la constante de deflexión por cortante

Si bien el problema es enfocado como una armadura en 2D para el análisis de cargas y momentos, el análisis de las deflexiones en cada uno de los elementos toma en cuenta las secciones estructurales de las vigas. Esto es válido y es la manera en que se acostumbra realizar los análisis de los marcos estructurales que constituyen los edificios. En este caso se están manejando secciones cuadradas de 1×1 micrómetros de lado. Los elementos viga tienen 10 micrómetros de largo. Cada línea de nodo a nodo en el modelo es un elemento viga y es equivalente a un lazo de la red.

Para modelar las redes con diferentes porosidades se utilizó lenguaje de programación APDL. Las consideraciones para desarrollar el algoritmo fueron las siguientes: se construyeron las redes de lazos con un solo material con las propiedades elásticas de los granos de arcilla, al cual se llamó "Material 1". Se introdujo otro material con propiedades idénticas al que se denominó "Material 2" y que en el futuro representará a los huecos del suelo. Mediante un algoritmo para llevar a cabo una selección aleatoria, se escoge cierta cantidad de lazos de acuerdo a la p, los cuales se cambiaron de "Material 1" a "Material 2" (Figura 3). El "Material 2" constituye (1–p) y representa a la masa "muerta" o "hueca" del suelo–compósito. Todos los lazos son elementos viga dentro del modelo de elementos finitos y tienen propiedades geométricas idénticas.

Para minimizar la participación de la parte "muerta" en el modelo, se redujo su rigidez en 1×10^–8. Esto significa que el "Material 2" se vuelve una parte extremadamente "suave" dentro de una estructura rígida cuyo esqueleto resistente son todos los lazos de "Material 1". Los elementos "muertos" no pueden ser simplemente removidos de la red porque la estructura se vuelve inestable. La rigidez del "Material 2" se elige muy cercana a cero para que se comporte como una red de poros del suelo –aunque estructuralmente esos lazos existen– ya que de otra manera no podría llevarse al cabo el análisis. Estos lazos fueron removidos de las imágenes para hacer más fácil la interpretación. Una vez que fueron modeladas las dos partes del compósito, se procedió a realizar el análisis del sistema con el método del elemento finito. Las redes fueron sometidas a esfuerzos de compresión que variaron desde cero hasta el necesario para hacer fallar la estructura por deformación excesiva. Sólo los desplazamientos en los nodos inferiores fueron restringidos en todos los sentidos. Al resto de los nodos se les permitió el libre desplazamiento, tal como ocurre con los especímenes ensayados en el laboratorio. Los resultados fueron analizados mediante otra subrutina para obtener el valor medio de las deformaciones verticales, las deformaciones horizontales y la relación de Poisson. Todos los procesos de carga–deformación fueron analizados a fin de obtener el módulo elástico de cada una de las redes modeladas. Así, se obtuvo el módulo elástico medio para cada valor de porosidad en cada red.

La Figura 3 muestra un ejemplo para la estructuración de la red cuadrada. Los lazos en color negro representan los lazos ocupados con probabilidad p. Los lazos ausentes representan la porosidad de la red. Las restricciones nodales están representadas por flechas en los nodos de la base de la red, en la cual fueron restringidas las deformaciones en x e y. Las cargas externas verticales sólo fueron aplicadas en la parte superior de la red. El origen del plano de referencia se encuentra en la esquina inferior izquierda. La Figura 4 es equivalente a la Figura 3, pero para el caso de una red hexagonal de lazos.

Las Figuras 5 y 6 muestran las redes deformadas después de que se les han aplicado las cargas verticales externas, rectangular y hexagonal, respectivamente y una vez que ha sido hecho el análisis con el FEM.

Se realizaron gráficas de los desplazamientos nodales en el sentido vertical para localizar los puntos con mayores deformaciones, las cuales se ilustran en las Figuras 7 y 8. Estas gráficas fueron escogidas de entre todas las existentes para ilustrar el procedimiento de análisis que se consigue con el modelo. En cada una de las redes y para diferentes valores de p se hicieron múltiples corridas, con el fin de garantizar la representatividad de los resultados obtenidos con el modelo.

Hecho así el planteamiento teórico, tomando en cuenta las propiedades de la arcilla y tratando de representar su estructura mediante redes sencillas como la cuadrada y la hexagonal, se procedió a realizar todo el trabajo de cómputo para modelar el comportamiento mecánico.

RESULTADOS

Los resultados obtenidos se pueden observar en la Figura 9. En esta figura se muestran los valores del módulo elástico medio obtenido para los dos tipos de redes que se modelaron: rectangular y panaloide, en comparación con los valores medidos en el laboratorio en pruebas de micromecánica.

Los datos experimentales obtenidos en el laboratorio fueron determinados en los intervalos de porosidad que se pueden fabricar con los equipos disponibles. Estas compacidades existen en la naturaleza y pueden ser reproducidas en el laboratorio (Pérez–Rea, 2005). Los modelos para las redes cuadrada y hexagonal fueron llevados más allá de los límites de porosidades máxima y mínima que se lograron en especímenes de suelo compactado en el laboratorio, para que los valores de p_cr fueran más fácilmente identificables en los gráficos.

El valor de p_cr determinado experimentalmente fue de 0.3943, menor que el umbral de percolación geométrico para una red cuadrada en 2D y mucho menor que el umbral de percolación geométrico en una red hexagonal también en 2D. La porosidad del suelo viene representada por el término (1–p), siendo p la porción de la red ocupada por los sólidos del suelo. (1–p_cr) representa la porosidad crítica del suelo para la rigidez mínima (E).

Un aspecto importante del espacio poroso es la interconectividad en la muestra. Si se observaran imágenes de la muestra en 2D y 3D, necesariamente la interconectividad será diferente. Las representaciones en 2D son una proyección de una red conectada donde algunos poros aparecen "no conectados"; de manera que se necesitará una mayor probabilidad de ocupación para encontrar racimos percolantes en 2D que en 3D. Lazos que están conectados en 3D aparecen como no conectados en 2D. Sin embargo, la distribución de tamaños de poros es de valores muy similares en 2D y 3D (Moreau et al. 1999).

CONCLUSIONES

Diferentes modelos de redes de percolación se han utilizado en diferentes áreas científicas y su uso ha sido aplicado para la predicción de propiedades de diferentes materiales. Aunque estos métodos aparecieron en la década de 1950, su uso en ciertas áreas de la ciencia como la Mecánica de Suelos es muy reciente y escaso, de acuerdo con la bibliografía existente. Hasta ahora, los modelos de percolación han sido únicamente aplicados a medios porosos para modelar propiedades hidráulicas y de conductividad eléctrica; sin embargo, pueden tener un gran potencial en el campo del modelado de comportamiento esfuerzo–deformación de materiales.

En este trabajo se realizaron modelos de redes regulares (cuadrada y hexagonal) para representar la microestructura de suelos arcillosos con el objeto de hacer un análisis del comportamiento esfuerzo–deformación y predecir así el módulo elástico que tendrá el suelo macroscópicamente. Los análisis mediante elemento finito de las redes arrojaron resultados muy cercanos a los determinados experimentalmente en suelos remoldeados. El módulo elástico del suelo es un parámetro sumamente importante en el diseño de cimentaciones de estructuras civiles y la dificultad para medirlo en campo hace que tales diseños se efectúen con valores deducidos de pruebas empíricas.

La estructura de un medio poroso depende de su heterogeneidad y la escala de longitud a la cual el medio es inspeccionado. La caracterización aquí mencionada se refiere a medios porosos macroscópicamente homogéneos (aunque microscópicamente sean desordenados) y cuyas propiedades son, para una muestra lo suficientemente grande, independientes de su dimensión lineal.

La universalidad del comportamiento crítico ha propiciado la aplicación con éxito de la teoría de percolación a un número de problemas que parecían estar muy lejos de ser solucionados, pero trabajos recientes han demostrado que en ciertos casos pueden ser muy acertados.

Quizás los modelos en 2D no sean todavía muy acertados, pero sí proporcionan un acercamiento suficiente al problema. En trabajos posteriores, modelos en 3D y configuraciones aleatorias podrían proporcionar mejores aproximaciones. La teoría de percolación proporciona una herramienta para la interpretación de datos experimentales y la comprensión de la estructura de medios porosos, como suelos de diversos tipos. La combinación de la teoría de percolación con herramientas numéricas poderosas, como el método del elemento finito, permite modelar el comportamiento de materiales que, si bien no son macroscópicamente continuos, pueden ser analizados como un medio continuo a nivel microscópico.

Barton, C.C., Hsieh, P.A., 1989, Physical and hydrological flow properties of fractures, en 28th International Geological Congress, Field trip Guidebook T385: Washington D.C, American Geophysic Union, 36 p.        [ Links ]

Bolander, J.E., Saitot, S., 1997, Discrete modeling of short–fiber reinforcement in cementitious composites: Advanced Cement Basic Materials, 6, 76–86.         [ Links ]

Broadbent, S.R., Hammersley, J.M., 1957, Percolation processes I. Crystals and mazes: Proceedings of Cambridge Philosophy Society, 53, 629–641.        [ Links ]

Clough, R.W., 1960, The finite element method in plane stress analysis, in Proceedings of the 2nd Conference on Electronic Computation: Pittsburgh, Pennsylvania, American Society of Civil Engineers, 23, 345–378.        [ Links ]

Collis, K., McGown, A., 1974, The form and function of microfabric features in a variety of natural soils: Geotechnique, 24, 223–254.         [ Links ]

Courant, R., 1943, Variational methods for the solution of problems of equilibrium and vibrations: Bulletin of the American Mathematical Society 49, 1–23.        [ Links ]

Chatzis, I., Dullien, F.A.L., 1977, Modelling pore structure by 2–D and 3–D networks with application to sandstones: Journal of Canadian Petroleum Technology, 16, 97–108.        [ Links ]

Davis, R.O., Selvadurai, A.P.S., 1996, Elasticity and Geomechanics: New York, USA, Cambridge University Press, 256 p.        [ Links ]

De Gennes, P.G., 1976, On a relation between percolation theory and the elasticity of gels: Journal of Physique (Paris), Letter 37, L1–L2.         [ Links ]

Deshpande, V.S., Fleck, N.A., Ashby, M.F., 2001, Effective properties of the octet–truss lattice material: Journal of the Mechanics and Physics of Solids, 49, 1747–1769.        [ Links ]

Flory, P.J., 1941, Molecular size distribution in three dimensional polymers. I. Gelation: Journal of American Chemical Society, 63, 3083–3090.         [ Links ]

Guéguen, Y., Chelidze, T., Le Ravalec, M., 1997, Microestructures, percolation thresholds, and rock physical properties: Tectonophysics, 279, 23–35.        [ Links ]

Hrenikoff, A., 1941, Solution of problems in elasticity by the framework method: Journal of Applied Mechanics, 8, 169–175.        [ Links ]

Isichenko, M.B., 1992, Percolation, statistical topography, and transport in random media: Review of Modern Physics, American Physical Society, 64 (4), 961–1043.        [ Links ]

Juárez–Badillo, E., Rico–Rodríguez, A., 1984, Mecánica de Suelos, Tomo I: México, Limusa, 642 p.        [ Links ]

Kantor, Y., Webman, I., 1984, Elastic properties of random percolating systems: Physical Review Letters, 52(21), 1891–1894.        [ Links ]

Keedwell, M.J., 1984, Reology and Soil Mechanics: New York, USA, Elsevier Applied Science Publishers, 340 p.        [ Links ]

Lee, J., 2006, Synthesis and characterization of polyethylene–based nanocomposites: Ontario, Canadá, Queen's University, Ph.D thesis.         [ Links ]

Liem, C., Jan, N., 1988, Fractal properties of the percolating backbone in three dimensions: Journal of Physics A, Mathematical and General, 21, L243–L245.        [ Links ]

Mitchell, J.K., Geenberg, J.A., Witherspoon, P.A., 1973, Chemico–osmotic effects in fine–grained soils: Journal of the Soil Mechanics and Foundation Division, American Society of Civil Enginneers (ASCE), 99(4), 307–322         [ Links ]

Moreau, E., Velde, B., Terribile, F., 1999, Comparison of 2D and 3D images of fractures in a vertisol: Geoderma, 92, 55–72.        [ Links ]

Murad, M.A., Moyne, C., 2002, Micromechanical computational modeling of expansive porous media: Comptes Rendus Mecanique, 330, 865–870.        [ Links ]

Pabst, W., Gregorová, E., Tichá, G., 2007, Effective properties of suspensions, composites and porous materials: Journal of the European Ceramic Society, 27, 479–482.        [ Links ]

Pérez–Rea, M.L., 2005, Predicción de propiedades mecánicas y reológicas de suelos usando teoría de percolación: Atizapán, Estado de México, México, Instituto Tecnológico y de Estudios Superiores de Monterrey, Campus Estado de México, Tesis Doctoral, 136 pp.        [ Links ]

Pérez–Rea, M.L., Rojas, F., Castaño, V.M., 2003, Nitrogen adsorption method for determining the meso slit–like size distribution of expansive soils: Materials Research Innovations, 7, 341–352.         [ Links ]

Poulos, H.G., Davis, E.H., 1974, Elastic Solutions for Soil and Rock Mechanics: New York, USA, John Wiley & Sons, 424 p.        [ Links ]

Przemieniecki, J.S., 1968, Theory of Matrix Structural Analysis: New York, USA, McGraw Hill, 480 p.        [ Links ]

Rastellini, F., Oller, S., Salomón, O., Oñate, E., 2008, Composite materials non–linear modelling for long fibre–reinforced laminates. Continuum basis, computational aspects and validations: Computers and Structures, 86, 879–896.        [ Links ]

Sahimi, M., 1994, Applications of Percolation Theory: Bristol Pennsylvania, USA, Taylor & Francis, 276 pp.         [ Links ]

Sahimi, M., 1995, Flow and Transport in Porous Media and Fractured Rock: Germany, VCH, 482 pp.         [ Links ]

Sahimi, M., Gavalas, G.R. Tsotsis, T.T., 1990, Statistical and continuum models of fluid solid reactions in porous media: Chemical Engineering Science, 45, 1443–1502.        [ Links ]

Shan, Z., Gokhale, A.M., 2002, Representative volume element for non–uniform micro–structure: Computational Materials Science, 24, 361–379.        [ Links ]

Scher, H., Zallen, R., 1970, Critical density in percolation processes: Journal of Chemical Physics, 53(9), 3759 – 3761         [ Links ]

Schmidt, V., McDonald, D.A., 1979, Texture and recognition of secondary porosity in sandstones, en Scholle, P.A., Schluger, P.R. (eds.), Aspects of Diagenesis: Society of Economic Paleontologists and Mineralogists, Special Publication 26, 209–225.        [ Links ]

Stauffer, D., Coniglio, A., Adam, M., 1982, Gelation and critical phenomena: Advanced Polymer Science, 44, 103–158.        [ Links ]

Stockmayer, W.H., 1943, Theory of molecular size distribution and gel formation in branched–chain plymers: Journal of Chemical Physical, 11(2), 45–55.        [ Links ]

Taya, M., 1999, Micromechanics modeling of smart composites: Composites, Part A, Applied science and manufacturing 30, 531–536.        [ Links ]

Turner, M.J., Clough, R.W., Martin, H.C., Topp, L.J., 1956, Stiffness and deflection analysis of complex structures: Journal of Aeronautical Science 23(9), 805–824.        [ Links ] ]]> 1989 Washington D.C 36 1997 6 76-86 1957 53 629-641 1960 23 Pittsburgh Pennsylvania 345-378 1974 24 223-254 1943 49 1-23 1977 16 97-108 1996 256 1976 37 2001 49 1747-1769 1941 63 3083-3090 1997 279 23-35 1941 8 169-175 1992 64 4 4 961-1043 1984 I 642 1984 52 21 21 1891-1894 1984 340 2006 1988 21 ^sMathematical and General ^sMathematical and General Mathematical and General 1973 99 4 4 307-322 1999 92 55-72 2002 330 865-870 2007 27 479-482 2005 136 2003 7 341-352 1974 424 1968 480 2008 86 879-896 1994 276 1995 482 1990 45 1443-1502 2002 24 361-379 1970 53 9 9 3759 - 3761 1979 26 209-225 1982 44 103-158 1943 11 2 2 45-55 1999 30 531-536 1956 23 9 9 805-824