sábado, 6 de julio de 2013

Herramienta numérica de análisis para losas de hormigón armado sometidas a aceleraciones verticales sísmicas

Debido a la magnitud momento de 8.8 del terremoto que sacudió la zona centro sur de Chile el 27 de febrero del año 2010, muchas estructuras de hormigón armado sufrieron deformaciones que sobrepasaron el límite elástico de varios de sus elementos estructurales. Deformaciones irreversibles se evidenciaron por medio de la aparición de fisuras, grietas e incluso desprendimiento del hormigón y falla del acero de refuerzo. Esto condujo a la falla de importantes elementos estructurales e incluso el colapso de algunas estructuras (Betanzo 2010). En general las fallas más importantes suceden en muros, vigas y columnas de hormigón armado, en cambio en las losas los daños son menores. Tal vez debido a esto las losas de hormigón armado son menos estudiadas. Sin embargo, las fotos de la Figura 1 muestran ejemplos de daños en losas ocurridas en el terremoto del 27/02/2010. Por lo tanto, también se hace necesario estudiar el comportamiento sísmico de las losas y en particular el efecto de la componente vertical del movimiento, el cual es aún menos estudiado.


El terremoto del 27/02/2010 se caracterizó por presentar altas aceleraciones verticales, alcanzando en algunos registros valores muy similares a las aceleraciones horizontales (Barrientos 2010; Boroschek et al. 2010). Sin embargo, la componente vertical de la aceleración es raramente considerada en el diseño sismorresistente de edificios y menos aún en el diseño sismorresistente de losas de hormigón armado.
Estos valores significativos de aceleración vertical merecen un estudio especial respecto a la respuesta que puedan originar en las superestructuras. Uno de los objetivos de este trabajo apunta a poder contribuir a las recomendaciones dadas por los códigos de diseño sismorresistentes respecto a la componente vertical del terremoto de diseño. Por ejemplo, la norma Chilena NCh 433Of96, la cual rige el diseño sismorresistente de edificios, no considera la componente vertical del terremoto de diseño, ya que solo se indica que las solicitaciones sísmicas relevantes provienen de aceleraciones horizontales.
Ecuación diferencial de placas
La expresión que permite estimar la deformación vertical de una placa y que puede ser usada para calcular la flecha de una losa de hormigón a partir de las propiedades geométricas y resistentes de este material y las cargas aplicadas viene dada por (Timoshenko, 1937; Chakraverty 2009):

La expresión (1) es una ecuación diferencial de cuarto orden que relaciona la flecha w con la carga repartida q y las propiedades del material. Donde w es el desplazamiento vertical o flecha, es la masa por unidad de superficie y q es la carga distribuida sobre la losa. El parámetro representa la rigidez flexural de la placa y su expresión está dada por:

donde E es el módulo de Young o de elasticidad del material, ν es la relación de Poisson y t es el espesor de la losa.
Discretizacion de la losa por medio del MEF
Para resolver (1) se ha utilizado el Método de los Elementos Finitos MEF. Para ello se han escogido dos elementos en particular: el elemento rectangular MZC (Melosh 1963, Zienkievicz y Cheung 1967) y el elemento triangular DKT (Discrete Kirchhoff Triangular, ver Batoz et al. 1980). La diferencia que existe entre estos dos elementos, es que el elemento triangular proporciona más versatilidad geométrica que el elemento rectangular, por ejemplo se pueden realizar mallados con geometrías de losas circulares.
Al discretizar la losa y al aplicar el MEF sobre (1), se obtienen matrices de rigidez y masa, tanto para el elemento MZC como para el elemento DKT. Una vez que han sido calculadas las matrices de rigidez y masa tanto para el elemento MZC como para el elemento DKT, la ecuación que permite encontrar los desplazamientos de la losa es,

