Una formulación general y eficiente de las fracturas en el MEF: II
Contents
1. 6 7 8 9 10 Figura 17 Mapa de concentraciones a tiempo T 2000 d as para el Modelo 2D con 18149 nudos M2D_20K Figura 18 Mapa de concentraciones a tiempo T 10000 d as para el Modelo 2D con 18149 nudos M2D_20K O i 2 3 4 5 6 7 8 9 10 Figura 19 Mapa de concentraciones a tiempo T 2000 d as para el Modelo 1D con 2201 nudos M1D_2K Figura 20 Mapa de concentraciones a tiempo T 10000 d as para el Modelo 1D con 2201 nudos M1D_2K 74 R Juanes y J Samper 0 1 2 3 4 5 6 7 8 10 Figura 21 el Modelo 1D con 10 Figura 22 Mapa de concentraciones a tiempo T 10000 d as para el Modelo 1D con 6501 nudos M1D_6K yk 0 1 2 3 4 5 6 7 8 9 10 Figura 23 Mapa de concentraciones a tiempo T 2000 d as para el Modelo 1D con 19705 nudos Figura 24 Mapa de concentraciones a tiempo T 10000 d as para el Modelo 1D con 19705 nudos M1D_20K Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 75 de solutos es necesario recurrir a una discretizaci n muy fina para representar de forma apropiada la heterogeneidad del medio La explicaci n matem tica de este efecto es objeto de futura investigaci n Los resultados obtenidos con el Modelo 1D muestran un comportamiento muy distinto Con una malla de 2201 nudos M1D_2K se obtienen los resultados que se muestran en las Figuras 19 y 20 mien2. Water Resour Res Vol 22 N 3 pp 361 374 1986 82 R Juanes y J Samper 10 11 12 13 14 15 16 G Gambolati F Sartoretto y F Uliana A Conjugate Gradient Finite Element Model for Large Multiaquifer Systemsp Water Resour Res Vol 22 N 7 pp 1003 1015 1986 A Medina G Galarza y J Carrera TRANSIN II Simulaci n y estimaci n de par metros en flujo y transporte en medios porosos Manual del usuario Barcelona Departamento de Ingenier a del Terreno ETSICCP UPC pp 118 1993 A Larabi y F de Smedt Solving three dimensional hexahedral finite element groundwater models by preconditioned conjugate gradient methods Water Resour Res Vol 30 N 2 pp 509 521 1994 ENRESA Empresa Nacional de Residuos Radiactivos S A El Berrocal Project Characterization and validation of natural radionuclide migration processes under real conditions on the fissured granitic environment European Commission Contract N F12W CT91 0080 Vol 4 1996 SKB Swedish Nuclear Fuel and Waste Management Co sp Hard Rock Laboratory Annual Report 1995 Technical Report TR 96 06 1996 G A Galarza A Medina y J Carrera TRANSIN IIT Applications to 3D media and non linear problems El Berrocal Project Vol 4 Hydrogeological modelling and code development s l ENRESA pp 423 680 1996 J Molinero J Samper R Juanes L Buj n y G Zhang N
3. e 2001 Revista Internacional de M todos Num ricos para C lculo y Dise o en Ingenier a Una formulaci n general y eficiente de las fracturas en el MEF IT Aplicaci n a casos sint ticos Ruben Juanes y Javier Samper Grupo de Hidrolog a Subterr nea ETS de Ingenieros de Caminos Canales y Puertos Universidad de La Coru a Campus de Elvi a 15192 La Coru a Espa a Tel 34 981 167 000 Fax 34 981 167 170 e mail juanesQOce berkeley edu e mail samperQiccp udc es Resumen Este trabajo forma parte de una serie de dos art culos en los que se presenta una formulaci n num rica eficiente y completamente general para el tratamiento de fracturas embebidas en un medio poroso y de las condiciones de contorno en el m todo de elementos finitos En el primero de ellos se presentan las expresiones para el c lculo de las derivadas Cartesianas y para la evaluaci n de las integrales sobre hipersuperficies m dimensionales en espacios Euclidianos n dimensionales Este tratamiento conduce a una formulaci n compacta que es aplicable a la integraci n num rica en l neas superficies y vol menes en dominios tridimensionales evit ndose de esta forma los c lculos farragosos de la formulaci n tradicional En este segundo art culo se presenta la aplicaci n de esta formulaci n a tres casos sint ticos de flujo de agua y transporte de solutos a trav s de medios porosos y fracturados 1 bloque de roca rectangular con una fractur
4. a dispersividad longitudinal en la fractura a 0 025 m para evitar un n mero de Pecletl excesivamente elevado que podr a ocasionar oscilaciones esp reas Por contra la soluci n anal tica considera que en la fractura s lo existe advecci n siendo nulo el transporte dispersivo t El n mero de Peclet es un n mero adimensional que representa el cociente entre el transporte advectivo y el transporte difusivo dispersivo En el caso muy frecuente de que la difusi n molecular sea despreciable frente a la dispersi n mec nica el n mero de Peclet de la malla se expresa como Pe Az a donde Az es la longitud del elemento en la direcci n del flujo y aL es la dispersividad longitudinal 68 R Juanes y J Samper X Analitica 0 9 A Modelo 2D O Modelo 1D CONCENTRACION mg L o o p o o po N w A 4 o N e a 5 eS 0 X m Figura 6 Evoluci n de las concentraciones en la fractura Comparaci n de la soluci n anal tica l nea continua Modelo 2D cruces y Modelo 1D c rculos CASO 2 BLOQUE RECTANGULAR DE ROCA CON FRACTURAS ORTOGONALES INTERCONECTADAS Descripci n del problema El segundo caso considera el flujo de agua y el transporte de un soluto a trav s de un bloque de roca que contiene una red de fracturas ortogonales conectadas entre s tal como se muestra en la Figura 7 Al igual que en el caso anterior se considera que los contornos superior e inferior son impermeables En los c
5. 10m Figura 7 Esquema del problema para el Caso 2 Par metro Descripci n Matriz Fractura e espesor Im 1m 2b apertura 0 04m K conductividad hidr ulica 1075 md 107 md S coef almacenamiento espec f 10m 10 m gt O porosidad 0 01 0 1 aL dispersividad longitudinal 107 m 5 x 107 m ar dispersividad transversal 2 x 10m 107 m D coef difusi n efectivo 107 m2 d7 107 m d Tabla I Par metros de flujo y transporte del Caso 2 Para ambos modelos se han realizado simulaciones con distintas mallas de elementos En la Tabla II se recogen las caracter sticas fundamentales de las mallas empleadas as como la nomenclatura utilizada posteriormente en la discusi n de resultados finitos progresivamente m s refinadas Identificador N mero de N mero de Tipo de Ancho de Ancho de de la malla nudos elementos elementos banda inicial banda final M2D_2K 2112 2015 2D lineales 451 53 M2D_6K 6105 5940 2D lineales 5664 90 M2D_20K 18149 5940 2D cuadr ticos 17067 269 M1D_2K 2201 2314 2D y 1D lineales 538 51 M1D_6K 6501 2314 2D y 1D cuadr ticos 4802 152 M1D_20K 19705 6962 2D y 1D cuadr ticos 18830 292 Tabla II Caracter sticas de las mallas de elementos finitos empleadas para las distintas simulaciones del Caso 2 70 R Juanes y J Samper En la Figura 8 se muestra la malla de elementos finitos correspondiente al Modelo 2D que contiene 2112 nudos Se puede observar
6. la Ingenier a Civil de Galicia Proyectos FEBEX c d 703231 y de Validaci n de C digos c d 703334 El proyecto FEBEX ha sido financiado en parte por la Comisi n Europea a trav s del Programa de Fisi n Nuclear Proyecto F14W C T95 0006 Tambi n se ha contado con financiaci n de un proyecto de investigaci n de la CICYT del Proyecto Nacional de Recursos H dricos Proyecto HID98 282 REFERENCIAS 1 R Juanes Un c digo para la modelizaci n tridimensional de flujo y transporte Proyecto T cnico ETSI Caminos Universidad de La Coru a 1997 2 R Juanes y J Samper Una formulaci n general y eficiente de fracturas y condiciones de contorno en el MEF I Aspectos te ricos Revista Internacional de M todos Num ricos para C lculo y Dise o en Ingenier a Vol 16 N 4 pp 471 491 2000 3 G Segol A Three dimensional Galerkin Finite Element Model for the analysis of contaminant transport in saturated unsaturated porous media en Finite Elements in Water Resources W G Gray G F Pinder y C A Brebbia Eds Pentech London pp 2123 2144 1977 4 P S Huyakorn B H Lester y C R Faust Finite element techniques for modeling groundwater flow in factured aquifers Water Resour Res Vol 19 N 4 pp 1019 1035 1983 5 P S Huyakorn B G Jones y P F Andersen Finite element algorithms for simulating three dimensional groundwater flow and solute transport in multilayer systems
7. medio Por los motivos expuestos en el primer art culo se ha optado por un modelo de fracturas embebidas en el medio poroso que resuelve simult neamente las ecuaciones de flujo y transporte en el medio poroso y en las fracturas Aunque este modelo est especialmente dise ado para la representaci n de zonas de fractura de gran extensi n superficial y grandes fracturas regionales tambi n es aplicable a fracturas de menor escala y a medios fracturados generados estoc sticamente Asimismo la formulaci n matem tica es completamente general de modo que es aplicable a sistemas de flujo bidimensional y tridimensional con fracturas 1D y 2D que se intersecan entre s Es factible la aplicaci n de este modelo a problemas definidos por otras ecuaciones como las de la elasticidad la elastoplasticidad y campos acoplados por ejemplo temperatura y deformaci n Primeramente se presenta una breve descripci n del programa TRANMEF 3 un c digo en elementos finitos para la resoluci n de las ecuaciones de flujo subterr neo transporte de varios solutos reactivos y transporte de calor a trav s de medios porosos y fracturados A continuaci n se muestran varios ejemplos de aplicaci n del c digo a tres casos sint ticos que ilustran el potencial la aplicabilidad y las ventajas num ricas de la formulaci n propuesta para el tratamiento de fracturas y contornos EL C DIGO TRANMEF 3 TRANMEF 3 es un programa de c lculo que resuelve la ecuaci n
8. que la malla se ha refinado en el entorno de las fracturas debido a la mayor velocidad del agua a lo largo de las zonas fracturadas y al cambio dr stico de las propiedades existente en las interfaces entre las fracturas y la matriz rocosa Las mallas para el resto de los casos son similares y se han procesado con los mismos criterios Se ha empleado la misma discretizaci n temporal para todas las simulaciones basada en un esquema de Crank Nicolson 0 0 5 con un incremento de tiempo que aumenta progresivamente desde At 0 1 d hasta At 10d Todas las simulaciones se han realizado en el ordenador DIGITAL AlphaServer 4000 5 466 4MB bajo sistema operativo OpenVMS V7 1 del Centro de C lculo de la Escuela de Ingenieros de Caminos de la Universidad de La Coru a Figura 8 Malla de elementos finitos del Caso 2 para el Modelo 2D con 2112 nudos Resultados En este apartado se presentan los resultados obtenidos en las seis simulaciones realizadas En primer lugar se comentan los resultados de flujo Por su mayor relevancia se dedica un an lisis m s detallado a los resultados del transporte de solutos Las soluciones num ricas de la ecuaci n de flujo estacionario obtenida con el Modelo 2D con las mallas M2D_2K y M2D_6K se muestran en las Figuras 9 y 10 respectivamente No se presenta la soluci n obtenida en la pasada M2D_20K por ser pr cticamente igual a la de la simulaci n M2D_6K Las peque as diferencias entre las Figuras 9 y 10 indi
9. y 18 Este hecho indica que cuando se emplea el Modelo 2D para la resoluci n de la ecuaci n de transporte Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 71 0 1 2 3 4 5 6 7 Figura 9 Mapa de niveles piezom tricos flujo estacionario para el Modelo 2D con 2112 nudos M2D_2K 7 Figura 10 es de niveles piezom tricos n P para el Modelo 2D con 6105 nudos M2D_6K L E A o y 7 8 Figura 11 Mapa de niveles piezom tricos flujo estacionario para el Modelo 1D con 2201 nudos M1D_2K CIN Figura 12 Mapa de niveles piezom tricos flujo estacionario para el Modelo 1D con 6501 nudos M1D_6K y m oF 72 R Juanes y J Samper 0 1 2 3 4 5 lt o 10 Figura 13 Mapa de concentraciones a tiempo T 2000 d as para el Modelo 2D con 2112 nudos M2D_2K 0 1 2 3 4 5 6 7 8 9 10 Figura 14 Mapa de concentraciones a tiempo T 10000 d as para el Modelo 2D con 2112 nudos M2D_2K 4 AE r i i E ER j 0 1 2 3 4 5 6 7 8 9 10 Mapa de concentraciones a tiempo T 2000 d as para el Modelo 2D con 6105 nudos M2D_6K Figura 15 7 9 10 Mapa de concentraciones a tiempo T 10000 d as para el Modelo 2D con 6105 nudos M2D_6K Figura 16 Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 73
10. I Aplicaci n a casos sint ticos TT El mapa de errores para el Modelo 2D con una malla de 6105 nudos se muestra en la Figura 27 Las caracter sticas son muy parecidas a las de la pasada anterior con errores generalizados que reproducen el patr n de las isol neas de concentraci n comp rense las Figuras 16 y 27 y con un error m ximo cercano al 25 En la Figura 28 se representan los errores correspondientes al Modelo 2D con una malla de 18149 nudos Son aplicables los mismos comentarios que para las pasadas anteriores aunque el error m ximo es en este caso del 8 La Figura 29 muestra los errores en concentraciones de la soluci n num rica obtenida con el Modelo 1D con una malla de 2201 nudos Se aprecia un comportamiento completamente diferente al de las simulaciones con el Modelo 2D El error m ximo es tan s lo del 2 y los errores se encuentran localizados en el entorno de las fracturas a pesar de que la discretizaci n empleada es relativamente grosera Es evidente que estos resultados son de una calidad excepcional especialmente si se comparan con los que proporciona el Modelo 2D para un coste computacional similar M2D_2K el error m ximo es 20 veces menor y los errores se concentran en el entorno de las fracturas Los errores correspondientes al Modelo 1D con 6501 nudos son virtualmente nulos En la Tabla III se recogen a modo de conclusi n los resultados m s relevantes de las distintas simulaciones realizadas De estos res
11. a 2 bloque de roca rectangular con una red de fracturas ortogonales interconectadas y 3 bloque de roca tridimensional con una fractura vertical Las simulaciones realizadas ilustran claramente el potencial y la aplicabilidad de la formulaci n propuesta y ponen de manifiesto que el tratamiento de las fracturas con elementos de menor dimensi n que el dominio por ejemplo elementos 1D en medios 2D redunda en soluciones num ricas m s precisas y menos costosas que las obtenidas con la formulaci n tradicional de elementos finitos A GENERAL AND EFFICIENT FORMULATION FOR FRACTURES IN THE F E M II APPLICATION TO SYNTHETIC CASES Summary In this series of two papers an efficient numerical formulation for the implementation of boundaries and fractures embedded in porous media by the finite element method is presented In the first of these papers the expressions for the computation of Cartesian derivatives and integrals over m dimensional surfaces within n dimensional Euclidean spaces are derived This leads to a compact formulation which allows to integrate numerically over lines surfaces and volumes in three dimensional domains avoiding the convoluted computations of classical formulations The second paper includes three synthetic cases of solute transport through porous and fractured media 1 rectangular block with a fracture 2 rectangular block with a network of connected orthogonal fractures and 3 three dimensional block with a vertical f
12. ak y Pickens asume que la fractura tiene porosidad unidad f 1 66 R Juanes y J Samper donde q es la componente horizontal de la velocidad de Darcy e es el gradiente de nivel piezom trico en la direcci n zx En la Figura 2 se muestra la malla de elementos finitos empleada en este caso En el Modelo 2D la fractura se ha discretizado con elementos cuadrangulares lineales 2 elementos en la direcci n y para un total de 1 100 elementos y 1173 nudos En el Modelo 1D la fractura se discretiza con elementos unidimensionales lineales para un total de 1050 elementos y 1071 nudos Ambas mallas coinciden excepto en la discretizaci n de la fractura Para la discretizaci n temporal se emplea en todos los casos un esquema de Crank Nicolson factor de ponderaci n temporal 0 0 5 y un incremento de tiempo que aumenta progresivamente desde At 0 1 d hasta At 20 d Figura 2 Malla de elementos finitos del Caso 1 para el Modelo 2D fractura discretizada con elementos 2D Resultados En la Figura 3 se muestran las concentraciones en la matriz al cabo de 1000 d as obtenidas a partir soluci n anal tica dada por la Ecuaci n 7 Las calculadas con los Modelos 2D y 1D se representan en las Figuras 4 y 5 respectivamente Las nicas discrepancias entre la soluci n anal tica y las num ricas se encuentran localizadas en el contorno izquierdo debido a que en el modelo num rico existe una peque a entrada de agua por dicho contorno q
13. ar metro Descripci n Matriz Fractura Pozo e espesor 0 1m A rea transversal 0 01 m K conductividad hidr ulica 107 md 107 md 10 md S5 coef almacenamiento espec f 1076 m t 107 m 1m Q porosidad 0 01 0 1 1 QL dispersividad longitudinal 0 5m 1m 1m aT dispersividad transversal 0 1m 0 2 m 1m D coef difusi n efectivo 107 m d 110 m 4 1074m dt Tabla IV Par metros de flujo y transporte del Caso 3 Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 79 La Figura 31 muestra la malla de elementos finitos empleada que consta de 16410 elementos y 18491 nudos Como se puede apreciar la malla se ha refinado en el entorno de la fractura ya que es en esa zona donde se producen los cambios m s bruscos en los par metros de flujo y transporte en el campo de velocidades y los mayores gradientes de concentraciones En mallas de este tipo en las que el n mero de nudos es elevado es crucial realizar una renumeraci n de los nudos de la malla que permite reducir dr sticamente el coste computacional En este caso el m ximo semiancho de banda se reduce de 11072 a 739 Figura 31 Malla de elementos finitos del Caso 3 Resultados La distribuci n espacial tridimensional de las concentraciones de soluto al cabo de 200 1000 y 5000 d as se muestra en las Figuras 32 33 y 34 respectivamente La evoluci n del transporte del so
14. can que es necesaria una malla del orden de los 6000 nudos para obtener una soluci n precisa de la ecuaci n de flujo cuando se emplea el Modelo 2D Las Figuras 11 y 12 muestran la soluci n de la ecuaci n de flujo obtenida con el Modelo 1D con las mallas M1D_2K y M1D_6K respectivamente Los valores num ricos de ambas pasadas difieren en la cuarta cifra significativa de modo que una malla de 2000 nudos es suficiente para este caso Las simulaciones M2D_6K y M1D_6K proporcionan pr cticamente los mismos resultados Las Figuras 13 a 24 muestran los mapas de concentraciones correspondientes a las seis pasadas al cabo de 2000 y 10000 d as Las concentraciones obtenidas en la pasada M2D_2K con el Modelo 2D y una malla de 2 112 nudos Figuras 13 y 14 difieren claramente de las obtenidas con el mismo modelo pero con una malla de 6105 nudos v anse las Figuras 15 y 16 Estas diferencias ponen de manifiesto que la discretizaci n empleada no es suficiente para obtener una soluci n num rica aceptable Este efecto persiste incluso con una discretizaci n todav a m s fina Las Figuras 17 y 18 ilustran que los resultados correspondientes al Modelo 2D con una malla de 18149 nudos M2D_20K son visualmente distintos a los de las Figuras 15 y 16 Adem s la tendencia no es err tica cuanto m s fina es la discretizaci n espacial mayor es la concentraci n de soluto en todo el dominio para un tiempo dado se pueden comparar por ejemplo las Figuras 14 16
15. del flujo subterr neo la del transporte de varios solutos y la del transporte de calor en medios porosos y fracturados completamente tridimensionales por el M todo de Elementos Finitos El c digo desarrollado tiene el m rito de ampliar el rango de aplicabilidad de los c digos comerciales y de investigaci n actualmente disponibles No existe restricci n alguna respecto al dominio del problema pudiendo ser 1 Unidimensional Puede ser conveniente para tratar problemas puramente unidimensionales o bidimensionales con simetr a radial 2 Bidimensional horizontal Este ha sido el dominio tradicionalmente empleado para el flujo en acu feros en los que la variaci n de las variables de estado en la direcci n vertical es despreciable frente a la variaci n en la direcci n horizontal 3 Bidimensional vertical Se emplea cuando la variaci n de las variables de estudio es mucho mayor en la direcci n vertical que en la lateral 4 Cuasi tridimensional Generalmente consiste en un dominio que se discretiza en diversas capas bidimensionales horizontales que representan los acu feros conectados entre s mediante elementos unidimensionales verticales que representan los acuitardos Se puede emplear para simular medios fracturados o como simplificaci n de medios tridimensionales Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 63 5 Tridimensional Tratamiento completamente 3D S
16. e de soluci n anal tica y ha sido estudiado con anterioridad por Grisak y Pickens Por otra parte ilustra el papel de la difusi n en la matriz en el transporte de solutos en medios fracturados La difusi n en la matriz es el mecanismo de transporte por el cual existe una transferencia de soluto desde las zonas permeables con flujo de agua fracturas hacia las zonas donde apenas hay flujo de agua matriz Carrera et al proporcionan una exhaustiva revisi n sobre distintas formulaciones y efectos cualitativos del proceso de difusi n en la matriz El dominio de flujo consiste en un medio poroso rectangular con una fractura en el sentido longitudinal tal y como se muestra en la Figura 1 Se considera flujo de agua de izquierda a derecha a lo largo de la fractura cuya conductividad hidr ulica es mucho mayor que la de la matriz Inicialmente la concentraci n del soluto es nula en todo el dominio El agua entra de forma continua en la fractura con una concentraci n Co Se supone que el soluto se transporta en la fractura exclusivamente por advecci n y que se difunde en la matriz s lo en la direcci n y e 4m i 0 04 m Advecci n en la fractura 10 m 1 Figura 1 Esquema del Caso 1 bloque de roca rectangular con una fractura en su contorno inferior Formulaci n matem tica y soluci n anal tica La formulaci n matem tica del problema se basa en las siguientes condiciones La concentraci n en la fractura s l
17. e debe emplear para representar adecuadamente situaciones en que el flujo y o transporte es claramente tridimensional En la bibliograf a existen ejemplos sint ticos y casos reales de este tipo de situaciones 6 Multidimensional Se denomina as a la combinaci n de los dominios anteriores para representar formaciones tridimensionales que contengan estructuras unidimensionales o bidimensionales embebidas en ellas Tiene su campo de aplicaci n m s importante en la simulaci n del flujo y transporte en medios fracturados En la bibliograf a se describe un c digo con caracter sticas similares a TRANMEF 3 y aplicable a este tipo de dominios El usuario puede emplear hasta diez tipos de elementos diferentes que pueden emplearse conjuntamente en cualquier simulaci n Estos elementos pueden no ser necesariamente de la misma dimensi n pueden ser 1D 2D 3D lo cual permite simular redes de fracuras bidimensionales y sus intersecciones unidimensionales en un medio tridimensional El programa realiza un tratamiento riguroso de las condiciones de borde mediante una integraci n completa sobre los contornos del dominio El r gimen de flujo y transporte puede ser estacionario o transitorio Las condiciones iniciales se pueden definir expl citamente o consistir en el resultado de una simulaci n estacionaria El programa resuelve las siguientes ecuaciones 1 Flujo subterr neo Se resuelve la ecuaci n de flujo en su forma lineal Esto im
18. ente para el tratamiento de las fracturas descrita en detalle en el primer art culo El primer caso ha servido para verificar la formulaci n propuesta mediante su comparaci n con la formulaci n tradicional y la soluci n anal tica Con el segundo caso se ha demostrado que el Modelo 1D fracturas discretizadas con elementos unidimensionales proporciona resultados mucho m s precisos que el Modelo 2D fracturas discretizadas con elementos bidimensionales para la soluci n num rica de la ecuaci n de transporte de solutos En particular la soluci n que se obtiene con el Modelo 1D con una malla de 2201 nudos presenta un error m ximo del 2 frente al 50 del Modelo 2D con 2 112 nudos y al 8 con 18 149 nudos Desde un punto de vista de eficiencia num rica la formulaci n propuesta permite obtener soluciones m s precisas con un coste computacional memoria utilizada y tiempo de CPU varios rdenes de magnitud inferior ver Tabla III Con el tercer caso se ha ilustrado la generalidad y flexibilidad de la formulaci n propuesta para la integraci n en fracturas y condiciones de contorno mostrando su aplicaci n a la simulaci n de un ensayo de trazadores en un medio tridimensional fracturado AGRADECIMIENTOS Este trabajo se ha realizado en el marco de proyectos de investigaci n financiados por la Empresa Nacional de Residuos Radiactivos ENRESA dentro del Convenio Marco con la Universidad de La Coru a y la Fundaci n de
19. es absolutos en concentraciones al cabo de 10000 d as normalizados respecto a la concentraci n Co del agua de entrada Se ha tomado como soluci n exacta la obtenida con el Modelo 1D con una malla de 19705 nudos El mapa de errores para el Modelo 2D y una malla con 2 112 nudos se representa en la Figura 26 Se observa que el error m ximo es pr ximo al 50 y que las isol neas de errores presentan un aspecto parecido a las isol neas de concentraciones comp rense las Figuras 14 y 26 Se puede comprobar que los errores son grandes en todo el dominio aunque los m ximos relativos se alcanzan en torno a las fracturas y especialmente cerca del punto de inyecci n 76 R Juanes y J Samper O 1 2 3 4 5 6 7 8 9 10 Figura 26 Mapa de errores de concentraciones a tiempo T 10 000 d as para el Modelo 2D con 2112 nudos M2D_2K 0 1 2 3 4 5 6 7 8 9 10 Figura 27 Mapa de errores de concentraciones a tiempo T 10 000 d as para el Modelo 2D con 6105 nudos M2D_6K 4 i 0 1 2 3 4 5 6 7 8 9 10 Figura 28 Mapa de errores de concentraciones a tiempo T 10 000 d as para el Modelo 2D con 18149 nudos M2D_20K A A e a r y DIO o OQ o 00 00 0000 o00000 oc000o 00 o 0000 e 2 3 4 5 6 7 8 99 0 Figura 29 Mapa de errores de concentraciones a tiempo 1 1U UUU dias para el Modelo 1D con 2201 nudos M1D_2K Una formulaci n general y eficiente de las fracturas en el MEF I
20. luto sigue el patr n que cabr a esperar El soluto se transporta desde el pozo de inyecci n a trav s de la matriz rocosa por advecci n y dispersi n hidrodin mica hasta que alcanza la fractura A partir de ese momento el soluto se canaliza por la fractura debido a que la conductividad hidr ulica de este material es mucho mayor Se debe destacar la gran calidad de la soluci n num rica en el entorno de la fractura donde se reproduce una variaci n espacial de la concentraci n en todo su rango de 0 a C de forma muy localizada ver Figura 32 A medida que avanza el tiempo en la parte anterior del macizo y lt 0 y cerca del contorno derecho x 80m el soluto migra hacia el punto de bombeo mientras que cerca del contorno izquierdo x 0 el proceso dominante es la difusi n en la matriz 80 R Juanes y J Samper Figura 32 Mapa de concentraciones al cabo de 200 d as Figura 33 Mapa de concentraciones al cabo de 1000 d as Figura 34 Mapa de concentraciones al cabo de 5000 d as Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 81 CONCLUSIONES En el segundo de esta serie de dos art culos se han presentado las caracter sticas del c digo TRANMEF 3 as como su aplicaci n a varios casos sint ticos de flujo y transporte de solutos en medios porosos y fracturados En este c digo se ha implementado una formulaci n matem tica compacta general y efici
21. ntaria de error Soluci n num rica Se han considerado dos modelos distintos para la discretizaci n en elementos finitos de la fractura El Modelo 2D emplea la formulaci n tradicional de elementos finitos y discretiza la fractura con elementos 2D El Modelo 1D utiliza la formulaci n propuesta en esta serie de art culos discretizando la fractura con elementos 1D Los valores num ricos de los par metros empleados en la soluci n anal tica son los siguientes b 0 04m vs 0 01md pp 0 1 m 0 01 Dm 107 m d y Co lmgL Para reproducir las condiciones de la soluci n anal tica en la soluci n num rica los contornos inferior y 0 04m y superior y 4m se consideran impermeables condici n de contorno de tipo Neumann y en los contornos izquierdo x 0m y derecho x 10m se impone el nivel piezom trico h 1m y h 0 m respectivamente condici n de contorno de tipo Dirichlet La conductividad hidr ulica de la fractura K 107 md es cuatro rdenes de magnitud mayor que la de la matriz Km 107 md De este modo se asegura que s lo existe flujo de agua en la direcci n x y que la velocidad del agua en la matriz sea muy peque a Para la matriz se ha empleado una dispersividad longitudinal a Om y una dispersividad transversal r 100m para reproducir un coeficiente de difusi n efectivo en la matriz Dm AT qs ar Km i 107 m d 7 La soluci n anal tica dada por Gris
22. o var a en la direcci n zx El transporte del soluto en la fractura es debido nicamente a la advecci n Existe equilibrio en la interfaz matriz fractura S lo existe difusi n en la matriz en la direcci n y La distancia a la que penetra el soluto es peque a en comparaci n con el espesor de la matriz O NR E LLE AU AS Para estas condiciones las ecuaciones que gobiernan las concentraciones en la fractura cf y en la matriz Cm son respectivamente 0 Cm brb ERE or Dn ee 2 gt 0 4 0 4 gt 0 1 y Cm cr mar Pm gt E 2 a D ay gt 0 y gt 0 t gt 0 2 Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 65 con las condiciones iniciales te 0 c lg 0 x gt 0 y gt 0 t 0 3 y las condiciones de contorno cs 0 t cm 0 0 t Co x 0 y 0 t gt 0 4 cplx t Cale 0 t D y 0 t gt 0 5 donde b es la semiapertura de la fractura fs es la porosidad de la fractura m es la porosidad de la matriz vs es la velocidad del agua en la fractura Dm es el coeficiente de difusi n efectivo en la matriz y x y denotan las coordenadas globales del dominio El origen de coordenadas se encuentra en la esquina inferior izquierda de la matriz La soluci n anal tica a este problema viene dada adaptada de por es Co 3 Din 06 b x E r vr 6 0 t lt x vs Cin erfc Oi EE Lao C 21Dm m t 2 0 14 0 t lt x vs donde erfc es la funci n compleme
23. ontornos laterales se impone el nivel piezom trico h 1m en el contorno izquierdo y h Om en el contorno derecho de tal forma que el flujo discurre de izquierda a derecha preferentemente por la red de fracturas debido a su mayor conductividad hidr ulica El flujo de agua es estacionario Inicialmente la concentraci n del soluto es cero y el agua que entra por la fractura lo hace con una concentraci n Co Este caso considera aspectos m s realistas del flujo y transporte en medios fracturados Por ello sirve para ilustrar la aplicabilidad de la formulaci n a condiciones reales Adem s se ha utilizado este caso para realizar una comparaci n exhaustiva entre el modelo num rio propuesto y el tratamiento convencional de elementos finitos Soluci n num rica Empleando la misma nomenclatura que en el caso anterior el Modelo 2D corresponde a la discretizaci n de las fracturas mediante elementos bidimensionales mientras que el Modelo 1D se basa en la formulaci n propuesta en el primer art culo de esta serie en el que se discretizan las fracturas con elementos unidimensionales En la Tabla I se recogen los valores num ricos de los par metros de flujo y transporte de la matriz y de la zona fracturada La concentraci n del agua que entra por la fractura es Co 1mgL Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 69 4m C Red de fracturas
24. plica que el programa estrictamente no es aplicable a problemas con superficie libre simulaciones en la zona no saturada o con efectos de densidad variable En la ecuaci n se consideran el flujo de tipo Darcy el almacenamiento recarga volum trica y t rminos fuente sumidero Las condiciones de contorno son generales siempre que sean lineales Neumann Dirichlet y Mixta distribuidas a lo largo de los contornos exteriores o interiores 2 Transporte de solutos Despu s de resolver la ecuaci n de flujo bien en r gimen estacionario o transitorio y una vez que se han calculado las velocidades de Darcy en los puntos de Gauss de los elementos del dominio TRANMEF 3 permite la simulaci n del transporte de solutos incluyendo los siguientes t rminos difusi n molecular dispersi n hidrodin mica advecci n almacenamiento recarga fuentes sumideros desintegraci n radiactiva y adsorci n La difusi n en la matriz a lo largo de las fracturas se considera autom ticamente dado que el programa resuelve las ecuaciones de flujo y transporte en la formaci n porosa as como en las fracturas Se pueden simular simult neamente tantas especies como se desee cada una de ellas con su propia constante de desintegraci n A y su coeficiente de distribuci n para adsorci n Ka Se incluyen los mismos tres tipos de condiciones de contorno Neuman Dirichlet y Mixta 3 Transporte de calor Son aplicables las mismas consideraciones que para el transpo
25. racture These simulations not only illustrate the potential and applicability of the proposed formulation but also demonstrate its numerical advantages over traditional formulations Representing fractures with elements of lower dimension than the domain e g 1D elements in a 2D medium leads to more accurate and numerically efficient solutions OUniversitat Polit cnica de Catalunya Espa a ISSN 0213 1315 Recibido Julio 1999 62 R Juanes y J Samper INTRODUCCI N En esta serie de dos art culos se presenta una formulaci n general y eficiente para el tratamiento num rico de fracturas y condiciones de contorno en modelos de elementos finitos En este segundo art culo de la serie se presentan una serie de ejemplos sint ticos de flujo y transporte de solutos en medios porosos y fracturados que ilustran el potencial y la eficiencia de la formulaci n propuesta La incorporaci n de fracturas y v as preferentes de flujo es crucial a la hora de representar de forma adecuada el comportamiento hidromec nico de las formaciones geol gicas La existencia de zonas fracturadas cuya conductividad hidr ulica es varios rdenes de magnitud mayor que la del macizo rocoso suele dominar el esquema de flujo regional Dado que el campo de velocidades del agua est controlado en gran medida por la presencia de fracturas la correcta predicci n de los fen menos de transporte de solutos depende de una representaci n adecuada de la fracturaci n del
26. rte de solutos ya que existe una analog a perfecta entre ambas ecuaciones ver la Tabla I en el primer art culo de la serie El c digo incorpora una rutina para la renumeraci n de los nudos de la malla de elementos finitos basada en la teor a de grafos y que es completamente transparente al usuario El objeto de la renumeraci n nodal es reducir el perfil de las matrices de coeficientes que resultan de la discretizaci n de las ecuaciones diferenciales Se debe enfatizar la importancia de dicha renumeraci n ya que los requisitos de almacenamiento en memoria y tiempo de CPU disminuyen varios rdenes de magnitud para simulaciones con un n mero de nudos elevado El programa principal y todas sus subrutinas constituyen un c digo de m s de 25000 l neas implementado de forma eficiente en lenguaje FORTRAN y verificado en numerosos casos con soluci n anal tica En resumen TRANMEF 3 es un programa que constituye un modelo num rico general para la resoluci n de problemas de flujo transporte de solutos y transporte de calor en formaciones subterr neas heterog neas medios porosos y fracturados 64 R Juanes y J Samper CASO 1 BLOQUE RECTANGULAR DE ROCA CON UNA FRACTURA Descripci n del problema El primer caso corresponde a un problema cl sico de transporte de un soluto conservativo a trav s de un medio fracturado El objeto de este caso es doble Por una parte permite la verificaci n de la formulaci n adoptada ya que dispon
27. tras que los obtenidos con una malla de 6501 nudos M1D_6K se representa en las Figuras 21 y 22 Se puede comprobar que ambas simulaciones proporcionan resultados cuyas diferencias son gr ficamente casi imperceptibles Los valores num ricos s lo difieren en la segunda o tercera cifra significativa Los resultados con una malla de 19705 nudos M1D_20K que se muestran en las Figuras 23 y 24 son virtualmente id nticos a los de la pasada M1D_6K Las curvas de paso del soluto concentraci n en funci n del tiempo en el punto de observaci n ver Figura 7 para las seis simulaciones realizadas se muestran en la Figura 25 El punto de observaci n est situado en el interior de una fractura en el centro geom trico del dominio x 5m y 2m Las curvas correspondientes al Modelo 1D pasadas M1D_2K M1D_6K y M1D_20K son coincidentes Asimismo puede observarse que las soluciones obtenidas con el Modelo 2D se aproximan a las del Modelo 1D a medida que se utiliza una discretizaci n m s fina 0 8 o u Concentracion mg L o o w gt 0 2 Modelo 1D M2D_20K 0 1 M2D_6K M2D_2K o 1000 2000 3000 4000 5000 6000 7000 8000 9000 10000 Tiempo dias Figura 25 Curvas de paso del soluto en el punto de observaci n indicado en la Figura 7 para los Modelos 1D y 2D y las distintas mallas empleadas Con el fin de cuantificar el mejor comportamiento num rico del Modelo 1D frente al Modelo 2D se presentan los mapas de error
28. ue la soluci n anal tica no considera Se observa que las isol neas de concentraciones calculadas con ambos modelos coinciden con las de la soluci n anal tica comprob ndose por tanto que el c digo resuelve correctamente la ecuaci n de transporte de solutos en estas condiciones TIEMPO 1000 dias CAS 2 3 4 5 6 7 8 9 10 X m Figura 3 Concentraciones en la matriz al cabo de 1000 d as Soluci n anal tica Una formulaci n general y eficiente de las fracturas en el MEF II Aplicaci n a casos sint ticos 67 TIEMPO 1000 dias X m Figura 4 Concentraciones en la matriz al cabo de 1000 d as Soluci n num rica del Modelo 2D TIEMPO 1000 dias 5 10 X m Figura 5 Concentraciones en la matriz al cabo de 1000 d as Soluci n num rica del Modelo 1D La Figura 6 muestra la distribuci n de la concentraci n del soluto a lo largo de la fractura y 0 al cabo de 200 600 y 1000 d as Se representan conjuntamente la soluci n anal tica Ecuaci n 6 y la soluci n num rica para ambos modelos Cabe destacar que la soluci n del Modelo 1D fracturas con elementos 1D es m s precisa que la soluci n del Modelo 2D fracturas con elementos 2D Para una misma discretizaci n espacial y temporal el Modelo 1D presenta menos dispersi n num rica Las peque as discrepancias entre la soluci n del Modelo 1D y la soluci n anal tica se deben a que el modelo num rico considera una peque
29. ultados se concluye que el Modelo 1D en el que las fracturas se discretizan mediante elementos unidimensionales propuesto en esta serie de art culos proporciona resultados mucho mejores que los del Modelo 2D fracturas discretizadas con elementos bidimensionales para la soluci n num rica de la ecuaci n del transporte de solutos en medios fracturados Identificador N mero Memoria Tiempo CPU Error m x Error medio de la malla de nudos MB s M2D_2K 2112 4 2 46 9 49 7 20 4 M2D_6K 6105 19 0 374 8 24 2 9 5 M2D_20K 18149 159 6 3438 6 8 0 2 3 M1D_2K 2201 4 7 50 5 2 48 0 16 M1D_6K 6501 33 0 450 3 0 21 0 00 M1D_20K 19705 196 4 4320 2 Tabla III Caracter sticas y resultados de las simulaciones realizadas para el Caso 2 CASO 3 BLOQUE DE ROCA TRIDIMENSIONAL CON UNA FRACTURA VERTICAL Descripci n del problema El tercer caso sint tico corresponde a un ensayo de trazadores en un bloque de roca tridimensional de dimensiones 80 x 40 x 20m que contiene una fractura vertical en el plano y 0m Figura 30 Todos los contornos son impermeables condici n de contorno de tipo Neumann con caudal nulo excepto la arista vertical x 0m y 20m en la que se fija el nivel piezom trico h 0m Esta arista coincide con un pozo de inyecci n por el que se inyecta agua con una concentraci n Co 1mgL Se bombea un caudal puntual Q 1m d en el v rtice x 80m
30. umerical modeling of the impact of the tunnel construction on the groundwater system at Asp Asp HRL Task Force 45 Progress Report 1999 G E Grisak y J F Pickens An analytical solution for solute transport through fractured media with matrix diffusion J Hydrol Vol 52 pp 47 57 1981 J Carrera X S nchez Vila I Benet A Medina G Galarza y J Guimer On matrix diffusion formulations solution methods and qualitative effects Hydrogeology Journal Vol 6 pp 178 190 1998 H S Carslaw y J C Jaeger Conduction of Heat in Solids Oxford University Press London 2nd ed pp 510 1959 R Juanes J Samper y J Molinero Modelizaci n num rica tridimensional de flujo de agua y transporte de solutos y calor en medios porosos y fracturados I Congreso Ib rico de Geoqu mica Soria Espa a pp 111 118 1997
31. y 20m z 0m Se considera flujo estacionario y transporte transitorio Inicialmente la concentraci n de soluto es nula en todo el dominio El error medio se ha calculado como medio 4 ie e z dz utilizando la regla de Simpson sobre una superficie krigeada con 101 x 41 puntos 78 R Juanes y J Samper Con este caso se pretende ilustrar el potencial del c digo TRANMEF 3 para la simulaci n de problemas complejos de flujo y transporte de solutos en medios porosos y fracturados tridimensionales Bombeo Inyecci n Figura 30 Esquema del Caso 3 correspondiente a una inyecci n de soluto por una de las aristas y bombeo en la esquina opuesta Soluci n num rica Al igual que en el Caso 2 se han utilizado dos esquemas num ricos distintos el Modelo 3D en el que las fracturas se discretizan mediante elementos tridimensionales y el Modelo 2D que emplea la formulaci n propuesta en esta serie de art culos y en el que las fracturas se discretizan con elementos bidimensionales Dado que en el caso sint tico anterior ya se ha realizado un an lisis detallado de las ventajas num ricas de la formulaci n propuesta se muestran s lamente los resultados correspondientes al Modelo 2D Para la discretizaci n del pozo de inyecci n se han empleado elementos unidimensionales Los par metros de flujo y transporte de solutos utilizados para la matriz rocosa la fractura y el pozo se recogen en la Tabla IV P
Download Pdf Manuals
Related Search