Análisis mecánico de albañilería arqueológica de adobe bajo cargas de compresión uniaxial : El caso de Huaca de la Luna

Análisis mecánico de albañilería arqueológica de adobe bajo cargas de compresión uniaxial: El caso de Huaca de la Luna Eduardo Ramírez1, Mijail Montesinos1, Rui Marques2, Ricardo Morales3, Santiago Uceda3, Paulo B. Lourenço2, Rafael Aguilar1 1Department of Engineering, Civil Engineering Division. Pontificia Universidad Católica del Perú-PUCP (PERU) {jeramirezc; mijail.montesinos; raguilar}@pucp.pe 2Department of Civil Engineering. ISISE, University of Minho (PORTUGAL) marques[email protected], [email protected] 3Facultad de Ciencias Sociales, Universidad Nacional de Trujillo (PERU) [email protected], [email protected] ABSTRACT This article presents an experimental and numerical study for the mechanical characterization under uniaxial compressive loading of the adobe masonry of one of the most emblematic archaeological complex in Peru, ‘Huaca de la Luna’ (100-650AD). Compression tests of prisms were carried out with original material brought to the laboratory. For measuring local deformations in the tests, displacement transducers were used which were complemented by a digital image correlation system which allowed a better understanding of the failure mechanism. The tests were then numerically simulated by modelling the masonry as a continuum media. Several approaches were considered concerning the geometrical modelling, namely 2D and 3D simplified models, and 3D refined models based on a photogrammetric reconstruction. The results showed a good approximation between the numerical prediction and the experimental response in all cases. However, the 3D models with irregular geometries seem to reproduce better the cracking pattern observed in the tests. Keywords: archaeological constructions, adobe masonry, mechanical testing, digital image correlation, numerical modelling 1. INTRODUCCIÓN El adobe es una forma tradicional de mampostería no reforzada que ha sido usada desde tiempos ancestrales [1]. Sin embargo, las construcciones de adobe presentan elevada vulnerabilidad sísmica debido a su baja resistencia a tracción, a la elevada fuerza inercial movilizada y a su comportamiento frágil [2]. En efecto, la observación post-sismo muestra que las construcciones de adobe sufren repetitivamente grandes daños, ej. [3]. Teniendo en cuenta que los ensayos experimentales de mampostería son costosos, una buena forma de analizar su comportamiento es a través del modelamiento numérico [4-6]. Sin embargo, una cantidad mínima de experimentos debe ser realizada para permitir la calibración y validación de los modelos. Por otro lado, la mampostería de adobe se 77 Mechanical analysis of adobe archaeological masonry under uniaxial compressive loading Third International Conference on Mechanical Models in Structural Engineering University of Seville. 24-26 june 2015. caracteriza por presentar una elevada variabilidad en sus propiedades mecánicas, aspecto que tiene que ser considerado tanto en los experimentos como en la simulación numérica. En el caso de construcciones arqueológicas la variabilidad geométrica se incrementa, la cual puede ser original, o asociada a procesos de degradación del material por su exposición en el tiempo. Las construcciones masivas de adobe en Perú son un ejemplo evidente de esta variabilidad. El presente artículo estudia la caracterización del material que constituye el sistema estructural de uno de los complejos arqueológicos de tierra más emblemáticos del Perú: “La Huaca de la Luna”. La Huaca se construyó entre los años 100 y 650 D.C. y se encuentra localizada en la costa norte del país, a 8 kilómetros del centro de la ciudad de Trujillo. Este sitio arqueológico fue uno de los templos más importantes de la cultura pre-inca ‘Mochica’ (I-IX D.C.) [7] y es particularmente importante en la actualidad por su fachada decorada muy bien preservada y su sistema de construcción con templos superpuestos en sucesivos periodos. Las construcciones patrimoniales de adobe en Perú, y en particular la Huaca de la Luna por su ubicación en la costa norte, están altamente expuestas a la acción de los sismos así como a fenómenos meteorológicos, ej. el Niño. Este estudio es parte de un análisis integral de la seguridad símica del monumento y contribuye con la caracterización del sistema de albañilería de adobe arqueológico. Asimismo, en este estudio se proponen propiedades mecánicas del sistema de albañilería que serán de utilidad para los futuros trabajos de análisis numérico predictivo en este monumento. 2. ENSAYOS EXPERIMENTALES La campaña experimental se desarrolló en el laboratorio de estructuras de la Pontificia Universidad Católica del Perú - PUCP. La caracterización mecánica se realizó utilizando componentes originales de la mampostería de adobe en la Huaca de la Luna. Estos materiales fueron recolectados exclusivamente para el presente estudio durante un trabajo especial de excavación arqueológica, conforme se muestra en la Figura 1. (a) (b) (c) Figura 1. Vista general de Huaca de la Luna y sus paredes: (a) Fachada Norte, (b) bloques de adobe usados para la caracterización mecánica y (c) proceso de extracción del mortero. Para la caracterización de la mampostería se construyeron 3 prismas de 210x320x770 mm3 siguiendo las recomendaciones establecidas en el RNE – E080 [8] y en HB 195-2002 [9]. En cada prisma se emplearon 5 unidades de adobe que se unieron con mortero de 2.5 cm de espesor 78 E. Ramírez1, M. Montesinos1, R. Marques2, R. Morales3, S. Uceda3, P.B. Lourenço2, R. Aguilar1 aproximadamente. Debido a que los adobes originales de la Huaca presentan una gran variabilidad dimensional, los prismas que se construyeron en laboratorio presentaron una geometría irregular (ver Figura 2). (a) (b) (c) Figura 2. Preparación de los prismas para la caracterización mecánica: (a) colocación del mortero, (b) construcción de prisma de adobe y (c) aspecto de los prismas antes del ensayo. El ensayo de compresión uniaxial se realizó por control de desplazamientos a una velocidad de 0.25 mm/min de forma tal que se alcance la falla en aproximadamente 30 minutos, conforme se sugiere en EN 1052-1[10]. Para evitar concentraciones de esfuerzos en la zona de contacto del equipo con el prisma se utilizó una capa de arena y una plancha de madera como material de transición [9]. Se midieron las deformaciones locales utilizando LVDTs fijados en las dos caras laterales de los prismas, como se puede observar en la Figura 3a-b. Adicionalmente, para el ensayo del tercer prisma se empleó un equipo de medición óptica de deformaciones a través de correlación de imagen digital (ARAMIS) con una precisión de hasta 0.005% [11]. Para este efecto, se analizó con este equipo una superficie de 200x160 mm2 pintada en una de las caras frontales (Figura 3c). (a) (b) (c) Figura 3. Instrumentación de prismas de adobe para el ensayo de compresión: (a) LVDTs en la cara frontal, (b) LVDTs en la cara posterior y (c) área con pintura en escala de grises para el análisis con el ARAMIS. En los ensayos se observa la aparición de grietas paralelas a la dirección de la carga. Se observó que las grietas aparecen primero en los extremos que están en contacto con los cabezales del equipo, lo cual es probablemente debido a que en esta zona se haya producido una concentración de esfuerzos, a pesar de la colocación de la capa de arena y plancha de madera. En la Figura 4 se observa la Prisma #01 Prisma #02 Prisma #03 79 Mechanical analysis of adobe archaeological masonry under uniaxial compressive loading Third International Conference on Mechanical Models in Structural Engineering University of Seville. 24-26 june 2015. evolución de las grietas durante el ensayo de compresión en el prisma #3 (los números indican la secuencia de aparición). Las grietas se originan en la zona de contacto con el equipo (grietas 1, 2, 3 y 5). Una grieta en la zona central (grieta 4) empieza en la junta de mortero y luego se propaga a las unidades conforme avanza el ensayo. La grieta 6 se genera en una cara lateral del prisma y se propaga verticalmente. Luego, se siguen generando grietas en la junta de mortero (grieta 7). Debido al progreso de la grieta 4, se produce la grieta 8 que presenta una longitud considerable. No se logró registrar la aparición y evolución del resto de grietas debido a que cuando se alcanzó la carga máxima estas se propagaron súbitamente en todo el espécimen. (a) (b) (c) (d) Figura 4. Evolución de grietas en el prisma 03: (a) vista frontal, (b) vista posterior y (c-d) caras laterales. El análisis de los resultados del ensayo con control óptico de deformaciones confirma lo observado y muestra que las grietas se originan en el mortero cuando el sistema aún presenta un comportamiento elástico y se encuentra aproximadamente al 20% de la carga máxima (Figura 4b). Por otro lado, la aparición de grietas en los adobes coincidió con el inicio de la etapa no lineal en el prisma, al 95% de la carga máxima aproximadamente (Figura 4c). Durante el transcurso de la etapa no lineal, las grietas se propagan tanto en el mortero como en el adobe (Figura 4d). Se observó además que existe una interacción entre mortero y adobe durante el ensayo que se manifiesta en que el mortero sufre grandes deformaciones longitudinales y transversales (similar a lo reportado en [12]). Las deformaciones transversales (perpendiculares a la carga) producen tracción en las unidades de adobe, dando origen a las grietas verticales. El proceso descrito de formación de grietas se puede observar en la Figura 5. 2 7 4 8 1 3 5 4 6 80 E. Ramírez1, M. Montesinos1, R. Marques2, R. Morales3, S. Uceda3, P.B. Lourenço2, R. Aguilar1 (a) (b) (c) (d) Figura 5. Formación de grietas: (a) estado inicial, (b) primeras grietas en la junta de mortero,, (c) primera grieta en ladrillo de adobe y (d) propagación de las grietas durante la etapa no-lineal. (El color azul indica mínima deformación en el espécimen mientras que el rojo representa una grieta completamente abierta) A partir del registro de la fuerza aplicada y de las deformaciones medidas mediante los LVDTs se obtuvieron las curvas esfuerzo deformación, que junto a la forma de falla de los prismas permiten evaluar el comportamiento a compresión uniaxial de la albañilería (Figura 6). El esfuerzo máximo a compresión promedio (σmáx) fue de 0.36 MPa con un coeficiente de variación (CV) de 15%, y la deformación unitaria promedio para el esfuerzo máximo fue de 0.0035 mm/mm con un CV de 12%. Utilizando el criterio propuesto por [13] se obtuvo un módulo de elasticidad promedio (E) de 107 MPa con un CV de 28%. La relación entre el módulo de elasticidad y la resistencia a compresión E/σmáx de 297 se encuentra dentro del rango propuesto por Tomaževič [14], que varía entre 200 y 1000. No fue posible obtener la respuesta completa post-pico, debido a que cuando las grietas atravesaban la superficie de apoyo de los LVDTs lo que ocasionó su desprendimiento de los prismas. Los parámetros mecánicos obtenidos en los ensayos se muestran en la Tabla 1. Tabla 1. Propiedades mecánicas de los prismas de adobe Prisma σ max (MPa) e max (x10-3) E (MPa) 01 0.31 3.85 84 02 0.35 3.03 102 03 0.41 3.47 136 Promedio 0.36 3.45 107 σ max: máximo esfuerzo a compresión e max: deformación unitaria local para el máximo esfuerzo a compresión E: Módulo de elasticidad 5.0 3.0 2.0 1.0 0.0 4.0 [%] Mortar joint Adobe unit Adobe unit 5.0 3.0 2.0 1.0 0.0 4.0 [%] Cracks on mortar 5.0 3.0 2.0 1.0 0.0 4.0 [%] First crack on adobe 5.0 3.0 2.0 1.0 0.0 4.0 [%] 81 Mechanical analysis of adobe archaeological masonry under uniaxial compressive loading Third International Conference on Mechanical Models in Structural Engineering University of Seville. 24-26 june 2015. (a) (b) (c) (d) Figura 6. Curvas esfuerzo-deformación del (a) prisma 01, (b)prisma 02, (c) prisma 03 y (d) envolvente, El análisis de los componentes de la albañilería a partir de las imágenes captadas por el equipo de video-correlación muestra que el mortero alcanza deformaciones hasta 12 veces más que el adobe, como se observa en la Figura 7b, lo que podría indicar que el comportamiento global del sistema depende fuertemente de la calidad del mortero. La influencia del mortero en el comportamiento resulta notoria al comparar la deformación medida en la albañilería con el LVDT y con el sistema de correlación de imagen digital, una vez que el LVDT cubre una menor extensión de mortero. La longitud inicial del LVDT era de 320 mm, a lo largo de la cual cerca del 80% era adobe y el resto mortero. En el caso de la longitud analizada con el ARAMIS (100 mm) el 70% aproximadamente era adobe. A pesar de la buena correspondencia de resultados, dicha diferencia podría ser el parámetro que genere que la curva de esfuerzo-deformación obtenida del ARAMIS tenga una menor pendiente y aparentemente mayores deformaciones, como se observa en la Figura 7. Aunque el comportamiento post-pico sea semejante, el mortero presenta mayores ductilidad y tasa de deformación plástica que las unidades de adobe. La ductilidad calculada en el mortero es igual a 2.08 mm, mientras que en el adobe es sólo de 0.08 mm, por otro lado, en la albañilería se estimó una ductilidad de 1.45 mm. El coeficiente de Poisson se calculó también en base a las mediciones realizadas con el equipo de correlación de imagen digital. En este caso se consideró estimar la componente horizontal de deformación a lo larga de una línea diagonal de forma tal que se consiga medir la deformación de la albañilería. 0.00 0.09 0.18 0.27 0.36 0.45 0.000 0.002 0.004 0.006 Stress (MPa) Strain (mm/mm) LVDT 01 LVDT 02 LVDT 03 Mean LVDTs 0.00 0.09 0.18 0.27 0.36 0.45 0.000 0.002 0.004 0.006 Stress (MPa) Strain (mm/mm) LVDT 01 LVDT 02 Mean LVDTs 0.00 0.09 0.18 0.27 0.36 0.45 0.000 0.002 0.004 0.006 Stress (MPa) Strain (mm/mm) LVDT 02 LVDT 03 Mean LVDTs 0.00 0.09 0.18 0.27 0.36 0.45 0.000 0.002 0.004 0.006 Stress (MPa) Strain (mm/mm) Experimental Prism #1 Prism #2 Prism #3 envelope 82 E. Ramírez1, M. Montesinos1, R. Marques2, R. Morales3, S. Uceda3, P.B. Lourenço2, R. Aguilar1 (a) (b) Figura 7. Comportamiento de la albañilería y sus componentes bajo cargas de compresión: (a) desfase entre el sistema de correlación óptico y los LVDTs y (b) componentes de la albañilería durante el ensayo de compresión. La influencia del mortero en el comportamiento resulta evidente al comparar los resultados experimentales obtenidos a partir de los LVDTs y del sistema ARAMIS, como se mencionó anteriormente. Además, es importante recordar que los LVDTs no lograron registrar el comportamiento post-pico completamente, lo que por un lado produce que el desfase mostrado en la Figura 7a sea aún más evidente, pero por otro lado indica también la utilidad del sistema de medición de deformaciones sin contacto. 3. MODELACIÓN NUMÉRICA La modelación numérica se desarrolló considerando modelos continuos de elementos finitos y asumiendo un material homogéneo. El comportamiento del material se representó usando el modelo total strain crack model disponible en DIANA [15]. El comportamiento a compresión fue descrito por una tendencia parabólica, mientras que en tracción se consideró una ley de decrecimiento exponencial. Parte de las propiedades empleadas en el modelo se obtuvieron de los ensayos experimentales. Sin embargo, debido a que no se realizaron ensayos de tracción, fue necesario recurrir a la literatura cuando no se disponía de información experimental. Se asumió la resistencia a tracción ft como fc/10, según lo propuesto en [16, 17]. Para el cálculo de la energía de fractura en compresión se consideró un factor de ductilidad (relación entre energía de fractura y resistencia a la compresión) igual a 1.6 mm [18]: fc = 0.576 N/mm. Para la energía de fractura en tracción se consideró igual a 0.029 mm de acuerdo a recomendaciones mostradas en [19]. Adicionalmente, se realizó un análisis paramétrico de las energías de fractura en compresión y tracción, para verificar que los valores empleados sean los más apropiados. Los valores de energía de fractura a compresión y tensión a emplear se seleccionaron considerando el criterio de obtener una respuesta más estable y con mayor similitud al comportamiento observado en los ensayos de compresión de prismas. Las propiedades utilizadas para la mampostería se resumen en la Tabla 2. La densidad utilizada fue la reportada en un informe técnico del proyecto arqueológico Huaca de la Luna [20]. Para el ajuste de la respuesta a compresión en función de la energía de fractura a tracción (Figura 8a), el mayor valor se tomó de Tarque [4]. En el resto de casos se tomaron los valores límites 0 0.09 0.18 0.27 0.36 0.45 0 0.005 0.01 0.015 0.02 Stress (MPa) Strain (mm/mm) Masonry (ARAMIS) Experimental envelope 0 0.09 0.18 0.27 0.36 0.45 0 0.01 0.02 0.03 0.04 0.05 Stress (MPa) Strain (mm/mm) Masonry (ARAMIS) Adobe (ARAMIS) Mortar (ARAMIS) 83 Mechanical analysis of adobe archaeological masonry under uniaxial compressive loading Third International Conference on Mechanical Models in Structural Engineering University of Seville. 24-26 june 2015. del rango propuesto por Angelillo [19] y su valor recomendado con base en la ductilidad de las unidades a tracción. En cuanto a la energía de fractura a compresión, el mayor valor corresponde a los datos experimentales en Almeida [21]. Los valores intermedios se encuentran al multiplicar la resistencia a compresión por un factor de ductilidad igual a 1.6 mm [18] y 1.0 mm [17]. Finalmente, el menor de los valores se determinó haciendo uso de la regresión lineal propuesta en [21], para la determinación de la energía de fractura a compresión con base en la resistencia a compresión. (a) (b) Figura 8. Comportamiento a compresión del prisma #2 con análisis paramétrico de la (a) energía de fractura de tensión y (b) energía de fractura a compresión Tabla 2. Propiedades del adobe usadas para los modelos numéricos Basic properties Density (kg/m3) 1735 [20] Elasticity modulus (MPa) 115 Poisson ratio 0.2 Inelastic properties Compression Compressive strength fc (MPa) 0.36 Compressive fracture energy Gc (N/m) 576 Tension Tensile strength ft (MPa) 0.036 Tensile fracture energy Gt (N/m) 1.044 En el análisis se consideraron cuatro tipos de modelos: dos usando elementos tipo lámina y dos utilizando elementos sólidos. En el primer caso se utilizaron elementos tipo shell de 8 nodos en estado plano de esfuerzo, CQ16M [15]. Con referencia al prisma #1 en la Figura 7a, se construyó un primer modelo con elementos cuadriláteros regulares. Después, a partir de un modelo 3D del prisma obtenido mediante reconstrucción fotogramétrica (Figura 9b), se tomó una sección representativa del prisma para desarrollar un modelo con elementos lámina, ver Figura 7c. En una fase posterior se emplearon elementos sólidos con 20 nodos, CHX60 [15]. También en este caso, se modeló un prisma regular, y otro considerando el modelo obtenido de la reconstrucción fotogramétrica (Figura 9d). La Figura 10 muestra todos los modelos geométricos considerados. En ninguno de los casos se consideró 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.000 0.001 0.002 0.003 0.004 0.005 Stress (MPa) Strain (mm/mm) Gt = 10 N/m [4] Gt = 1.44 N/m [19] Gt = 1.044 N/m [19] Gt = 0.648 N/m [19] [4] [19] [19] [19] 0.00 0.05 0.10 0.15 0.20 0.25 0.30 0.35 0.000 0.001 0.002 0.003 0.004 0.005 Stress (MPa) Strain (mm/mm) Gc = 1070 N/m [20] Gc = 576 N/m [18] Gc = 360 N/m [17] Gc = 120 N/m [20] [21] [18] [17] [21] 84 E. Ramírez1, M. Montesinos1, R. Marques2, R. Morales3, S. Uceda3, P.B. Lourenço2, R. Aguilar1 modelar la interacción entre el prisma y el equipo de compresión, se restringieron los desplazamientos en la base del prisma. (a) (b) (c) (d) Figura 9. Generación de mallas: (a) fotografía del prisma #1, (b) modelo 3D obtenido de la reconstrucción fotogramétrica, (c) modelo de elementos finitos 2D y (d) modelo 3D usando elementos sólidos. El ensayo de compresión fue simulado aplicando un desplazamiento en la extremidad superior del modelo. El método de solución utilizado fue el Newton-Raphson Modificado en conjugación con la técnica de line-search, de modo a obtener un compromiso entre la carga computacional y la precisión de los resultados. Los resultados de la simulación numérica de los ensayos de compresión se muestran en la Figura 11. Se observa, en general, una buena aproximación a la respuesta general de los varios prismas ensayados. Es evidente que al considerar una geometría más detallada con elementos sólidos se consiguen resultados más cercanos a los experimentales, en especial con respecto a la predicción del modo de falla. Por otro lado, el costo computacional al emplear elementos sólidos es mucho más alto comparativamente a los modelos 2D y la respuesta obtenida de ambos no difiere notablemente. Cuando se simplifica la geometría, la predicción sigue aproximando mínimamente la respuesta experimental (Figura 11 d). (a) (b) Figura 10. Modelos con elementos finitos: (a) 2D en estado plano de esfuerzo y (b) sólidos. (Prisma #1 al #3, y geometría simplificada, de izquierda a derecha en ambos casos) 85