La expresión (3) es una ecuación de movimiento para sistemas de múltiples grados de libertad. Esta ecuación matricial de movimiento posee la característica de involucrar las matrices de masa [M], amortiguación [C] y rigidez [K] de la losa. Cabe destacar que los coeficientes que componen la matriz de amortiguación no pueden ser obtenidos explícitamente, debido a esto se incorporan en las ecuaciones de movimiento solo a nivel modal (Paz 1992). También se debe mencionar que (3) es un sistema de ecuaciones acophdo, por lo que para resolverlo se hace necesario desacoplarloa través del Método de Superposición Monal (Paz 1992). Al aplicar el Método de Superposición Modal sobre (3) se obtiene la siguiente ecuación modal de movimiento,


La expresión (4) ha sido resuelta por el método Newmark (1959), donde el lado derecho corresponde al registro verticañ de aceleraciones del terremoto 27/02/2010 obtenido en el Colegio Concepción de San Pedro de la Paz (Barrientos 2010).
Matriz de rigidez y masa del elemento MZC
La Figura 2 muestra el elemento rectangular de cuatro nodos denominado MZC, el cual se basa en la siguiente aproximación,



La aproximación (5) ha sido resuelta por Melosh (1963) y Zienkiewicz y Cheung (1967), entregando la siguiente expresión para la flecha w,


donde N son las matrices de funciones de forma y a es el vector de movimiento del elemento e de un nodo i,respectivamente. Las expresiones analíticas de las funciones de forma en coordenadas naturales son,

La matriz de deformación generalizada de flexión se obtiene a partir del vector de deformaciones generalizadas de flexión de la siguiente manera,

con

y

La matriz de rigidez elemental viene dada por:

donde

es la matriz constitutiva de flexión para un material elástico isótropo. La matriz de masa elemental es,

La expresión del vector de fuerzas nodales equivalentes para una carga q uniformemente distribuida sobre el elemento e es:

Las integrales (13) y (15) se calculan por el método de integración numérica de Gauss para dos dimensiones (Oñate 1992), obteniendo así la matriz de rigidez y de masa elemental. Las matrices en coordenadas globales de rigidez y de masa, al igual que el vector en coordenadas globales de fuerzas nodales, se obtienen por medio de un proceso de ensamblaje a través de los vectores de coordenadas y conectividades de cada elemento que compone la malla de elementos rectangulares. Los esfuerzos en puntos nodales se calculan con la siguiente expresión,

Debido a que las tensiones se calculan en los nodos de cada elemento, puede ocurrir que en nodos comunes sus valores varíen considerablemente. Esto se puede solucionar promediando los valores de las tensiones que comparten el mismo nodo, lo cual se denomina alisado de tensiones. Una vez que se ha resuelto (1) y se han utilizado (3) y (4) se pueden obtener las tensiones en la losa a través de (17). Cabe destacar que se utiliza (17) para el cálculo de tensiones del caso estático y sísmico del elemento rectangular.
Matriz de rigidez y masa del elemento DKT
La formulación del elemento DKT se basa en la teoría discreta de Kirchhoff para placas delgadas, la que da origen al elemento triangular discreto de Kirchhoff (Figura 3). El elemento DKT ha demostrado tener un excelente comportamiento con respecto a su convergencia (Batoz et al. 1980).


Para placas delgadas la deformación transversal por corte y, por lo tanto su energía de deformación Us, es despreciable comparado con la energía de deformación por flexión U. Teniendo en cuenta lo anterior, la formulación de la matriz de rigidez del elemento DKT se basa en la siguiente expresión,

donde U es la energía de deformación por flexión, A es (17) el área del plano medio del elemento, εf es el vector de deformaciones por flexión y es la matriz constitutiva. No se considera la deformación por corte ya que para placas delgadas el aporte energético de la deformación por corte transversal es despreciable. Además (18) contiene solamente las primeras derivadas de θx y θy por lo tanto es posible establecer funciones de interpolación que satisfagan los requerimientos de compatibilidad, relacionando la rotación de la normal (θx y θy) a la superficie media con los desplazamientos transversales w. Este objetivo se alcanza con las siguientes consideraciones:
a.
El elemento triangular debe tener solo 9 grados de libertad, los cuales son el desplazamiento w y los giros θxy θy y de los tres nodos esquineros.
b.
Los giros nodales deben ser , de modo que se satisfagan las condiciones de borde del elemento.
c.
Ya que los modelos de elementos de placas delgadas son gobernados por la teoría de Kirchhoff, ésta se puede aplicar sobre cualquier punto discreto del elemento.
d.
La compatibilidad de los giros θx y θy no debe perderse.
El punto de partida para obtener la matriz de rigidez del elemento de 3 nodos DKT, es el elemento triangular de placa de Reissner-Mindlin de 6 nodos mostrado en la Figura 4, sometido a varias condiciones y/o consideraciones (Batoz et al. 1980).


Utilizando las consideraciones anteriores se puede establecer una interpolación del desplazamiento w (flecha) y los giros nodales θx y θy a través de,

donde u es el vector de movimientos y a es el vector de movimientos de cada grado de libertad del elemento,

N es la matriz de funciones de forma del elemento DKT que aproximan o interpolan los desplazamientos y giros nodales. La matriz de deformación generalizada de flexión se obtiene a partir del vector εf y (21) como,

donde Bf es la matriz de deformación generalizada del elemento definida por:

con 2A = x31y12 - xl2y3l y teniendo en cuenta que A es el área del elemento triangular. Las derivadas de Nθx yNNθy con respecto a ξ y η se pueden encontrar explícitamente en Batoz et al. (1980), al igual de los coeficientesxij 1ij La expresión para encontrar la matriz de rigidez del elemento DKT es la siguiente,

y para la matriz de masa se tiene de acuerdo a Luo y Hutton (2002) que,

Finalmente, la expresión del vector de fuerzas nodales equivalente para una carga q uniformemente distribuida sobre el elemento es,

donde q es la carga repartida sobre el elemento triangular y A es el área del elemento triangular. Al igual que para el elemento rectangular, las ecuaciones (25) y (26) se calculan por el método de integración numérica de Gauss para dos dimensiones (Oñate 1992), obteniendo así la matriz de rigidez y de masa elemental. Las matrices en coordenadas globales de rigidez y de masa, al igual que el vector en coordenadas globales de fuerzas nodales, se obtienen con el proceso de ensamblaje a través de la matriz de transformación de coordenadas, y también, por los vectores de coordenadas y conectividades de cada elemento que compone la malla de elementos triangulares. Los esfuerzos en los puntos nodales se calculan con la siguiente expresión,

Nuevamente se tiene que las tensiones calculadas en los nodos de cada elemento puede causar que en nodos comunes sus valores varíen. Esto se puede solucionar usando alisado de tensiones. Como también se dijo anteriormente, una vez que se ha resuelto (1) y se han utilizado (3) y (4) se pueden obtener las tensiones sobre la losa a través de (28). Cabe destacar que (28) se ha utilizado para el cálculo de tensiones en el caso estático y sísmico del elemento triangular.
Validación numérica
Para poder utilizar el código programado, primero se ha validado, es decir, se comparan los resultados del programa con las soluciones analíticas o exactas conocidas, las cuales se han basado en la teoría de la elasticidad lineal, y específicamente en la teoría de placas y en la teoría de dinámica estructural. Debido a que el programa involucra el comportamiento sísmico de una estructura (respuesta de una losa bajo la acción de la componente vertical de aceleraciones), ha sido necesario obtener las matrices de rigidez y masa para la resolución de la ecuación de movimiento. La matriz de rigidez de la losa se ha validado con el caso de una losa empotrada en todos sus bordes y sometida a una carga repartida q. Y la matriz de masa, mediante la comparación de las frecuencias naturales entregadas por una expresión analítica proveniente de la teoría de vibración libre de placas con las que entrega el programa elaborado con el MEF. Finalmente, la integración en el tiempo del registro de aceleraciones, se ha validado comparando la respuesta de un oscilador de un grado de libertad, obtenida mediante su solución analítica con la respuesta que otorga el método numérico de Newmark (1959).
Validación del caso estático
En este caso se estudia la convergencia de la flecha central de una losa cuadrada con todos sus bordes empotrados. La solución analítica del problema, según la teoría de Kirchhoff para placas delgadas, está dada por la siguiente expresión (Novoa 2005),

donde w es la flecha central, q es la carga repartida, L es el largo, E es el módulo de elasticidad estático, ν es la razón de Poisson y t es el espesor de la losa. Las propiedades de la losa con las que se ha realizado el análisis y la comparación de los resultados se muestran la Tabla 1 y en la Figura 5.




De la Figura 5 se puede deducir que tanto el elemento MZC como DKT entregan aproximaciones razonables de la solución analítica de la flecha central de la losa. En las Figuras 6 y 7 se muestran los resultados de las deformadas de la losa con el elemento MZC ν DKT.




Validación de caso dinámico
La validación del caso dinámico se ha llevado acabo a través de la comparación de las 6 primeras frecuencias naturales obtenidas de la teoría de vibración libre de placas y las obtenidas a través del programa. Cabe destacar que la losa que se ha utilizado para la validación es una losa cuadrada simplemente apoyada. Existe una expresión teórica para las frecuencias naturales válida para estos tipos de losas dada por (Timoshenko, 1937; Chakraverty 2009):

donde es la rigidez flexural de una placa, ρ es la densidad del material de la placa, t es el espesor de la placa,a y b son las longitudes de los lados de la placa, y finalmente, m indica la numeración de los modos de vibrar de la placa. La Tabla 1 presenta los valores de las propiedades de la losa usados en los cálculos numéricos.
Las Tablas 2, 3 y 4 presentan los resultados obtenidos de las frecuencias naturales, en donde la letra m indica el número de frecuencia natural calculada.




De las Tablas 2, 3 y 4 es posible apreciar que los valores de las frecuencias naturales de las placas calculadas con el elemento MZC y DKT son aproximaciones razonables del valor teórico otorgado por (30).
Validación de la integración en el tiempo del registro de aceleraciones
Se ha analizado la respuesta de un oscilador de 1 grado de libertad por medio de la solución exacta y posteriormente, mediante la implementación del método numérico de Newmark (1959). El oscilador que se ha considerado posee las propiedades y valores mostrados en la Tabla 1.
La solución analítica de este oscilador está dada por la siguiente expresión:

donde ω0 es la frecuencia de la fuerza externa. Cabe destacar que esta expresión corresponde a la respuesta de un oscilador de 1 grado de libertad sometido a una carga sinusoidal monotónica, la cual no se puede extender a un registro aleatorio como lo es un registro sísmico. Finalmente, en la Figura 8 se muestra la comparación de ambas respuestas, en donde se aprecia que el resultado numérico es prácticamente igual a la solución exacta.


Análisis comparativo
El análisis comparativo consiste en la obtención de las tensiones σχ, σy y τχy, producto de la aplicación de fuerzas estáticas y sísmicas sobre 2 losas distintas. En el caso sísmico, la carga sobre la losa corresponde al registro de aceleraciones verticales del terremoto del 27/02/2010 (Colegio Concepción, San Pedro de la Paz) y para el caso estático es la carga que especifica la norma chilena NCh 1537Of.86, la cual en su tabla 3, especifica que la sobrecarga de uso uniformemente distribuida para pisos es de 2 kPa.
Para la comparación de ambos casos en el análisis sísmico de la losa, se ha considerado el punto 5.5.1 de la norma NCh 433Of96, la cual establece que para el cálculo de masas se deben considerar las cargas permanentes más un porcentaje de la sobrecarga de uso, que no podrá ser inferior a 25 % en construcciones destinadas a la habitación o al uso público. La Tabla 6 resume las características del material y condiciones de apoyo de la losa. Los valores adoptados permanecen constantes durante el análisis. Cabe destacar que se consideró en el análisis el módulo de elasticidad estático del hormigón. Las características geométricas de las losas que se han analizado se muestran en las Figuras 9 y 10. La Tabla 7 entrega las cargas estáticas uniformemente repartida sobre la losa, sobrecarga, peso propio y carga total.








La carga sísmica sobre la losa corresponde al registro de aceleraciones verticales del terremoto del 27 de febrero del 2010 mostrado en la Figura 11.


Se debe tener en cuenta que en este trabajo, la masa de la losa corresponde a su masa como tal más un 100 % de la sobrecarga de uso (según NCh433Of96), la cual es de 2 kPa. En este caso se ha utilizado el 100 % del valor de la sobrecarga debido a que es el caso más desfavorable.
Las Figuras 12, 13 y 14 muestran las distribuciones de tensiones a partir del elemento MZC y DKT para la losa 1 bajo carga estática.


Las Figuras 15, 16 y 17 muestran las envolventes de tensiones a partir del elemento MZC y DKT para la losa 1 bajo carga sísmica. Notar que la escala en la dimensión horizontal χ ha sido reducida y que se han separado las tensiones positivas y negativas.


Las tensiones máximas positivas sobre la losa 1 para el caso estático son mayores a las del caso sísmico. Sin embargo, se producen tensiones máximas negativas que son mayores a las del caso estático (-σx, -σy, -τxy).
Las Figuras 18, 19 y 20 muestran las distribuciones de tensiones a partir del elemento MZC y DKT para la losa 2 bajo carga estática.


Las Figuras 21, 22 y 23 muestran las envolventes de tensiones a partir del elemento MZC y DKT para la losa 2 bajo carga sísmica. Notar que nuevamente se ha cambiado la escala y se han incluído las tensiones negativas.



En este caso tanto las tensiones máximas positivas como negativas que produce la aplicación del registro vertical del terremoto son siempre mayores a las que produce la aplicación de una carga estática vertical. Excepto en el caso de σx positivo con elementos triangulares.
Los resultados con tensiones negativas de -4 MPa podrían explicar las grietas y desprendimientos de hormigón especialmente en las uniones de las losas con los muros o vigas tal como se mostró en la Figura 1.
Las periodos y frecuencias fundamentales calculados para la losa 1 y 2 son 0.031 s (198.9 rad/s) y 0.104 s (60.4 rad/s) respectivamente. De esta manera el periodo fundamental de la losa 2 conlleva a aceleraciones espectrales mayores que la losa 1. Esto es porque la losa 2 es más flexible, ya que posee un periodo fundamental mayor. Por lo tanto, este aumento en la aceleración espectral se debe a que este terremoto en particular (en ese rango de periodos la aceleración espectral es ascendente), afecta mayormente a la losa 2 que a la losa 1. Esto explica que las tensiones provocadas por el registro vertical de aceleraciones aplicadas sobre la losa 2 superen a las que origina la carga estática aplicada sobre la misma losa, y para la losa 1 resulten menores.
Conclusiones y recomendaciones
Se han utilizado las teorías de placas basadas en la elasticidad lineal y la del movimiento de estructuras de 1 GDL para caracterizar el comportamiento sísmico de una losa de hormigón armado. Para ello se ha implementado un programa computacional en código MATLAB, el cual permite obtener desplazamientos verticales y tensiones ±σx, ±σy y ±τxy, producto de la aplicación de cargas estáticas y sísmicas. Se utilizaron los elementos finitos MZC y DKT, encontrándose no mayor diferencia en los resultados obtenidos al usar uno u otro elemento.
Las tensiones máximas en losas cuadradas y rectangulares empotradas surgen en los bordes de éstas, tanto para el caso estático como para el sísmico. La forma general de las distribuciones de tensiones en losas cuadradas y rectangulares son muy similares entre si, tanto para el caso estático como para el sísmico. Por otro lado, las fibras de una losa que habitualmente se encuentran en compresión (o tracción) bajo carga estática vertical, cambian temporalmente a tracción y compresión (compresión y tracción) durante la aplicación del registro de aceleraciones vertical del terremoto.
Los registros que posean una componente de aceleración vertical importante, como el del 27/02/2010 (San Pedro de la Paz), pueden provocar tensiones superiores a las que provocan cargas estáticas, como por ejemplo las cargas de diseño de la norma NCh 1537Of.86. Es por ello que se debería considerar la componente vertical en los análisis sísmicos, ya que sismos con aceleraciones verticales mayores a la del terremoto del 27/02/2010 podrían, eventualmente, provocar tensiones superiores a las que se consideran en los análisis para el diseño, tal como ocurre en las losas analizadas en este trabajo. Se recomienda en futuras investigaciones incluir a los muros en el análisis de aceleración vertical para así evaluar el efecto de interacción losa muro.























































































































No hay comentarios:

Publicar un comentario