Source: https://www.scribd.com/doc/146403662/No-lineal
Timestamp: 2017-10-23 19:57:49+00:00

Document:
Juan Tomás Celigüeta Departamento de Ingeniería Mecánica
7 Matriz de rigidez tangente 13 RESOLUCIÓN DE LAS ECUACIONES INCREMENTALES 13. Barra apoyada .4 Métodos restringidos 13.2 Deformaciones unitarias 12.3 Deformaciones virtuales 11.5 Trabajo virtual interior 12. FORMULACIÓN CO-ROTACIONAL 10. FORMULACIÓN CO-ROTACIONAL 11.3 Matriz de rigidez tangente 11 ELEMENTO VIGA PLANA.2 Vector de fuerzas interiores 9.2 Método de Newton-Raphson 13.3 8.2 Vector de fuerzas interiores 10.4 Formulación isoparamétrica 10 ELEMENTO BIARTICULADO. FORMULACIÓN LAGRANGIANA TOTAL 9.8.deslizante .1 Ejemplo 1.1 Deformación unitaria 9.5 8.1 Campo de deformaciones 12.6 Trabajo virtual Ecuación de equilibrio Linealización de las ecuaciones de equilibrio Ecuaciones de equilibrio incrementales Formulación isoparamétrica Fuerzas nodales equivalentes a las fuerzas exteriores 39 40 40 43 43 45 46 46 48 48 50 52 53 54 54 57 57 58 59 60 60 63 63 63 65 67 68 68 68 71 72 73 74 75 77 79 84 84 9 ELEMENTO BIARTICULADO.1 Deformación unitaria 10.1 8.3 Matriz de rigidez tangente 9.2 8.5 Matriz de rigidez tangente 12 FLEXIÓN DE PLACAS. FORMULACIÓN LAGRANGIANA TOTAL 12.5 Criterios de convergencia 14 MÉTODO DE LA LONGITUD DEL ARCO 15 EJEMPLOS ESTÁTICOS 15.4 8.3 Variación de la deformación unitaria 12.4 Trabajo virtual 11.4 Deformaciones unitarias de cortadura 12.6 Vector de fuerzas interiores 12.3 Método de Newton modificado 13.2 Deformación y momentos de flexión 11.1 Deformación axial y esfuerzo axial 11.1 Método incremental puro 13.
Formulación lagrangiana total 16.5 Métodos implícitos de integración de paso simple 16.3 Ejemplo 3.2 Traza 20.5 Teoremas de integración 21 ANEJO 3.1 Principio del trabajo virtual en dinámica 16. PROCEDIMIENTOS MATLAB . Celosía Ejemplo 5. Voladizo muy flexible.3 Método explícito basado en diferencias centrales 16. Voladizo muy flexible Ejemplo 4. NOTACIÓN 20 ANEJO 2.3 Gradiente 20. Pórtico biarticulado 87 88 90 91 94 94 94 95 96 99 102 103 103 104 106 108 109 110 110 111 111 112 113 115 16 DINÁMICA 16.6 Criterios de convergencia 17 EJEMPLOS DINÁMICOS 17.2 15.4 Divergencia 20. PRELIMINARES MATEMÁTICOS 20. 17.4 15. Cable pretensado 18 BIBLIOGRAFÍA 19 ANEJO 1.4 Estabilidad del método de diferencias centrales 16. Barra apoyada – deslizante 17.5 Ejemplo 2.2 Ecuaciones de equilibrio.3 15.1 Resumen de álgebra de vectores y tensores 20.2 Ejemplo 2.1 Ejemplo 1.15. Barra deslizante apoyada elásticamente Ejemplo 3.
Es necesario seguir un proceso de carga incremental. La naturaleza no lineal del fenómeno hace que no pueda calcularse en general la situación deformada final en un sólo paso. la ecuación anterior proporciona su trayectoria temporal. de tal manera que no puede aceptarse la hipótesis de que la posición final deformada coincide con la posición inicial. t ) (1) Para cada instante de tiempo. la respuesta final depende del orden de aplicación de las mismas y se hace necesario el proceso de carga paso a paso. y determinando la respuesta para cada uno de esos incrementos. La figura 1 muestra el sólido referido a un sistema de coordenadas cartesiano. Cada partícula queda definida por unas coordenadas x i agrupadas en un vector X . por incrementos. A esta no linealidad de origen geométrico se puede añadir la no linealidad debida a la ecuación constitutiva del material si es que este fenómeno se pone de manifiesto en el proceso. La necesidad de un proceso de carga incremental y de un parámetro al cual referir el mismo es importante asimismo cuando existen condiciones de carga diversas que pueden aplicarse en diferente orden. Al ser el sistema no lineal. aplicando las cargas finales paso a paso. ni siquiera siguiendo un proceso iterativo. Configuración en un instante cualquiera t. 1 . aplicando la totalidad de la carga. Por lo tanto la única diferencia entre los casos estático y dinámico está en la consideración o no de las fuerzas de inercia. En este contexto la respuesta del sólido es altamente no lineal pues por una parte no se conoce la posición deformada final y por otra la presencia de grandes deformaciones implica el uso de medidas de la deformación adecuadas que son esencialmente no lineales. y en ella se identifican:   Configuración inicial en t=0.Introducción al análisis de estructuras con no linealidad geométrica 1 INTRODUCCIÓN Se estudia la deformación de un medio continuo deformable bajo la acción de cargas exteriores que provocan grandes deformaciones. aunque no se trate más que de un parámetro arbitrario. La posición de cada partícula queda ahora definida por unas coordenadas x El movimiento entre t=0 y t se puede representar matemáticamente como una función: x = φ(X. pero por comodidad se le denominará así. Para identificar los distintos pasos del proceso se empleará un parámetro de tiempo t. la función φ define una transformación de coordenadas entre las dos configuraciones espaciales. En principio no se considerará aquí esta no linealidad debida al material. al cual se referirán todos los incrementos de carga y las distintas configuraciones deformadas. Para una partícula cualquiera. En el caso de que las cargas sean estáticas no tiene sentido hablar del parámetro tiempo en el sentido que tiene en dinámica. inicial y final.
El planteamiento Lagrangiano es adecuado al estudio de la mecánica de sólidos. Un cambio en el tiempo en las ecuaciones anteriores implica que una determinada posición está ocupada por 2 . t 1. t ) (3) En este planteamiento se persigue una determinada posición en el espacio y se determina la posición inicial que tenían las partículas que pasan por dicha posición. t ) Las deformaciones se obtienen con respecto a esa posición deformada: u(X. es decir: X = φ−1(x. Configuración inicial y deformada. se trata de expresar las coordenadas finales de una partícula en función de sus coordenadas iniciales: x = φ(X. el comportamiento se refiere a la posición final x ocupada por una partícula y se trata de obtener la posición inicial que ocupaba dicha partícula. o material.Introducción al análisis de estructuras con no linealidad geométrica (X. t ) La deformación del sólido se trata de expresar asimismo en función de dichas coordenadas: u(X.t) u t=0 X x x2 x3 x1 Figura 1. t ) = φ(X. • En el planteamiento Euleriano.1 Planteamientos material y espacial La resolución de un problema no lineal puede abordarse genéricamente de dos maneras distintas. o espacial. dependiendo de a qué sistema de coordenadas se refieran las magnitudes fundamentales involucradas en el proceso. Un cambio en el tiempo en las ecuaciones anteriores implica que la misma partícula ocupa una posición diferente. • En el planteamiento Lagrangiano. t ) − X (2) En este planteamiento por lo tanto se persigue el movimiento de una misma partícula material cuya posición inicial X se conoce y se trata de obtener su posición final. en el que es necesario incluir alguna ecuación constitutiva del comportamiento de las partículas del material. t ) = x − φ−1(x.
t ) Sin embargo. t ) Figura 2. como la temperatura T. su valor se describe con respecto a la posición inicialmente ocupada por una partícula. si consideramos una magnitud escalar cualquiera. se describe con respecto a una posición actual en el espacio. • Por ejemplo. el valor la temperatura T. Planteamientos material y espacial. en el que no interesa tanto la evolución de las partículas sino la distribución espacial de las magnitudes. en el planteamiento Euleriano espacial. en el planteamiento Lagrangiano.Introducción al análisis de estructuras con no linealidad geométrica una partícula diferente. y se obtiene como la temperatura en esa misma partícula a medida que pasa el tiempo (y posiblemente también cambie la posición de la partícula): T = T (X. El planteamiento Euleriano es adecuado a problemas de mecánica de fluidos. 3 . y se obtiene como la temperatura en esa posición a medida que pasa el tiempo (y posiblemente también pasen distintas partículas por ese punto): T = T (x.
y establece la relación entre los elementos diferenciales de coordenadas en dos configuraciones. es: dx = F dX Por lo tanto el tensor gradiente de la deformación se define como el gradiente de las coordenadas deformadas x respecto de las coordenadas iniciales X. Así entre la configuración inicial X y deformada x. Denominamos ∇0 a dicho gradiente: F = ∇0 x = Cada uno de sus términos es: Fij = ∂x ∂X (4) ∂x i ∂X j Es un tensor definido positivo. Consideremos el diferencial de volumen en el 4 . como se demuestra al estudiar la transformación del volumen.1 GRADIENTES DE DEFORMACIÓN Tensor gradiente de deformación Es la magnitud fundamental en la definición de la deformación de un sólido. También se puede poner en notación de matrices como: F = ∇0 xT ( ) T (5) Siendo ∇0 la representación como vector del operador gradiente respecto de las coordenadas iniciales: ⎡ ∂ ⎤ ⎢ ⎥ ⎢ ∂X 1 ⎥ ⎢ ⎥ ⎢ ∂ ⎥ ⎥ ∇0 = ⎢ ⎢ ∂X ⎥ 2 ⎢ ⎥ ⎢ ∂ ⎥ ⎢ ⎥ ⎢ ∂X ⎥ ⎣ 3⎦ • La relación inversa entre ambos elementos diferenciales permite introducir el tensor gradiente de deformación inverso F-1: dX = F−1 dx Cuyo valor es: (6) F−1 = T ∂X = (∇XT ) ∂x Siendo ∇ el operador gradiente respecto de las coordenadas deformadas x. • El determinante del tensor gradiente de deformación F=|F| establece la relación entre los diferenciales de volumen en los estado 0 y t.Introducción al análisis de estructuras con no linealidad geométrica 2 2.
dX2. y cuya dirección n0 viene dada por su producto vectorial: dA0 = n 0 dA0 = dX1 × dX2 El diferencial de área en el estado deformado es: dA = n dA = dx1 × dx2 = (FdX1 ) × (FdX2 ) Empleando la propiedad del producto vectorial (A a) × (A b) = A A−T (a × b) se obtiene: dA = F F−T (dX1 × dX2 ) = F F−T dA0 dA = J F−T dA0 2.Introducción al análisis de estructuras con no linealidad geométrica estado inicial dv0 formado por un paralelepípedo cuyos lados están definidos por tres vectores diferenciales dX1. 3 dv = dx1 ⋅ (dx2 × dx 3 ) = FdX1 ⋅ (FdX2 × FdX3 ) Empleando la siguiente propiedad del producto mixto: Aa ⋅ (Ab × Ac) = A a ⋅ (b × c) se puede poner: dv = FdX1 ⋅ (FdX2 × FdX3 ) = F dX1 ⋅ (dX2 × dX3 ) = F dv 0 dv = F dv 0 (7) Como ambos diferenciales de volumen siempre deben ser positivos.2 Tensores gradiente de desplazamientos El tensor gradiente de desplazamientos lagrangiano (o material) H se define como el gradiente de los desplazamientos u respecto de las coordenadas iniciales: H = ∇0 u = ∂u ∂X (8) En su representación como matriz. el determinante del tensor gradiente de deformación también lo es. los términos de este tensor son: 5 . dX3. Consideremos un elemento diferencial de área en el estado inicial dA0: es una cantidad vectorial cuyo módulo dA0 es igual al área de un paralelogramo definido por dos vectores dX1 y dX2. El valor de este diferencial de volumen es el producto mixto de los tres vectores: dv 0 = dX1 ⋅ (dX2 × dX3 ) Los vectores se transforman en el estado deformado a: dxi = F dXi El diferencial de volumen en el estado deformado es: i = 1. • La transformación que sufre el diferencial de área a consecuencia de la deformación se establece mediante el tensor gradiente inverso.
3 ∑λ α n α NT α Esto indica la naturaleza de F en el sentido de que involucra a las direcciones principales en los estados inicial y deformado. con lo que λα representa el alargamiento en la dirección α: λα = dlα dLα • Es instructivo expresar el tensor F en función de las direcciones principales en ambos estados. • Considerando un vector dXα orientado en la dirección de un eje cualquiera α.Introducción al análisis de estructuras con no linealidad geométrica dl 1 1d L1 1d L1 dL 1 Figura 3. Para ello sustituimos el valor de U: F = RU = Considerando que n = R N. de módulo dLα: dXα = dLα Nα En el estado deformado será: dxα = FdXα = R U Nα dLα Pero se cumple que U Nα = λα Nα luego: dxα = λαdLα R Nα dxα = λαdLα n α Luego el módulo del vector deformado es dlα = λαdLα . Descomposición polar del tensor gradiente de deformación.3 ∑λ α R Nα NT α F= α=1. α =1. 9 .
• La variación de los tensores gradiente de desplazamientos al variar los desplazamientos es sencillamente: δH = ∂δ u ∂X δH ij = ∂δui ∂X j ∂δui ∂x j Si se emplea la definición en función de las coordenadas deformadas x la variación es: δ Ht = ∂δ u ∂x δH tij = • La variación del tensor gradiente de la deformación es inmediata: δF = δH Esta expresión se puede poner en función de las coordenadas deformadas x efectuando la derivación en cadena: δFij = ∂δ x i ∂δui ∂δui ∂x 1 ∂δui ∂x 2 ∂δui ∂x 3 = = + + ∂X j ∂X j ∂x 1 ∂X j ∂ x 2 ∂X j ∂x 3 ∂ X j δ F = δ Ht F Por lo tanto la variación del gradiente de deformación se puede poner: 10 .6 Variación de los tensores gradiente Consideramos un cuerpo en un estado deformado t.Introducción al análisis de estructuras con no linealidad geométrica 2. sometido a unos desplazamientos u=x-X. Se aplica una variación virtual δu a dichos desplazamientos. la cual produce una variación de los gradientes de deformaciones.
que corresponde al tensor de deformaciones unitarias empleado en el análisis con pequeñas deformaciones. Es igual a la parte simétrica del tensor gradiente de desplazamientos Ht evaluado respecto de las coordenadas deformadas: ε= 1 (Ht + HT t ) 2 (16) Se puede representar en forma vectorial como: ⎧ ε11 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε ⎪ ⎪ 22 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε ⎪ ⎪ 33 ⎪ ⎪ ε =⎨ ⎬ ⎪ 2ε12 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε 2 13 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε 2 ⎪ ⎪ 23 ⎪ ⎪ ⎩ ⎭ Dado su valor.1 DEFORMACIONES UNITARIAS Tensor infinitesimal de deformaciones unitarias Sea u el campo de deformaciones existente en el sólido en el instante t.Introducción al análisis de estructuras con no linealidad geométrica 3 3. la expresión detallada para el caso de dos dimensiones es: ⎧ ∂u1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎧ ⎫ ⎪∂ ⎪ x ⎪ ⎪ ∂ u 1⎪ 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂u1 ⎪ ∂x 1 ⎡ 1 0 0 0⎤ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎢ ⎥ ⎪ ⎪ ⎪ ⎪ ⎪ ∂x 2 ⎪ ⎪ ⎥⎪ ⎪ ∂u 2 ⎪ ⎢ ε =⎨ ⎬ = ⎢ 0 0 0 1⎥ ⎨ ⎬ ⎪ ⎪ ∂u 2 ⎪ ∂x 2 ⎢ ⎥⎪ ⎪ ⎪ ⎪ ⎪ ⎢ ⎥ ⎪ ⎪ ⎪ ⎪ 0 1 1 0 ⎪ ⎪ ⎪ ⎢ ⎥ ∂ x ∂u1 ∂u2 ⎪ ⎣ ⎦⎪ 1⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ + ⎪ ⎪ ⎪ ⎪ ∂ u ⎪ ⎪ ⎪ ⎪ ∂ ∂ x x 2 1⎭ ⎪ ⎪ ⎩ 2 ⎪ ⎪ ⎪ ⎪ ⎪ ∂ x ⎪ 2⎪ ⎪ ⎩ ⎭ (17) • La variación del tensor infinitesimal de deformaciones unitarias al variar los desplazamientos es: δε = 1 (δ Ht + δ HT t ) 2 11 . se puede expresar siempre en la forma: ⎧ε ⎪ ⎫ ⎪ 11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε =⎪ ⎨ ε22 ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ε 2 ⎪ 12 ⎪ ⎪ ⎪ ⎩ ⎭ ε = AC Ht Siendo AC una matriz constante. Se define el tensor de deformaciones unitarias infinitesimales como: ∂u j ⎞ 1 ⎛ ∂u i ⎟ ⎟ ⎜ + εij = ⎜ ⎟ ⎜ ⎟ 2 ⎜∂ ∂x i ⎠ ⎝ xj (15) Es un tensor lineal. Por ejemplo.
12 .2 Tensor de deformaciones unitarias de Green-Lagrange (18) El cuadrado de la distancia entre dos puntos infinitamente próximos en los estados inicial y deformado es: (ds 0 ) 2 = dXT dX (ds ) = dxT dx El cambio en esta distancia al cuadrado se puede expresar como: 2 (ds ) − (ds 0 ) = dxT dx − dXT dX = dXT FT FdX − dXT dX = dXT FT F − I dX 2 2 ( ) Este cambio de la distancia al cuadrado referido a la distancia inicial define el tensor de deformación unitaria de Green – Lagrange entre los estados 0 y t: (ds ) − (ds 0 ) = 2 dXT E dX 2 2 E= 1 T 1 F F − I = (C − I) 2 2 ( ) (19) Este tensor se puede expresar en función de los desplazamientos sustituyendo el valor del tensor gradiente de deformaciones F: E= 1⎡ T (I + H) (I + H) − I⎤⎥⎦ ⎢ ⎣ 2 1⎡ H + HT + HT H⎤⎥ ⎢ ⎣ ⎦ 2 (20) E= Sustituyendo el valor del gradiente de desplazamientos H: T T ⎛ ∂u ⎞ 1 ⎡⎢ ∂u ⎛ ∂u ⎞ ∂u ⎤⎥ ⎟ ⎟ ⎜ ⎜ E= +⎜ + ⎟ ⎟ ⎜ ⎟ ⎟ ∂X ⎥ ⎝ ∂X ⎠ 2 ⎢⎣ ∂X ⎝ ∂X ⎠ ⎦ (21) ⎞ ⎛ ∂u ∂ u ⎞ ∂u j 1 ⎛ ∂ui ⎟ ⎜ k k ⎟ ⎟ ⎟ ⎜ ⎜ Eij = ⎜ + + ⎟ ⎟ ∑ ⎜ ⎜ ⎟ ⎟⎠ 2⎜ ∂X i ⎟ ⎝ ∂X i ∂ X j ⎠ k ⎜ ⎝ ∂X j Los dos primeros sumandos corresponden al tensor infinitesimal lineal de pequeñas deformaciones. y el último término corresponde a los términos no lineales habituales en grandes deformaciones.Introducción al análisis de estructuras con no linealidad geométrica Cuyos términos valen: ∂δu j ⎞ 1 ⎛ ∂δui ⎟ ⎟ ⎜ δεij = ⎜ + ⎟ ⎜ ⎟ 2⎜ ∂x i ⎠ ⎝ ∂x j En forma vectorial esta variación se puede poner: δ ε = AC δ Ht 3.
en la forma: ⎧ ⎫ ⎡ ⎪ ⎪ ∂u1 ⎪ ⎪ ⎢ ∂u1 ⎪ ⎪ ⎢ ∂X 1 ⎪ ⎪ ∂X 1 ⎪ ⎪ ⎢ ⎪ ⎪ ⎪ ⎪ ∂u 2 ⎪ ⎪ 1 ⎢⎢ E=⎨ ⎬+ ⎢ 0 ⎪ ⎪ ∂X 2 2⎢ ⎪ ⎪ ⎪ ⎪ ⎢ ⎪ ∂u ⎪ ∂u 2 ⎪ ⎪ 1 ⎢ ∂u1 ⎪ ⎪ + ⎪ ⎪ ⎢ ∂X ⎪ ⎪ ∂X 2 ∂ X 1 ⎪ ⎪ ⎩ ⎭ ⎣ 2 ∂u 2 ∂X1 ⎧ ∂u1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎤ ⎪∂X1 ⎪ ⎪ ⎪ 0 ⎥⎪ ⎪ ⎪ ⎥⎪ ∂u1 ⎪ ⎪ ⎪ ⎥⎪ ⎪ ⎪ ⎥ ∂X 2 ⎪ ⎪ ∂u 2 ⎥ ⎪ ⎪ ⎪ ⎨ ⎬ ⎪ ∂ u ∂X 2 ⎥⎥ ⎪ 2 ⎪ ⎪ ⎪ ⎪ ∂X 1 ⎪ ⎪ ∂u2 ⎥⎥ ⎪ ⎪ ⎪ ⎪ ⎪ ⎥ ⎪ ∂u 2 ⎪ ∂X 1 ⎦ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂ X ⎪ ⎪ ⎩ 2⎪ ⎭ 0 ∂u1 ∂X 2 ∂u1 ∂X 1 0 ∂u 2 ∂X 2 La estructura del primer término. Agrupando los términos lineales y los no lineales se puede poner.2. Entre t y t+Δt se aplica una rotación de sólido rígido definida por una matriz R. El tensor gradiente de deformación en el nuevo estado es: Ft +Δt = R Ft El tensor de deformaciones unitarias de Green – Lagrange en el nuevo estado es: Et +Δt = T T 1 ⎡ t +Δt T t +Δt 1 1 − I⎤⎥ = ⎡⎢(Ft ) RT R Ft − I⎤⎥ = ⎡⎢(Ft ) Ft − I⎤⎥ = Et F ) F ( ⎢ ⎦ 2⎣ ⎦ ⎦ 2⎣ 2⎣ 3. sometido a un tensor de deformación Ft . permite expresarlo como: 13 .Introducción al análisis de estructuras con no linealidad geométrica • El tensor de Green – Lagrange es invariante ante rotaciones de sólido rígido. para el caso de 2 dimensiones. Sea un sistema en un instante t. para problemas de 2 y 3 dimensiones (la barra sobre el símbolo indica una representación como vector): ⎧E ⎪ ⎫ ⎪ 11 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ E=⎪ ⎨ E22 ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 2 E ⎪ 12 ⎪ ⎪ ⎪ ⎩ ⎭ ⎧ E11 ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ E 22 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ E 33 ⎪ ⎪ ⎪ ⎪ E=⎪ ⎨ ⎬ ⎪ 2 E 12 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 2 E 13 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 2E 23 ⎪ ⎪ ⎪ ⎩ ⎭ (22) Donde se han multiplicado por 2 los términos fuera de la diagonal para poder sustituir el producto contracto de tensores de orden 2 por el producto escalar de vectores. que es lineal.1 Expresión vectorial del tensor de Green – Lagrange El tensor de Green-Lagrange se puede expresar en forma de vector en la forma siguiente.
Se aplica una variación virtual δ u a dichos desplazamientos.3 0 ∂u 2 ∂X1 0 ∂u1 ∂X 2 ∂u1 ∂X 1 ∂u 2 ∂X 2 ⎤ 0 ⎥ ⎥ ⎥ ∂u2 ⎥⎥ ∂X 2 ⎥⎥ ∂u2 ⎥⎥ ∂X1 ⎥⎦ (26) Variación del tensor de Green – Lagrange Consideramos un cuerpo en un estado deformado t. Precisamente la dependencia de A de la deformación es el origen de la no-linealidad del problema.Introducción al análisis de estructuras con no linealidad geométrica ⎛ ⎡ ⎜ ⎢ ∂u1 ⎜ ⎜ ⎢ ∂X 1 ⎜⎡ ⎤ 1 0 0 0 ⎜ ⎢ ⎢ ⎥ ⎜ ⎜⎢ ⎥ 1 ⎢⎢ E=⎜ 0 0 0 1 ⎜ ⎢ ⎥+ ⎢ 0 ⎜⎢ 2⎢ ⎥ ⎜ ⎜ ⎢ 0 1 1 0⎥ ⎜ ⎢ ∂u ⎜⎢ ⎥⎦ ⎣ ⎜ ⎢ 1 ⎜ ⎜ ⎢ ⎜ ⎝ ⎣ ∂X 2 0 ∂u1 ∂X 2 ∂u1 ∂X 1 ∂u 2 ∂X 1 0 ∂u 2 ∂X 2 ⎧ ∂u1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎤⎞ ⎪∂ X ⎪ 1 ⎟ ⎪ ⎪ ⎥ ⎟ 0 ⎟⎪ ⎪ ⎪ ⎥⎟ ⎟ ⎪ ∂u1 ⎪ ⎪ ⎥ ⎟⎪ ⎪ ⎟ ⎪ ⎥ ⎟ ∂ X ⎪ 2⎪ ⎪ ∂u 2 ⎥ ⎟ ⎪ ⎪ ⎟ ⎨ ⎬ ⎟ ⎥ ⎟⎪ ∂X 2 ⎥ ⎟ ⎪ ∂u2 ⎪ ⎪ ⎟⎪ ⎪ ⎪ ⎪ ⎟ ∂ X ∂u2 ⎥⎥ ⎟ ⎪ 1⎪ ⎟ ⎪ ⎪ ⎟ ⎟⎪ ⎟ ⎪ ∂u 2 ⎪ ∂X1 ⎦⎥ ⎠ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎩ ⎭ (23) Esta expresión se puede poner en forma compacta definiendo dos matrices y recordando la forma vectorial del tensor de desplazamientos H: ⎛ 1 ⎞ E=⎜ A + A(H)⎟ ⎟H ⎜ ⎝ C 2 ⎠ (24) • La matriz AC es constante. la variación se puede poner: δE = 1 δ HT F + FT δ H 2 ( ) (27) 14 . Según (19). sus términos son sencillamente una reordenación de los términos de H. con tamaño 3x4 en 2 dimensiones y 6x9 en 3 dimensiones: ⎡ 1 0 0 0⎤ ⎢ ⎥ ⎢ ⎥ AC = ⎢ 0 0 0 1 ⎥ ⎢ ⎥ ⎢ 0 1 1 0⎥ ⎥⎦ ⎣⎢ (25) • La matriz A(H) depende del vector gradiente de desplazamientos H. sometido a unos desplazamientos ui = x i − X i . la variación del tensor de deformaciones unitarias de Green Lagrange es: δE = 1 δ FT F + FT δ F 2 ( ) 3.1 Expresión en función del gradiente de desplazamientos. Dado que se cumple que δ F = δ H . pues como puede comprobarse. ⎡ ∂u1 ⎢ ⎢ ∂X 1 ⎢ ⎢ A=⎢ 0 ⎢ ⎢ ⎢ ∂u ⎢ 1 ⎢ ∂X ⎣ 2 3.3.
con lo que finalmente la variación del tensor de deformaciones de Green – Lagrange es: δ E = FT δε F 3. (ds ) − (ds 0 ) = dxT dx − dXT dX = dxT dx − dxT F−T F−1dx 2 2 (ds ) − (ds 0 ) = dxT (I − F−T F−1 )dx 2 2 Con lo que se define este tensor euleriano.4 Tensor de deformaciones unitarias euleriano (29) De forma análoga al tensor de Green-Lagrange se puede definir el tensor de deformación unitaria en el planteamiento euleriano e o tensor de Almansi.2 Expresión en función del tensor infinitesimal de deformaciones unitarias La variación del tensor gradiente F es δ F = δ Ht F . también se comprueba fácilmente que: A(δ H) H = A(H) δ H De esta forma el último sumando de la expresión anterior se puede poner como: 1 1 1 δA H = A(δ H) H = A(H) δ H 2 2 2 Con lo que se obtiene que: 1 1 δ E = AC δ H + A δ H + A δ H = (AC + A) δ H 2 2 3.3. Por otra parte. Para ello se considera la diferencia entre los cuadrados de las distancias entre dos puntos infinitamente próximos. como: e= 1 1 (I − F−T F−1 ) = (I − B−1 ) 2 2 15 . desarrollando las expresiones.Introducción al análisis de estructuras con no linealidad geométrica 3. pero referidas al estado deformado final. Por lo tanto: δE = 1 T 1 T T F δ HT F δ HT t F + F δ Ht F = t + δ Ht F 2 2 ( ) ( ) En esta expresión se identifica la variación del tensor infinitesimal de deformaciones unitarias δε en el instante t.3 Variación del tensor de Green-Lagrange en forma de vector La variación del tensor de Green-Lagrange en su forma de vector es: (28) 1 1 δ E = AC δ H + A δ H + δ A H 2 2 Observando el valor de A se comprueba que se cumple que δ A(H) = A(δ H) .3. en función del tensor de Finger.
5. Para ello ambos tensores emplean la diferencia entre los cuadrados de las distancias entre dos puntos próximos. mientras que los dos últimos se deben a las deformaciones iniciales ya existentes en el material en el instante t. en los estados t y t+Δt. Para su empleo en las formulaciones incrementales.1 Tensor incremental de Green – Lagrange El tensor incremental de Green – Lagrange E 0 se define como la diferencia entre los cuadrados de las distancias entre dos puntos próximos. resulta útil estudiar el incremento que sufre el tensor de Green – Lagrange al pasar desde una configuración t a otra t+Δt. en relación al incremento de deformación u ⎞ ⎛ ∂u ∂ u ⎞ ⎛ ∂u ˆj ∂u ˆi ˆk ⎟ ˆk ∂uk ⎞ 1 ⎛ ∂u ⎟ ⎟ ⎜ ⎜ k ⎟ ⎟ ⎜ ˆij = ⎜ ⎜ ⎜ e + + + ⎟ ⎟ ∑⎜ ⎟⎟ ∑ ⎜ ⎜ ⎟ ⎜ ⎟ ⎟ ⎜ ∂X i ∂ X j ⎠ ⎜ ∂X i ∂X j ⎠⎠ 2 ⎝ ∂X j ∂X i ⎟ k ⎝ k ⎝ Los dos primeros términos son similares al tensor de deformaciones infinitesimales (aunque referidos a las coordenadas iniciales). los dos tensores se expresan en función de la deformación incremental entre los dos estados: ˆ = xt +Δt − xt = ut +Δt − ut u 3. el tensor incremental se puede expresar como suma de dos tensores: ˆ =ˆ ˆ E e+η ˆ: El primer tensor contiene los términos lineales. El tensor no lineal en el incremento de deformación vale: ˆij = η ˆk ∂u ˆk 1 ∂u 2 ∂X i ∂ X j 16 . pero referidas a estados de referencia distintos: uno emplea el estado inicial como referencia y el otro emplea el estado t como referencia. A estos efectos se emplean dos tensores incrementales que miden el incremento en el tensor de Green – Lagrange entre los estados t y t+Δt.Introducción al análisis de estructuras con no linealidad geométrica 3.5 Tensores incrementales Para el desarrollo de formulaciones incrementales. expresadas en los estados indicados. referidas al estado inicial (por eso se añade el subíndice 0) ˆ dX = (ds t +Δt )2 − (ds t )2 2 dXT E con lo cual coincide con la diferencia entre los tensores en los instantes t y t+Δt: ˆ dX = (ds t +Δt )2 − (ds 0 )2 − (ds t )2 − (ds 0 )2 2 dXT E ˆ dX = 2 dXT (Et +Δt − Et )dX 2 dXT E ( ) ( ) (30) ˆ = Et +Δt − Et E Sustituyendo el valor de los tensores en función de las deformaciones.
en los estados t y t+Δt.5. referidas al estado t: ˆ dx = (ds t +Δt )2 − (ds t )2 2 dxT E t (31) Sustituyendo el valor de los diferenciales en función de las coordenadas y de los incrementos de deformación.2 Tensor incremental actualizado de Green – Lagrange El tensor incremental actualizado de Green – Lagrange Et se define como la diferencia entre los cuadrados de las distancias entre dos puntos próximos. el tensor incremental actualizado se puede expresar como suma de dos tensores: ˆ =ˆ ˆt E et + η t ˆ: El primer tensor contiene los términos lineales en relación al incremento de deformación u ˆj ⎞ ∂u ˆi 1 ⎛ ∂u ⎟ ⎟ ˆtij = ⎜ ⎜ e + ⎟ ⎜ ⎟ x x 2 ⎜∂ ∂ ⎝ j i ⎠ El tensor no lineal en el incremento de deformación vale: ˆtij = η ˆk ∂u ˆk 1 ∂u 2 ∂x i ∂ x j 17 .Introducción al análisis de estructuras con no linealidad geométrica 3.
y está asociado al movimiento de la partícula que ocupa la posición x. De hecho se puede definir la velocidad como una función de la coordenada espacial x. t ) ∂σ(x. t ) ∂x(X. t ). a pesar de que se ha expresado en las coordenadas materiales de la partícula X. la derivada temporal material requiere efectuar la derivación en cadena: σ (x.t). representa la velocidad en sentido clásico. t ) + (∇σ ) v ∂t σ (x. su derivada respecto al tiempo se define como: σ (X.2 Tensor gradiente de velocidad Habiendo definido la velocidad v. Es importante notar que la velocidad espacial es la derivada temporal de ninguna función. pues mide el cambio de σ asociado con la partícula material que está situada inicialmente en X.Introducción al análisis de estructuras con no linealidad geométrica 4 VELOCIDAD x = φ(X. 4. t ) = ∂φ(X. t ) + ∂t ∂x ∂t ∂σ(x. t ) ∂t Por su definición resulta obvio que se trata de un campo vectorial material. su derivada respecto a las coordenadas espaciales x define el tensor gradiente de velocidad L: 18 . Si el campo σ está definido en función de la posición espacial x. Físicamente representa la velocidad de la partícula que ocupa la posición X en el instante t = 0. t ) = dt ∂t Y se conoce como la derivada temporal material de la magnitud σ. v (x.1 Derivada temporal material Sea un campo escalar. t ) = ∂σ(x. de la partícula que en el instante de tiempo t ocupa la posición x. vectorial o tensorial expresado en las coordenadas materiales σ(X. t ) = v (φ−1 (x. t ) = El segundo término de esta expresión se denomina convectivo o de transporte. t ) = dσ ∂σ(X. t ) Se trata de un campo espacial que físicamente. t ) Considerando el movimiento de una partícula cualquiera definida por su posición: Se define la velocidad material de la partícula como: v (X. 4. La velocidad es un campo espacial.
Introducción al análisis de estructuras con no linealidad geométrica L (x. otra expresión del tensor gradiente de velocidad: L = F F−1 4. t ) = ∇v ∂x Lij = ∂ vi ∂x j Este tensor proporciona la velocidad relativa entre dos partículas situadas en puntos próximos P y Q. a su vez. mediante el siguiente desarrollo algebraico: 1 1 (L + LT ) = F F−1 + F−T FT 2 2 1 D = F−T FT F + FT F F−1 = F−T E F−1 2 D= ( ) ( ) • La tasa de deformación D es una medida de la variación con el tiempo del cuadrado de la longitud de un elemento diferencial: 19 .4 Tasa de deformación El tensor gradiente de velocidad se puede descomponer en sus componentes simétrica y antisimétrica: L= 1 1 (L + LT ) + (L − LT ) 2 2 ∂v j ⎞ 1 ⎛ ∂v i ⎟ ⎟ ⎜ Dij = ⎜ + ⎟ ⎟ ⎜∂ 2⎜ x x ∂ ⎝ j i ⎠ • El tensor tasa de deformación D se define como la parte simétrica de L: D= 1 (L + LT ) 2 La parte antisimétrica de L define el tensor de giro W: W= 1 (L − LT ) 2 ∂v j ⎞ 1 ⎛ ∂vi ⎟ ⎟ ⎜ Wij = ⎜ − ⎟ ⎜ ⎟ 2 ⎜∂ x x ∂ ⎝ j i ⎠ • La tasa de deformación D se relaciona con la derivada temporal del tensor de Green. t ) = ∂v(x.3 Derivada temporal del tensor de Green Se obtiene fácilmente a partir de su definición: E= 1 1 C = (FT F − FT F) 2 2 4. El tensor gradiente de velocidad permite obtener una expresión útil de la derivada temporal del gradiente de deformación: F= d ⎛ ∂x ⎞ ∂ ⎛ ∂x ⎞ ∂v ∂ v ∂x ⎟ ⎟ ⎜ ⎜ = = = = LF ⎟ ⎟ ⎜ ⎜ ⎟ ⎟ ⎝ ⎠ ⎝ ⎠ dt ∂X ∂X ∂t ∂X ∂x ∂ X Esta expresión proporciona.
Introducción al análisis de estructuras con no linealidad geométrica d d d d (ds 2 ) = (dxTdx) = dXT FT FdX = dXT C dX dt dt dt dt ( ) ( ) d (ds 2 ) = dXT C dX = 2 dXT E dX = 2 dxT F−T E F−1dx dt 1d (ds 2 ) = dxT D dx 2 dt 20 .
mediante la fórmula de Cauchy. en el estado indeformado. pero el siguiente proceso explica su naturaleza. 5. definido por su vector normal n. Sea un elemento diferencial de área en el estado deformado de módulo dA y dirección n. El vector tensión en dicha área es la fuerza de tracción actuante por unidad de área: t= df dA La fórmula de Cauchy proporciona el valor de esta fuerza de tracción por unidad de área en función del tensor de tensiones en dicho punto y de la dirección normal: t = σn 5.Introducción al análisis de estructuras con no linealidad geométrica 5 5. pero referida a la superficie sin deformar. empleando un cierto tensor de tensiones que resulta ser el segundo tensor de Piola-Kirchhoff: df 0 = S n 0dA0 21 . Según la fórmula de Cauchy esta fuerza es: df = σ n dA Supongamos que transformamos esta fuerza diferencial al estado original indeformado. Se trata de un tensor simétrico. dando lugar al primer tensor de tensiones de Piola-Kirchhoff P: df = P dA0 = P n 0 dA0 Por lo tanto este tensor proporciona la fuerza en el estado deformado.3 Segundo tensor de tensiones de Piola-Kirchhoff Este tensor no tiene mucho significado físico. definidas como las fuerzas internas por unidad de área en la situación deformada. Es un tensor no simétrico. Sea un elemento diferencial de área de módulo dA en el estado deformado. utilizando el tensor gradiente de la deformación (que puede usarse para transformar cualquier vector diferencial): df 0 = F−1 df = F−1 σ n dA Este vector de fuerzas ya transformado al estado indeformado puede expresarse.1 TENSIONES Tensor de tensiones de Cauchy Este tensor σ representa las tensiones reales existentes en el material. La fuerza actuante sobre él es df.2 Primer tensor de tensiones de Piola-Kirchhoff (32) La fuerza df actuante sobre el diferencial de área deformada dA se puede referir al área inicial no deformada dA0 . como puede deducirse de las ecuaciones de equilibrio de momentos de un cubo diferencial. sobre el cual la fuerza actuante es df.
Segundo tensor de tensiones de Piola. Este tensor corresponde a la fuerza en el estado deformado.Kirchhoff • El segundo tensor de Piola – Kirchhoff es invariante ante rotaciones de sólido rígido. sometido a un tensor de deformación Ft . Sea un sistema en un instante t. El tensor gradiente de deformación en el nuevo estado es: Ft +Δt = R Ft La tensión de Piola – Kirchhoff en el estado t+Δt es: St +Δt = F t +Δt (Ft +Δt ) −1 σt +Δt (Ft +Δt ) −T Teniendo en cuenta que el determinante del gradiente de deformaciones es el mismo al aplicar la rotación de sólido rígido: St +Δt = F t (Ft ) −1 RT σt +Δt R (Ft ) −T 22 .Introducción al análisis de estructuras con no linealidad geométrica Igualando ambas expresiones: S n 0dA0 = F−1 σ n dA La relación entre las áreas inicial y final es: n dA = F F−T n 0 dA0 Por lo tanto S n 0dA0 = F F−1 σ F−T n 0 dA0 Como esta expresión debe satisfacer para cualquier diferencial de área. Como puede comprobarse en su expresión. es un tensor simétrico. Entre t y t+Δt se aplica una rotación de sólido rígido definida por una matriz R. pero transformada al estado inicial y referida a la unidad de área del estado inicial. df n dA n0 df 0 dA0 S t=0 F-1 Figura 4. se debe cumplir que: S = F F−1 σ F−T (33) Esta expresión permite transformar las tensiones reales en el estado deformado σ en el segundo tensor de tensiones de Piola – Kirchhoff S.
Introducción al análisis de estructuras con no linealidad geométrica Durante la rotación de sólido rígido las tensiones se mantienen constantes en el sistema móvil y sufren dicha rotación de sólido rígido. 23 . Ello es debido a que se aplica la misma matriz de rotación a las tensiones de Cauchy y al gradiente de deformaciones. por lo tanto su valor es: σt +Δt = R σt RT Sustituyendo este valor se obtiene St +Δt = F t (Ft ) RT R σt RT R (Ft ) St +Δt = F t (Ft ) −1 −1 −T σt (Ft ) −T ≡ St Es decir que las tensiones de Piola-Kirchhoff son las mismas que antes de efectuar la rotación de sólido rígido.
Introducción al análisis de estructuras con no linealidad geométrica 6 6. Sean qv las fuerzas de volumen aplicadas sobre él y t las fuerzas en su superficie circundante (al ser el trozo de sólido arbitrario. en el estado deformado en el instante t. ∫ r×q v v dv + ∫ r × t ds = 0 s Las fuerzas en la superficie se pueden sustituir por las tensiones σ en la superficie empleando la fórmula de Cauchy: ∫ r×q v v dv + ∫ r × (σ n)ds = 0 s La segunda integral se puede transformar en dos integrales de volumen empleando el teorema de integración y efectuando ciertos desarrollos algebraicos: ∫ r×q v v dv + ∫ r × div(σ) + ∫ e : σT dv = 0 v v 24 . el integrando tiene que ser nulo siempre: div(σ) + q v = 0 Esta es la ecuación de equilibrio estático del sólido. parte de estas fuerzas de superficie serán fuerzas aplicadas conocidas y otras serán fuerzas interiores desconocidas). de volumen v y área lateral s. y se aplica la ecuación de equilibrio estático de momentos respecto de un punto cualquiera. puesta en función de las tensiones de Cauchy.2 Equilibrio de momentos Se considera de nuevo un trozo cualquiera de sólido. ∫q v v dv + ∫ σ n ds = 0 s La segunda integral se puede transformar en una integral de volumen empleando el teorema de integración aplicado al tensor σ ∫q v v dv + ∫ div(σ) dv = 0 v Como el trozo de sólido es arbitrario. Sea r el vector que define la posición de la fuerza respecto del punto donde se toman momentos. en su forma más compacta.1 ECUACIONES DE EQUILIBRIO Equilibrio de fuerzas Se considera un trozo cualquiera de sólido. 6. El equilibrio estático de dichas fuerzas implica que: ∫q v v dv + ∫ t ds = 0 s Las fuerzas en la superficie t se pueden sustituir por las tensiones σ en la superficie empleando la fórmula de Cauchy.
Introducción al análisis de estructuras con no linealidad geométrica El símbolo e representa el tensor alternador de orden 3. en el que existe un campo de deformaciones u y en el que se aplica una variación virtual δ u a dicho campo de deformaciones. ∫ r × (q v v + div(σ))dv + ∫ e : σT dv = 0 v La primera integral es nula pues su integrando contiene la ecuación de equilibrio. como el trozo de sólido es arbitrario. el integrando de la segunda tiene que ser nulo siempre: e : σT = 0 Desarrollando el producto contracto se obtiene un tensor de orden 1 (se contraen 2 índices): ⎧ ⎪ ⎪ σ32 − σ23 ⎫ ⎪ ⎪ ⎪ ⎪ T ⎪ e : σ = ⎨σ13 − σ31 ⎪ ⎬=0 ⎪ ⎪ ⎪ ⎪ − σ σ ⎪ ⎪ 21 12 ⎪ ⎪ ⎩ ⎭ Esta condición indica que el tensor de tensiones de Cauchy σ es simétrico. empleando la fórmula de Cauchy q S = σ n δWE = ∫ δu ⋅ q s v ds + ∫ δ u ⋅ σ n ds = s ( ) ∫ δu ⋅ q v v dv + ∫ (δ u σ) ⋅ n ds s La integral a la superficie se puede transformar en integral al volumen empleando el teorema de la divergencia: δWE = ∫ δu ⋅ q v v dv + ∫ div(δ u σ) dv v El integrando de la segunda integral se puede desarrollar utilizando la propiedad indicada en el anejo: div(δ u σ) = div(σ δ u) = δ u ⋅ div(σ) + σ : grad(δ u) Sustituyendo y agrupando las integrales se obtiene: δWE = ∫ δ u ⋅ [q v v + div(σ)]dv + ∫ σ : ∇(δ u ) dv v 25 .3 Principio del trabajo virtual Sea un cuerpo en equilibro en un estado cualquiera t.k} es par. Además.j. definido como eijk = 1 si la permutación {i. El trabajo virtual de las fuerzas exteriores aplicadas sobre el volumen y sobre la superficie es: δWE = ∫ δu ⋅ q v v dv + ∫ δ u ⋅ q S ds s (34) Las fuerzas de superficie se pueden sustituir por las tensiones de Cauchy σ en la superficie. 6. -1 si la permutación es impar y 0 si hay índices repetidos.
El tensor de tensiones de Cauchy se puede poner en función del 2º tensor de Piola – Kirchhoff despejando de la ecuación (33) : σ = F −1 F S FT La variación del tensor de deformaciones unitarias se puede poner en función de la variación del tensor de Green-Lagrange según la ecuación (28) : 26 . deformaciones infinitesimales ε como el volumen de integración se refieren al sólido en el estado deformado en el instante t. que se refieren a un estado conocido. En la segunda integral. con lo que: δWE = ∫ σ : δ(∇u ) dv v El tensor gradiente de la deformación ∇u se puede descomponer como suma de sus componentes simétrica y hemisimétrica.Introducción al análisis de estructuras con no linealidad geométrica El integrando de la primera integral es nulo pues corresponde a la ecuación de equilibrio del sólido. Para resolver este problema se emplean las magnitudes de medida de tensión y deformación anteriormente definidas. su producto contracto con la parte hemisimétrica es nulo. compatible con las condiciones de ligadura: Este principio es la herramienta fundamental para el desarrollo de una formulación de análisis no lineal. el gradiente de la variación de deformación es igual a la variación del gradiente ∇(δ u) = δ(∇u) . que es desconocido. con lo que se llega a: ∫ σ : δ(∇u) dv = ∫ σ : δε dv + ∫ σ : δΗ dv t v v v = ∫ σ : δε dv v Recuperando la expresión inicial del trabajo virtual de las fuerzas exteriores se puede poner la siguiente expresión general del principio del trabajo virtual: ∫ δu ⋅ q v v dv + ∫ δ u ⋅ q S ds = s ∫ σ : δε dv v El término de la derecha corresponde al trabajo virtual de las fuerzas interiores: δWI = ∫ σ : δε dv v (35) Por lo que de forma compacta se pude poner el principio del trabajo virtual como: δWI = δWE (36) El principio del trabajo virtual aplicado a esa configuración indica que la condición necesaria y suficiente para que exista equilibrio es que el trabajo virtual de las fuerzas interiores δWI sea igual al trabajo virtual de las fuerzas exteriores δWE para cualquier variación virtual de las deformaciones δ u . siendo la componente simétrica el tensor infinitesimal de deformaciones unitarias ε: ∇u = Ht = ε + Ht Por ser el tensor σ simétrico. pues tanto los tensores de tensiones σ. Sin embargo su aplicación directa no es fácil.
Introducción al análisis de estructuras con no linealidad geométrica δε = F−T δ E F−1 Sustituyendo estas expresiones en (35) se obtiene el siguiente valor del trabajo virtual interior: δWI = Operando se llega a: Vt ∫ F (F S F ) : (F −1 T −T δ E F−1 dv t ) δWI = Vt ∫F −1 S : δ E dv t La relación entre los diferenciales de volumen es dv t = F dv 0 Con lo que finalmente se obtiene: δWI = ∫ S : δE dv V 0 0 (37) Empleando la representación como vectores de los tensores S y E. Esta expresión es la base de las formulaciones que se desarrollan a continuación 27 . la expresión final del trabajo virtual interior es: δWI = V0 ∫ δE T S dv 0 (38) Se ha obtenido así una expresión del trabajo virtual interior en el estado t. pero en la que se emplean magnitudes de tensión (S) y de deformación unitaria (E) referidas al estado inicial conocido del cuerpo. por lo que su empleo es mucho más sencillo.
28 . Su tamaño en problemas de dos dimensiones es de 4 filas y tantas columnas como grados de libertad tiene el elemento. U 2 { } T Con lo que: H = ∂ 0 u = ∂ 0 N U = G0 U (41) La matriz G0 contiene las derivadas de las funciones de interpolación con respecto a las coordenadas iniciales y no depende de las deformaciones.. en el caso plano (el primer superíndice especifica el nudo): 1 1 2 n U = U1 U2 U 12 U 2 .1 FORMULACIÓN LAGRANGIANA TOTAL Trabajo virtual interior El equilibrio en el instante t viene dado por el principio del trabajo virtual: δWI = V0 ∫ δE T S dv 0 = δWE La variación de la deformación unitaria viene dada por (29) δ E = (AC + A) δ H Las derivadas de los desplazamientos contenidas en H se pueden expresar en función del campo de deformaciones a través de un operador de derivación ∂ 0 . a través de unas funciones de interpolación N: u = NU (40) El vector U contiene las deformaciones de los nudos en el instante t y es. nos centramos en el estudio de un solo elemento e introducimos la hipótesis de discretización del método de los elementos finitos: el campo de deformaciones se aproxima por interpolación de las deformaciones de los nudos U. que en el caso plano es: ⎡ ∂ ⎢ ⎢ ∂X 1 ⎢ ⎢ ∂ ⎢ ⎢ ∂X 2 H=⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎣ ⎤ 0 ⎥ ⎥ ⎥ ⎥ 0 ⎥ ⎥⎪ ⎧u1 ⎪ ⎫ ⎪ ⎥⎪ ⎨ ⎬ = ∂0 u u ∂ ⎥⎪ 2⎪ ⎪ ⎪ ⎩ ⎭ ⎥ ∂X 1 ⎥ ⎥ ∂ ⎥ ⎥ ∂X 2 ⎥⎦ (39) A partir de este punto.Introducción al análisis de estructuras con no linealidad geométrica 7 7..
. la ecuación de equilibrio resulta: ∫B V 0 T S dv 0 = P (48) 29 ..2 Ecuaciones de equilibrio El principio del trabajo virtual queda: δWI = δ UT Q = δWE (47) El trabajo virtual de las fuerzas exteriores se sustituye por el trabajo virtual producido por las fuerzas nodales equivalentes a las fuerzas exteriores P. ..... cuya determinación se deja para más adelante: δWE = δ UT P Al ser arbitraria la variación de los desplazamientos se cumple que: Q=P Empleando el valor detallado de las fuerzas interiores.⎥ ⎥ ⎥ ⎥ .⎥ ⎥ ⎦ (42) La variación de las deformaciones unitarias queda por lo tanto: δ E = (AC + A) G0 δ U = B δ U Donde se ha definido la matriz B = (AC + A) G0 (43) (44) que proporciona la relación entre la variación de las deformaciones de los nudos y la variación de las deformaciones unitarias de Green.... ..⎥ ⎥ ⎥ ⎥ ..... .Introducción al análisis de estructuras con no linealidad geométrica ⎡ ∂N 1 ⎢ ⎢ ∂X 1 ⎢ ⎢ ∂N 1 ⎢ ⎢ ∂X G0 = ⎢ 2 ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎣ 0 0 ∂N 1 ∂X1 ∂N 1 ∂X 2 ∂N 2 ∂X 1 ∂N 2 ∂X 2 0 0 0 0 ∂N 2 ∂X 1 ∂N 2 ∂X 2 ⎤ .. El trabajo virtual interior vale: δWI = V0 ∫ δU T BT S dv 0 = δ UT ∫ BT S dv 0 V0 (45) En esta expresión se define el vector de fuerzas nodales equivalentes a los esfuerzos interiores en el elemento.⎥ ⎥ ⎥ ⎥ . .. Q≡ ∫B V 0 T S dv 0 (46) 7.
la matriz de rigidez tangente: Δ (δWI ) = δ UT ΔQ = δ UT ∂Q ˆ ˆ U ˆ U = δ UT K ∂U El incremento en el trabajo virtual interior es. Al estar empleando un método incremental. y el incremento de dicha deformación tiene una expresión similar. Esta aproximación es válida para materiales elásticos e incluso para otros comportamientos más sofisticados. por conveniencia para desarrollos posteriores.Kirchhoff. proporcional al incremento de la deformación de los nudos del elemento U empleando para ello una matriz de rigidez constante. por lo que procedemos a linealizarlo con respecto a un incremento de ˆ. la notación de vectores para el primer término y la de tensores para el segundo. pues es no lineal. tanto en S como en E. se obtiene el valor del primer término en el incremento del trabajo virtual interior: 30 . y ser los incrementos de deformaciones pequeños. que deseamos ponerlo en ˆ. • Para evaluar la primera integral de (49). en el que buscamos obtener el equilibrio en t+Δt a partir del equilibrio conocido en t. forma lineal. Donde el vector U Sustituyendo este incremento de la deformación unitaria. que viene dado por el principio del trabajo virtual en dicho instante: δWIt +Δt = δWEt +Δt La resolución directa de esta ecuación es muy difícil.Introducción al análisis de estructuras con no linealidad geométrica 7.3 Linealización de las ecuaciones de equilibrio Suponemos conocido el equilibrio en el instante t y buscamos el equilibrio en t + Δt . El término no lineal corresponde al trabajo virtual interior. Por otra parte. la deformación de los nudos del elemento finito que denominaremos U Desarrollando en serie alrededor del punto anterior t cuyo equilibrio se conoce: δWIt +Δt ≈ δWI + Δ (δWI ) El segundo sumando es el incremento en el trabajo virtual interior. así como su variación. resulta aceptable suponer que el incremento de dicha tensión es proporcional al incremento en las deformaciones de Green-Lagrange: ΔS = C ΔE La matriz C es constante y representa la ecuación constitutiva del material en términos incrementales. Para su resolución efectuamos un planteamiento incremental. dada por la ecuación (43). según (38): Δ (δWI ) = V0 ∫ δE T ΔS dv 0 + ∫ Δ (δ E) : S dv 0 V0 (49) En esta expresión se ha mantenido. la variación en la deformación de Green-Lagrange viene dada por (43). es necesario establecer un valor del incremento en la tensión de Piola . dada la similitud entre variaciones e incrementos: ˆ ΔE = B U ˆ contiene los incrementos en las deformaciones nodales del elemento.
ambos productos contractos son iguales: I2 = ∫ (δH T ΔH : S dv 0 ) La evaluación del integrando en esta forma resulta complicada para la implementación práctica. Al ser S simétrica. su variación δH e incremento ΔH tienen la misma expresión que ˆ respectivamente. Al ser H lineal en las deformaciones (ver (9)). Para el caso de 2 H. Apoyándose en la ecuación (27). por lo que se desarrolla en función del valor de los tensores H y S.Introducción al análisis de estructuras con no linealidad geométrica V0 ∫ δE T ˆ ΔS dv 0 = δ UT ∫ BT C B dv 0 U V0 • Para calcular la segunda integral de (49) es necesario establecer una expresión del incremento de la variación de la deformación unitaria de Green . pero sustituyendo u por δu y por su incremento u dimensiones la integral queda: I2 = ∫ ⎡ ∂δu1 ⎢ ⎢ ∂X1 ⎢ ⎢ ∂δu2 ⎢ ⎢⎣ ∂X1 ∂δu1 ⎤ ⎥ ∂X 2 ⎥ ⎥ ∂δu2 ⎥ ⎥ ∂X 2 ⎥⎦ T ⎡ ∂u ˆ ⎢ 1 ⎢ ∂X 1 ⎢ ⎢ ∂u ˆ ⎢ 2 ⎢⎣ ∂X1 ˆ1 ⎤ ∂u ⎥ ∂X 2 ⎥ ⎡⎢S11 S12 ⎤⎥ 0 ⎥: dv ˆ2 ⎥ ⎢⎢S12 S 22 ⎥⎥ ∂u ⎦ ⎥ ⎣ ∂X 2 ⎥⎦ Desarrollando los productos se puede demostrar que es posible poner el integrando en la forma siguiente: I2 = ∫ ⎧ ∂δu1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ X ⎪ ∂ ⎪ 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δ u ∂ ⎪ 1⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎪ ⎪ ⎨ ⎬ ⎪ ∂δu2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ X ∂ 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δ u ∂ ⎪ ⎪ 2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎩ ⎭ T ⎡S11 S12 0 ⎢ ⎢S ⎢ 12 S 22 0 ⎢ 0 S11 ⎢0 ⎢ ⎢0 0 S12 ⎣⎢ ⎧ ˆ1 ⎫ ∂u ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ X ⎪∂ ⎪ 1 ⎪ ⎪ ⎪ 0 ⎤⎪ ⎪ ˆ u ∂ ⎥⎪ ⎪ 1⎪ ⎪ ⎪ 0 ⎥⎥ ⎪ ⎪ ⎪ ∂X 2 ⎪ 0 ⎪ ⎪ ⎥⎨ ⎬ dv ⎪ ⎪ S12 ⎥ ⎪ ∂u ˆ2 ⎪ ⎥⎪ ⎪ ⎪ ∂X 1 ⎪ ⎪ S 22 ⎥⎥ ⎪ ⎪ ⎦⎪ ⎪ ˆ2 ⎪ u ∂ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎭ ⎩ 31 .Lagrange. este incremento es: Δ(δ E) = Dado que ΔF = ΔH queda: 1 δ HT ΔF + ΔFT δ H 2 ( ) ) ) Δ(δ E) = La segunda integral es: 1 δ HT ΔH + ΔHT δ H 2 ( I2 = ∫ Δ(δE) : S dv 0 = 1 δ HT ΔH + ΔHT δ H : S dv 0 ∫ 2 ( Transponiendo el segundo sumando I2 = T 1 1 0 T T : : S dv 0 δ H Δ H S dv + δ H Δ H ∫ ∫ 2 2 ( ) ( ) En ambas integrales la matriz que multiplica a S es la misma transpuesta.
sobre un material matriz K D dado. El último factor corresponde a las derivadas de los ˆ y define un nuevo vector incrementos de las deformaciones u ⎧ ∂u ⎫ ˆ1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪∂ ⎪ X 1 ⎪ ⎪ ⎪ ⎪ ⎪ ˆ1 ⎪ ∂u ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎪ ˆ ⎪ ⎪ H=⎨ ⎬ ⎪ ⎪ ˆ ∂ u 2 ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ˆ2 ⎪ ∂u ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂ X ⎪ 2⎪ ⎪ ⎩ ⎭ Siendo los incrementos de deformación: ⎧u ⎫ ˆ1 ⎪ ⎪ ⎪ ˆ =⎪ u ⎨ ⎬ ⎪ ˆ⎪ u ⎪ ⎩ 2⎪ ⎭ Con lo que la segunda integral queda: I2 = ∫ δH T ˆ dv 0 S H En esta expresión se ha definido la matriz S . Nótese su similitud con la matriz de rigidez en el análisis lineal. que tiene dos sumandos. La primera Esta expresión define la matriz de rigidez tangente K ˆ corresponde a la rigidez asociada al incremento de las tensiones.Introducción al análisis de estructuras con no linealidad geométrica El primer factor del integrando es la variación del vector H (representación como vector del tensor gradiente de la desplazamiento H). ⎡S11 S12 0 0⎤ ⎢ ⎥ ⎢S ⎥ 0 0 S 12 22 ⎢ ⎥ S=⎢ ⎥ 0 S11 S12 ⎥ ⎢0 ⎢ ⎥ ⎢0 ⎥ 0 S S 12 22 ⎥ ⎢⎣ ⎦ (50) ˆ se puede expresar en función de los incrementos de las deformaciones de los El vector H nudos efectuando el mismo desarrollo que para el vector H (ver (41)): ˆ =∂ u ˆ ˆ H 0 ˆ = ∂ 0 N U = G0 U Finalmente la segunda integral queda: ∫ Δ (δE) : S dv 0 = (δ U) T ∫G T 0 ˆ S G0 dv 0 U • El incremento del trabajo virtual queda por lo tanto: ˆ + δ UT GT S G dv 0 U Δ (δWI ) = δ UT ∫ BT C B dv 0 U ∫ 0 0 ˆ ˆ U ˆ + δ UT K ˆ U ˆ = δ UT K ˆ U ˆ Δ (δWI ) = δ UT K σ D ˆ . que consiste en una agrupación diagonal del 2º tensor de tensiones de Piola-Kirchhoff tantas veces como dimensiones tenga el problema. 32 . aunque ahora la matriz B es dependiente de las deformaciones existentes.
No depende de las propiedades del material sino sólo del estado de tensiones (a través de S) y de la geometría (a través de G).4 ∫G T 0 S G0 dv 0 (52) Ecuaciones de equilibrio incrementales Tras la linealización del trabajo virtual interior. en forma compacta ˆ U ˆ = Pt +Δt − Q K (54) Esta expresión es la ecuación incremental de equilibrio del elemento finito. se debe cumplir: ˆU ˆ = Pt +Δt − BT S dv 0 K ∫ V 0 (53) O también. mientras que el trabajo virtual de las fuerzas exteriores se sustituye por sus fuerzas nodales equivalentes Pt+Δt. El término de la izquierda representa el incremento aproximado de fuerza interior que se obtiene al aplicar un incremento a la deformación. el trabajo virtual de las fuerzas interiores en el instante conocido t viene dado por (45). de ahí su nombre. la ecuación de equilibrio en t + Δt queda: t +Δt δWIt +Δt ≈ δWI + Δ (δWI ) = δWE Sustituyendo el incremento linealizado del trabajo virtual y reordenando: ˆ U ˆ = δW t +Δt − δW δ UT K E I En el término de la derecha. El término de la derecha representa el desequilibrio entre las fuerzas exteriores aplicadas en t+Δt y la fuerzas interiores existentes en t. con lo que queda: ˆ U ˆ = δ UT Pt +Δt − δ UT BT S dv 0 δ UT K ∫ V0 Al ser arbitraria la variación de las deformaciones. ˆ ≡ K σ 7. 33 .Introducción al análisis de estructuras con no linealidad geométrica ˆ ≡ K D ∫ B C B dv T 0 (51) ˆ se denomina habitualmente matriz de rigidez geométrica y La segunda matriz K σ corresponde a la rigidez asociada al incremento de las deformaciones unitarias actuando sobre un estado de tensiones ya existente.
Introducción al análisis de estructuras con no linealidad geométrica Pt+ t Pt ^ u Qt ut ut+ t Figura 5. Para las deformaciones en el instante t la interpolación es: ui = ∑ N k U ik k u = NU (56) Para los incrementos de deformaciones se emplea una expresión similar: ˆk ˆi = ∑ N k U u i k ˆ ˆ = NU u (57) 7. para la derivada de una función de interpolación: 34 .5. Por ejemplo. 7.2 Interpolación de deformaciones.5 Formulación isoparamétrica Asumiendo una formulación isoparamétrica para el elemento.5. Sistema linealizado. que se interpolan con respecto a las de los nudos (el superíndice k indica el nudo): X i = ∑ N k X ik k (55) 7.1 Interpolación de coordenadas En principio sólo son necesarias las coordenadas en el estado inicial.3 Transformación de derivadas Las derivadas de las distintas magnitudes involucradas se transforman entre el sistema local normalizado y el general por medio de la matriz jacobiana habitual. Se supone un sistema de coordenadas normalizadas ξi local al elemento. 7. resulta sencillo desarrollar el proceso para obtener la matriz de rigidez tangente y el vector de fuerzas interiores.5. en el que se definen las funciones de interpolación.
cada uno de los cuales contiene las derivadas de una función de interpolación respecto de las coordenadas iniciales. Gn G2 0 0⎥ ⎦ Su cálculo es inmediato empleando la transformación de coordenadas inversa: 35 . ⎡ ∂N k ⎢ ⎢ ∂X 1 ⎢ ⎢ ∂N k ⎢ ⎢ ∂X 2 k G0 = ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎣ ⎤ 0 ⎥ ⎥ ⎥ ⎥ 0 ⎥ ⎥ ⎥ ∂N k ⎥ ⎥ ∂X1 ⎥ ⎥ ∂N k ⎥ ⎥ ∂X 2 ⎥⎦ G0 = ⎡⎢G1 ⎣ 0 ⎤ .Introducción al análisis de estructuras con no linealidad geométrica ⎡ ∂N k ⎤ ⎡ ∂X1 ⎢ ⎥ ⎢ ⎢ ∂ξ1 ⎥ ⎢ ∂ξ1 ⎢ ⎥=⎢ ⎢ ∂N k ⎥ ⎢ ∂X1 ⎢ ⎥ ⎢ ⎢⎣ ∂ξ2 ⎥⎦ ⎢⎣ ∂ξ2 ∂X 2 ⎤ ⎡ ∂N k ⎤ ⎥⎢ ⎥ ∂ξ1 ⎥ ⎢ ∂X1 ⎥ ⎥⎢ ⎥ ∂X 2 ⎥ ⎢ ∂N k ⎥ ⎥⎢ ⎥ ∂ξ2 ⎥⎦ ⎢⎣ ∂X 2 ⎥⎦ (58) ∂N k ∂N k = J0 ∂ξ ∂X Los distintos términos de la jacobiana se obtienen mediante la interpolación de coordenadas y las derivadas de las funciones de interpolación en coordenadas locales: J 0ij = ∂X j ∂ξi =∑ k ∂N k k Xi ∂ξi Es necesaria asimismo la transformación de coordenadas inversa: ∂N k 1 ∂N k = J− 0 ∂X ∂ξ ⎧ ⎧ ∂N k ⎫ ∂N k ∂N k ⎫ ⎪ ⎪ ⎪ ⎪ 1 1 ⎪ ⎪ ⎪ ⎪ + ( J− ( J− 0 )11 0 )12 ⎪ ⎪ ⎪ ⎪ ξ ξ X ⎪ ∂ ⎪ ⎪ ∂ ∂ ⎪ 1 1 2 ⎪ ⎪ ⎪ ⎪ ⎨ ⎬=⎨ ⎬ ⎪ ⎪ ⎪ ⎪ N N N ∂ ∂ ∂ 1 1 − − k⎪ k k⎪ ⎪ ⎪ + ( J0 )22 J0 )21 ( ⎪ ⎪ ⎪ ⎪ ⎪ ∂ξ1 ∂ξ2 ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ ⎩ ⎭ 7..5..4 Matriz G Está formada por una serie de tantos bloques como nudos tiene el elemento.
5 Vector de gradiente de los desplazamientos H Su cálculo es asimismo inmediato.7 Matriz de rigidez tangente Su expresión tiene dos sumandos: (61) ˆ = K ∫ B C B dv T 0 0 + ∫ GT 0 S G0 dv • En la primera integral.6 Matriz A En realidad esta matriz sólo contiene los términos del tensor H.5.5. ordenados de otra manera: ⎡H 0 ⎤⎥ ⎢ 11 0 H 21 ⎢ ⎥ 0 H 22 ⎥ A = ⎢ 0 H 12 ⎢ ⎥ ⎢H 12 H 11 H 22 H 21 ⎥ ⎥⎦ ⎣⎢ 7. la matriz B se puede descomponer en dos sumandos: B = (AC + A) G0 = AC G0 + A G0 = BL 0 + BN 0 El primer sumando proviene de los términos lineales en la deformación y da lugar a la matriz constante BL 0 . Para la matriz constante: 36 .Introducción al análisis de estructuras con no linealidad geométrica ⎡ −1 ∂N k ⎢( J0 )11 ⎢ ∂ξ1 ⎢ ⎢ −1 ∂N k ⎢( J0 )21 ⎢ ∂ξ1 k G0 = ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ∂N k ∂ξ2 ∂N k 1 + ( J− 0 )22 ∂ξ2 1 + ( J− 0 )12 0 0 ∂N k (J ) ∂ξ1 ∂N k 1 ( J− 0 )21 ∂ξ1 −1 0 11 ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ 0 ⎥ ⎥ ∂N k ⎥ −1 ⎥ + ( J0 )12 ∂ξ2 ⎥ ⎥ ∂N k ⎥ −1 ⎥ + ( J0 )22 ∂ξ2 ⎥⎦ 0 (59) 7. un bloque para cada nudo. El segundo es proporcional al estado de deformaciones existente a través de la matriz A y da lugar a la matriz no lineal BN 0 . Los valores de estas matrices se obtienen fácilmente a partir de las A. Tienen una estructura de bloques. similar a la de G. y requiere conocer las deformaciones de los nudos en el estado conocido t: ⎧ ∂u1 ⎪ ⎫ ⎪ ⎧ ∂N k k ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ U1 ⎪ ∑ ⎪ ⎪ ⎪ ⎪ X1 ⎪ ⎪ ⎪∂ ⎪ k ∂X 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎧H 11 ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ∂ u1 ⎪ ∂N k k ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ U ⎪ ⎪ ⎪ ∑ 1⎪ ⎪ ⎪ ⎪ ⎪ ⎪ H 12 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎪ k ∂X 2 ⎬ H=⎪ ⎨ ⎬=⎨ ⎬=⎨ ⎪ ⎪ ⎪ ⎪ ⎪ H u N ∂ ∂ k⎪ 21 ⎪ k 2 ⎪ ⎪ ⎪ ⎪ ⎪ U2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∑ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ X X ∂ ∂ k ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 1 1 H ⎪ ⎪ ⎪ ⎪ ⎩ 22 ⎪ ⎭ ⎪ ⎪ ⎪ ⎪ ⎪ u N ∂ ∂ ⎪ ⎪ ⎪ k⎪ k 2 ⎪ ⎪ ⎪ U2 ⎪ ∑ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂X 2 ⎪ ⎪ ⎪ ⎪ k ∂X 2 ⎪ ⎪ ⎩ ⎭ ⎩ ⎭ (60) 7. AC y G.5.
Introducción al análisis de estructuras con no linealidad geométrica BL 0 = AC G0 = ⎡⎢ B(1) ⎣ L0 n) ⎤ ... en base a la matriz B y al vector de tensiones de Piola-Kirchhoff S en el estado conocido t: Q= V0 ∫G T 0 (AC + A) S dv 0 T (63) 37 ... 7. B( B(2) N0 N0⎥ ⎦ k) B( N0 ⎡ ⎢ ⎢ ⎢ ⎢ k = A G0 = ⎢ ⎢ ⎢ ⎢ ∂N ⎢ k ⎢ ∂X 2 ⎣ ˆ = K D ˆ = K D ∂N k H 11 ∂X1 ∂N k H 12 ∂X 2 ∂N k H 11 + H 12 ∂X 1 ⎤ ∂N k ⎥ H 21 ⎥ ∂X 1 ⎥ ⎥ ∂N k ⎥ H 22 ⎥ ∂X 2 ⎥ ⎥ ∂N k ∂N k H 21 + H 22 ⎥⎥ ∂X 2 ∂X 1 ⎦ Sustituyendo en la primera integral se obtiene el primer sumando de la matriz tangente: ∫ (B L0 + BN 0 ) C (BL 0 + BN 0 )dv 0 T (62) ∫B +∫ B T L0 0 C BL 0 dv 0 + ∫ BT L 0 C BN 0 dv 0 C BL 0 dv 0 + ∫ BT N 0 C BN 0 dv T N0 ˆ El primer sumando corresponde a la matriz de rigidez lineal K D 0 y los 3 restantes a la componente no lineal. ˆ =K ˆ +K ˆ +K ˆT + K ˆ K D D0 D1 D1 D2 • El segundo sumando de la matriz tangente corresponde a la matriz de rigidez geométrica y se puede evaluar directamente..8 Vector de fuerzas interiores Su expresión general es: Q= V0 ∫B T 0 S dv 0 Es sencillo de evaluar.5. B( B(2) L0 L0 ⎥ ⎦ k) B( L0 ⎡ ∂N k ⎢ ⎢ ∂X 1 ⎢ ⎢ ⎢ 0 = AC Gk = 0 ⎢ ⎢ ⎢ ∂N ⎢ k ⎢ ∂X 2 ⎣ ⎤ 0 ⎥ ⎥ ⎥ ∂N k ⎥⎥ ∂X 2 ⎥⎥ ∂N k ⎥⎥ ∂X1 ⎥⎦ Para la matriz no lineal: BN 0 = A G0 = ⎡⎢B(1) ⎣ N0 n) ⎤ ..
).e. En caso contrario. Estas expresiones son válidas si las fuerzas no dependen de la deformación. empleando la fórmula de Nanson. Se obtiene la siguiente expresión: Pt +Δt = ∫N v0 T t t 0 q tv+Δ dv 0 + ∫ NT q tS+Δ 0 0 ds s0 (64) t t En esta expresión q tv+Δ y q tS+Δ son los valores de las fuerzas de volumen y superficie en el 0 0 instante t+Δt. peso propio. las variaciones de deformación se pueden poner en función de las variaciones de deformación nodales: δWEt +Δt = δ UT ∫ NT qtv+Δt dv t + δ UT ∫ NT q tS+Δt ds t = δ UT Pt +Δt vt st Por lo tanto las fuerzas nodales equivalentes son: vt +Δt ∫ NT qtv+Δt dv t +Δt + s t +Δt ∫ NT qtS+Δt ds t +Δt = Pt +Δt La evaluación de estas fuerzas no es posible pues no se conoce ni el volumen ni la superficie en t+Δt. concentradas.Introducción al análisis de estructuras con no linealidad geométrica 7. Para poderlas evaluar se transforman al estado inicial. El proceso requiere sustituir el diferencial de volumen empleando para ello el determinante del tensor F y el diferencial de área. etc. 38 . pero referidas al volumen y superficie iniciales. como es el caso de muchas fuerzas habitualmente (p. la presencia de fuerzas dependientes de la deformación origina nuevos términos en las ecuaciones de equilibrio que no han sido tenidos en cuenta.6 Fuerzas nodales equivalentes a las fuerzas exteriores Las fuerzas nodales P equivalentes a las fuerzas exteriores producen el mismo trabajo virtual que ellas. En el instante t+Δt su valor es: δWEt +Δt = ∫ δu vt T qtv+Δt dv t + ∫ δ uT q tS+Δt ds t = δ UT Pt +Δt st Introduciendo la interpolación de deformaciones.
Introducción al análisis de estructuras con no linealidad geométrica 8 8. las derivadas de los desplazamientos contenidas en Ht se pueden expresar en función del campo de deformaciones a través de un operador de derivación ∂ t . 39 . que empleando las magnitudes en forma de vectores. Su tamaño en problemas de dos dimensiones es de 4 filas y tantas columnas como grados de libertad tiene el elemento. puede escribirse: δWI = Vt ∫ δε T σ dv t = δWE La variación de la deformación unitaria se relaciona con la variación de las derivadas de los desplazamientos por la ecuación (18): δ ε = AC δ Ht A su vez.1 FORMULACIÓN LAGRANGIANA ACTUALIZADA Trabajo virtual El equilibrio en el instante t viene dado por el principio del trabajo virtual. que en el caso plano es: ⎡ ∂ ⎢ ⎢ ∂x 1 ⎢ ⎢ ∂ ⎢ ⎢ ∂x 2 Ht = ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎣ ⎤ 0 ⎥ ⎥ ⎥ ⎥ 0 ⎥ ⎥⎪ ⎧u1 ⎪ ⎫ ⎪=∂ u ⎥⎪ ⎨ t u2 ⎬ ∂ ⎥⎪ ⎪ ⎪ ⎪ ⎭ ⎥⎩ ∂x 1 ⎥ ⎥ ∂ ⎥ ⎥ ∂x 2 ⎥⎦ Se introduce la interpolación de deformaciones apoyándose en las deformaciones de los nudos U y en las funciones de interpolación N: u = NU Con lo que: Ht = ∂ t u = ∂ t N U = Gt U (65) La matriz Gt contiene las derivadas de las funciones de interpolación con respecto a las coordenadas deformadas x y no depende de las deformaciones.
2 Ecuación de equilibrio El trabajo virtual interior vale: δWIt = ∫ δU V t T t T T t T BT t σ dv = δ U ∫ Bt σ dv = δ U Q V t (67) En esta expresión se ha definido el vector de fuerzas nodales equivalentes a los esfuerzos interiores en el instante t... . .⎥ ⎥ ⎦ La variación de las deformaciones unitarias queda por lo tanto: δ ε = AC Gt δ U = Bt δ U (66) Donde se ha definido la matriz que relaciona la variación de las deformaciones de los nudos y la variación de las deformaciones unitarias: Bt = AC Gt 8.⎥ ⎥ ⎥ ⎥ . que viene dado por el principio del trabajo virtual en dicho instante: δWIt +Δt = δWEt +Δt De la misma manera que en la formulación anterior.. Q≡ Vt ∫B T t σ dv t (68) Finalmente.. la ecuación de equilibrio resulta: Q≡ Vt ∫B T t σ dv t = P 8. expresándolo interior en el punto t. al ser arbitraria la variación de los desplazamientos.. en la dirección de un incremento de la deformación u t ˆ mediante una matriz de rigidez tangente K : ˆ U ˆ δWIt +Δt ≈ δWI + Δ (δWI ) = δWI + δ UT K El trabajo virtual es....3 Linealización de las ecuaciones de equilibrio Suponemos conocido el equilibrio en el instante t y buscamos el equilibrio en t + Δt ..Introducción al análisis de estructuras con no linealidad geométrica ⎡ ∂N 1 ⎢ ⎢ ∂x 1 ⎢ ⎢ ∂N 1 ⎢ ⎢ ∂x Gt = ⎢ 2 ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎣ 0 0 ∂N 1 ∂x 1 ∂N 1 ∂x 2 ∂N 2 ∂x 1 ∂N 2 ∂x 2 0 0 0 0 ∂N 2 ∂x 1 ∂N 2 ∂x 2 ⎤ .⎥ ⎥ ⎥ ⎥ ..⎥ ⎥ ⎥ ⎥ ..... .. procedemos a linealizar el trabajo virtual ˆ . según (35): 40 . ..
es necesario establecer un valor del incremento en la tensión. y su variación se obtiene el valor del primer término en el incremento del trabajo virtual interior: Vt ∫ δε T t ˆ Δσ dv t = δ UT ∫ BT t C Bt dv U Vt • Para la segunda integral se emplea el incremento de la variación de la deformación unitaria: Δ(δε) = 1 T δ HT t ΔHt + ΔH t δ Ht 2 ( ) Sustituyendo en la segunda integral. Donde el vector U Sustituyendo el incremento de la deformación. y el incremento de dicha deformación tiene una expresión similar. Al estar empleando un método incremental. se obtiene: I2 = ∫ (δH T t ΔHt : σ dv t ) La evaluación del integrando en esta forma resulta complicado para la implementación práctica. La variación en la deformación unitaria viene dada por (66). trasponiendo el segundo sumando y considerando que el tensor de tensiones σ es simétrico. y ser los incrementos de deformaciones pequeños.Introducción al análisis de estructuras con no linealidad geométrica δWI = Vt ∫ σ : δε dv T t = Vt ∫ δε T σ dv t Su incremento es: Δ (δWI ) = Vt ∫ δε Δσ dv t + ∫ Δ (δε) : σ dv t Vt • Para evaluar la primera integral. resulta aceptable suponer que el incremento de la tensión es proporcional al incremento en las deformaciones unitarias: Δσ = C Δ ε La matriz C es constante y representa la ecuación constitutiva del material en términos incrementales. Esta aproximación es válida para materiales elásticos e incluso para otros comportamientos más sofisticados. Para el caso 2D: I2 = ∫ ⎡ ∂δu1 ⎢ ⎢ ∂x 1 ⎢ ⎢ ∂δu2 ⎢ ⎢⎣ ∂x 1 ∂δu1 ⎤ ⎥ ∂x 2 ⎥ ⎥ ∂δu2 ⎥ ⎥ ∂x 2 ⎥⎦ T ⎡ ∂u ˆ ⎢ 1 ⎢ ∂x 1 ⎢ ⎢ ∂u ˆ ⎢ 2 ⎣⎢ ∂x 1 ˆ1 ⎤ ∂u ⎥ ∂x 2 ⎥ ⎡σ11 ⎥:⎢ ˆ2 ⎥ ⎢⎣σ12 ∂u ⎥ ∂x 2 ⎦⎥ σ12 ⎤ ⎥ t σ22 ⎥ dv ⎦ Desarrollando los productos se puede demostrar que es posible poner el integrando en la forma siguiente: 41 . por lo que se desarrolla en función del valor de los tensores. dada la similitud entre variaciones e incrementos: ˆ Δ ε = Bt U ˆ contiene los incrementos en las deformaciones nodales del elemento.
puesta en forma de vector. ˆ se puede Efectuando el mismo desarrollo que para el vector Ht . 42 . que tiene dos sumandos. el valor del vector H t expresar en función de los incrementos de las deformaciones de los nudos (ver (65)): ˆ =∂ u ˆ ˆ H t t ˆ = ∂ t N U = Gt U Con lo que la segunda integral queda: ∫ Δ (δε) : σ t ˆ dv t = δ UT ∫ GT t σ Gt dv U • El incremento del trabajo virtual queda finalmente como: t ˆ T T t ˆ Δ (δWI ) = δ UT ∫ BT t C Bt dv U + δ U ∫ Gt σ Gt dv U ˆ +K ˆ )U ˆ = δ UT K ˆ U ˆ Δ (δWI ) = δ UT (K D σ (69) ˆ . pues ahora la matriz primera matriz K D Bt es constante y coincide con la matriz empleada en dicho análisis lineal. La En esta expresión se ha definido la matriz de rigidez tangente K ˆ coincide con la matriz de rigidez en el análisis lineal. El último factor del integrando corresponde al gradiente de los incrementos de las deformaciones: ⎧ ∂u ⎫ ˆ1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x1 ⎪ ⎪∂ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ˆ u ∂ ⎪ 1⎪ ⎪ ⎪ ⎪ ⎪ ∂ x ⎪ ⎪ 2 ˆ =⎪ ⎪ H ⎨ ⎬ t ⎪ ⎪ ˆ ∂ u 2⎪ ⎪ ⎪ ⎪ ⎪ ∂x 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ˆ2 ⎪ ∂ u ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂ x ⎪ 2⎪ ⎪ ⎩ ⎭ Con lo que: I2 = ∫ δH T t ˆ dv t σ H t Donde se ha definido la matriz σ que consiste en una agrupación diagonal del tensor de tensiones tantas veces como dimensiones tenga el problema.Introducción al análisis de estructuras con no linealidad geométrica I2 = ∫ ⎧ ∂δu1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ⎪ ∂ ⎪ 1 ⎪ ⎪ ⎪ ⎪ ⎪ ∂δu1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ∂ ⎪ ⎪ 2 ⎪ ⎪ ⎨ ⎬ ⎪ ⎪ δ u ∂ 2⎪ ⎪ ⎪ ⎪ ⎪ ∂x 1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δ u ∂ ⎪ ⎪ 2 ⎪ ⎪ ⎪ ⎪ ⎪ x ∂ ⎪ 2⎪ ⎪ ⎩ ⎭ T ⎡σ11 ⎢ ⎢σ12 ⎢ ⎢0 ⎢ ⎢ ⎢0 ⎣ σ12 σ22 0 0 0 0 σ11 σ12 ⎧ ˆ1 ⎫ ∂u ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ⎪∂ ⎪ 1⎪ ⎪ ⎪ ⎪ 0 ⎤⎪ ⎪ ˆ1 ⎪ ∂u ⎥⎪ ⎪ ⎪ ⎥ 0 ⎥⎪ ∂x 2 ⎪ ⎪ ⎪ t ⎪ ⎪ ⎨ ⎬dv ⎥ σ12 ⎥ ⎪ ˆ2 ⎪ ∂u ⎪ ⎪ ⎪ ⎥⎪ ⎪ ⎪ x ∂ σ22 ⎥ ⎪ 1⎪ ⎪ ⎦⎪ ⎪ ˆ2 ⎪ u ∂ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ∂ ⎪ ⎪ 2⎪ ⎭ ⎩ El primer factor de esta expresión es la variación del tensor Ht.
introduciendo el vector de fuerzas interiores: ˆ U ˆ = Pt +Δt − Q K 8.5. ˆU ˆ = Pt +Δt − BT σ dv t K ∫ t Vt (72) O también. Por ejemplo. ˆ ≡ K σ 8.Introducción al análisis de estructuras con no linealidad geométrica ˆ ≡ K D ∫B T t C Bt dv t (70) ˆ tiene una estructura muy similar a la correspondiente en La matriz de rigidez geométrica K σ al formulación Lagrangiana total. que tienen la misma expresión general que en aquel caso. a base de añadirles las deformaciones obtenidas en cada paso de carga. aunque empleando la tensión de Cauchy en lugar de la de Piola Kirchhoff. y sólo se diferencia en los valores de la matriz de rigidez tangente y del vector de fuerzas interiores.3 Transformación de derivadas Las derivadas de las distintas magnitudes involucradas se transforman entre el sistema local normalizado y el general por medio de la matriz jacobiana habitual. en forma compacta.2 Interpolación de deformaciones Para las deformaciones en el instante t la interpolación es: ui = ∑ N k U ik k u = NU Para los incrementos de deformaciones se emplea una expresión similar: ˆk ˆi = ∑ N k U u i k ˆ ˆ = NU u 8.1 Interpolación de coordenadas En principio sólo son necesarias las coordenadas en el estado t: x it = ∑ N k x ikt k Estas coordenadas se deben ir actualizando a medida que progresa el análisis incremental.5 Formulación isoparamétrica 8. y evaluando todos sus términos en el estado t en lugar de en el estado inicial. 8.5. para la derivada de una función de interpolación: 43 .4 ∫G T t σ Gt dv t (71) Ecuaciones de equilibrio incrementales Efectuando el mismo desarrollo que en la formulación lagrangiana total se llega a las ecuaciones de equilibrio incrementales.5.
5. cada uno de los cuales contiene las derivadas de la función de interpolación de ese nudo respecto de las coordenadas en el instante t.4 Matriz G Está formada por una serie de tantos bloques como nudos tiene el elemento. ⎡ ∂N k ⎢ ⎢ ∂x 1 ⎢ ⎢ ∂N k ⎢ ⎢ ∂x 2 k Gt = ⎢ ⎢ ⎢ 0 ⎢ ⎢ ⎢ ⎢ 0 ⎢ ⎣ ⎤ 0 ⎥ ⎥ ⎥ ⎥ 0 ⎥ ⎥ ⎥ ∂N k ⎥ ⎥ ∂x 1 ⎥ ⎥ ∂N k ⎥ ⎥ ∂x 2 ⎥⎦ Gt = ⎡⎢G1 Gt2 … Gtn ⎤⎥ ⎣ t ⎦ Su cálculo es inmediato empleando la transformación de coordenadas inversa: 44 .Introducción al análisis de estructuras con no linealidad geométrica ⎡ ∂N k ⎤ ⎡ ∂x 1 ⎢ ⎥ ⎢ ⎢ ∂ξ1 ⎥ ⎢ ∂ξ1 ⎢ ⎥=⎢ ⎢ ∂N k ⎥ ⎢ ∂x 1 ⎢ ⎥ ⎢ ⎢⎣ ∂ξ2 ⎥⎦ ⎢⎣ ∂ξ2 ∂x 2 ⎤ ⎡ ∂N k ⎤ ⎥⎢ ⎥ ∂ξ1 ⎥ ⎢ ∂x 1 ⎥ ⎥⎢ ⎥ ∂x 2 ⎥ ⎢ ∂N k ⎥ ⎥⎢ ⎥ ∂ξ2 ⎥⎦ ⎢⎣ ∂x 2 ⎥⎦ ∂N k ∂N k = Jt ∂ξ ∂x Los distintos términos de la jacobiana se obtienen mediante la interpolación de coordenadas y las derivadas de las funciones de interpolación en coordenadas locales: J tij = ∂x j ∂ξi =∑ k ∂N k k xi ∂ξi Es necesaria asimismo la transformación de coordenadas inversa: ∂N k ∂N k = Jt−1 ∂x ∂ξ ⎧ ⎧ ∂N k ⎫ ∂N ∂N ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ( Jt−1 )11 k + (Jt−1 )12 k ⎪ ⎪ ⎪ ⎪ ⎪ ξ ξ x ⎪ ∂ ⎪ ⎪ ∂ ∂ ⎪ 1 1 2 ⎪ ⎪ ⎪ ⎪ ⎨ ⎬=⎨ ⎬ ⎪ ⎪ ⎪ ⎪ N N N ∂ ∂ ∂ − − 1 1 k⎪ k k⎪ ⎪ ⎪ + ( Jt )22 Jt )21 ( ⎪ ⎪ ⎪ ⎪ ⎪ ∂ξ1 ∂ξ2 ⎪ ⎪ ∂x 2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ ⎩ ⎭ 8.
5.Introducción al análisis de estructuras con no linealidad geométrica ⎡ −1 ∂N k ⎢( Jt )11 ⎢ ∂ξ1 ⎢ ⎢ −1 ∂N k ⎢( Jt )21 ⎢ ∂ξ1 k Gt = ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ ∂N k ∂ξ2 ∂N k + ( Jt−1 )22 ∂ξ2 + ( Jt−1 )12 0 0 ∂N k (J ) ∂ξ1 ∂N ( Jt−1 )21 k ∂ξ1 −1 t 11 ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ 0 ⎥ ⎥ ∂N k ⎥ −1 ⎥ + ( Jt )12 ∂ξ2 ⎥ ⎥ ∂N k ⎥ −1 ⎥ + ( Jt )22 ∂ξ2 ⎥⎦ 0 (73) 8.. las fuerzas nodales equivalentes a las fuerzas exteriores se evalúan con los valores en el instante t+Δt pero referidas al área y volumen de la posición conocida t.5.5 Matriz tangente Su expresión tiene dos sumandos: ˆ = K ∫B T t t C Bt dv t + ∫ GT t σ Gt dv La matriz B es constante y tiene una estructura de bloques similar a la de G..6 Fuerzas nodales equivalentes a las fuerzas exteriores Al igual que en la formulación total.6 Vector de fuerzas interiores ⎤ 0 ⎥ ⎥ ⎥ ∂N k ⎥⎥ ∂x 2 ⎥⎥ ∂N k ⎥⎥ ∂x 1 ⎥⎦ Su expresión es sencilla de obtener. en base a la matriz B y al vector de tensiones σ en el estado conocido t: Q≡ ∫B V t T t σ dv t = ∫G V t T t t AT C σ dv (74) 8. con un bloque para cada nudo: Bt = AC Gt = ⎡⎢ B1 ⎣ t Bt2 . Btn ⎤⎥ ⎦ ⎡ ∂N k ⎢ ⎢ ∂x 1 ⎢ ⎢ k k Bt = AC Gt = ⎢ 0 ⎢ ⎢ ⎢ ∂N ⎢ k ⎢ ∂x ⎣ 2 8. Pt +Δt = ∫N vt T +Δt +Δt qtvt dv t + ∫ NT q tSt ds t st 45 .
9. Formulación lagrangiana.Introducción al análisis de estructuras con no linealidad geométrica 9 ELEMENTO BIARTICULADO. Elemento de celosía plana. FORMULACIÓN LAGRANGIANA TOTAL Consideramos un elemento biarticulado plano. definido en su posición inicial mediante las coordenadas de sus nudos extremos 1 y 2: Xe = X1 Y1 X 2 Y2 Las deformaciones de los nudos en el instante t son: { } T U = U 1 V1 U 2 U 2 { } T En su posición deformada en el instante t. las coordenadas de los nudos son: xe = x 1 y1 x 2 y2 xe = Xe + U { } T Figura 6.1 Deformación unitaria La longitud inicial del elemento es: 2 L2 Y2 −Y1 )2 = XT 0 = (X 2 − X 1 ) + ( 21 X21 Siendo el vector de diferencias entre las coordenadas: X21 = (X 2 − X1 ) La longitud final deformada del elemento es: { (Y2 −Y1 ) } T 46 .
se puede poner en función de las deformaciones de los nudos: E11 = BL U + Siendo BL una matriz constante que depende de las coordenadas iniciales de los nudos y A una matriz constante: BL = 1 L2 0 ⎡−X Y21 X 21 Y21 ⎤⎥ ⎢⎣ 21 − ⎦ ⎡+1 0 −1 0⎤ ⎢ ⎥ ⎢ 0 +1 ⎥ 0 1 − ⎢ ⎥ A=⎢ ⎥ 0 +1 0⎥ ⎢ −1 ⎢ ⎥ ⎢ 0 −1 0 +1⎥ ⎣ ⎦  La variación de la deformación unitaria de Green-Lagrange es: δE11 = BL δ U + 2 T U A δU 2L2 0 δE11 = BL δ U + BN δ U = B δ U Siendo BN una matriz proporcional a las deformaciones de los nudos: BN = 1 T 1 V1 −V2 ) (U 2 −U 1 ) ( V2 −V1 )⎤⎥ U A = 2 ⎡⎢(U 1 −U 2 ) ( 2 ⎣ ⎦ L0 L0 BN = 1 L2 0 ⎡− V21 U 21 V21 ⎤⎥ ⎢⎣ U 21 − ⎦ Sumando las expresiones de la parte lineal y no lineal.Lagrange es un escalar: E11 = 2 L2 (X21 + U21 )T (X21 + U21 ) − XT t − L0 21X21 = 2 2 2L0 2L0 E11 = ⎞ 1⎛ T 1 ⎟ ⎜X U + UT 21U21 ⎟ 2 ⎜ 21 21 ⎠ L0 ⎝ 2 1 T U AU 2L2 0 Desarrollando. los vectores de diferencias entre las coordenadas finales y entre las deformaciones: x21 = (x 2 − x 1 ) { (y2 − y1 ) } T = X21 + U21 U21 = (U 2 −U 1 )  { (V2 −V1 ) } T El tensor de deformación unitaria de Green . la matriz B resulta ser: B = BL + BN = 1 L2 0 ⎡(−X −U ) (− Y21 −V21 ) (X 21 + U 21 ) (Y21 + V21 )⎤⎥ 21 21 ⎢⎣ ⎦ 47 .Introducción al análisis de estructuras con no linealidad geométrica 2 2 T T L2 t = (x 2 − x 1 ) + (y 2 − y1 ) = x21 x21 = (X21 + U21 ) (X21 + U21 ) Siendo. análogamente.
pero evaluada en la situación deformada. Como B es constante dentro del elemento la integración es inmediata: ⎧ ⎫ ⎪−(X 21 + U 21 )⎪ ⎪ −x 21 ⎫ ⎧ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ y − −(Y21 + V21 )⎪ S A0 ⎪ 21 ⎪ ⎪ S A L ⎪ 0 0 ⎪ ⎪ ⎪ ⎪ Q = S A0 L0 BT = ⎨ ⎬= ⎨ x ⎬ 2 ⎪ ⎪ ⎪ ⎪ L0 ⎪ X 21 + U 21 ⎪ L0 ⎪ 21 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ y21 ⎪ ⎪ ⎪ ⎪ ⎪ Y V + ⎩ ⎭ 21 21 ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ 9. se emplea su definición como derivada del vector de fuerzas interiores: ˆ = ∂Q K ∂U Derivando en la expresión de las fuerzas interiores: ∂Q ∂BT T ∂S ˆ ˆ +K ˆ = A0 L0 B + SA0 L0 =K K= D σ ∂U ∂U ∂U  La primera matriz se debe a la variación de la tensión S al deformarse la barra. 9.3 Matriz de rigidez tangente (75) Para su obtención. Por lo tanto: ∂(E E11 ) ∂E11 ˆ = A L BT ∂S = A L BT K = EA0 L0 BT 0 0 0 0 D ∂U ∂U ∂U La derivada de la deformación unitaria de Green-Lagrange es: ∂E11 ∂U Luego: = BL + 2 T U A = BL + BN = B 2L2 0 ˆ = EA L BT B K D 0 0 (76) Desarrollando se obtiene: 48 . Se supone que dicha tensión de Piola es proporcional a la deformación unitaria de Green siendo la contante de proporcionalidad el módulo de elasticidad del material E.Introducción al análisis de estructuras con no linealidad geométrica B= 1 L2 0 ⎡−x ⎤ ⎢⎣ 21 −y21 x 21 y21 ⎥⎦ Que corresponde a la misma expresión que la componente lineal. en lugar de aplicar las expresiones detalladas ya obtenidas.2 Vector de fuerzas interiores Su expresión es: Q= V0 ∫BS T dv 0 Siendo S la tensión de Piola – Kirchhoff en la barra.
Si ambas longitudes son muy similares (lo cual puede suponerse siempre que los alargamientos sean pequeños). salvo por el empleo del área y longitud iniciales.Introducción al análisis de estructuras con no linealidad geométrica 2 2 ⎡x 21 x 21y21 −x 21 −x 21y21 ⎤ ⎢ ⎥ ⎢ 2 2 ⎥ y21 −x 21y21 −y21 ⎥ ⎢ ˆ = EA0 ⎢ ⎥ K D 2 ⎢ ⎥ L3 Sim x x y 0 21 21 21 ⎢ ⎥ ⎢ ⎥ 2 ⎢ y21 ⎥ ⎣ ⎦ Esta matriz coincide con la matriz de rigidez convencional de una barra biarticulada. Considerando que el seno y coseno del ángulo final de la barra θ valen sθ ≡ sin θ = y21 / L y cθ = cos θ = x 21 / L . ambas matrices coinciden. pero evaluada en la posición deformada. el primero de los cuales corresponde a la matriz de rigidez lineal de la barra en su posición inicial: 2 2 ⎡X 21 X 21Y21 −X 21 −X 21Y21 ⎤ ⎢ ⎥ ⎢ 2 2 ⎥ Y21 Y21 ⎥ −X 21Y21 − EA0 ⎢⎢ ⎥ = EA0L0 BT L BL = 2 L0 ⎢⎢ sim X 21 X 21Y21 ⎥⎥ ⎢ ⎥ 2 ⎢ ⎥ Y21 ⎣ ⎦ ˆ K D0 Los sumando 2 y 3 son traspuestos uno del otro y corresponden a los términos lineales en las deformaciones: ˆ K D1 ⎡ X 21U 21 X 21V21 −X 21U 21 −X 21V21 ⎤ ⎢ ⎥ ⎢YU ⎥ Y V Y U Y V − − 21 21 21 21 21 21 21 21 EA0 ⎢ ⎥ T = EA0L0 BL BN = ⎢ ⎥ X 21V21 ⎥ L0 ⎢−X 21U 21 −X 21V21 X 21U 21 ⎢ ⎥ ⎢− ⎥ − Y U Y V Y U Y V 21 21 21 21 21 21 ⎥ ⎢⎣ 21 21 ⎦ El último sumando es cuadrático en la deformación: 49 . la matriz anterior se puede poner: 2 2 ⎡cθ cθsθ −cθ ⎢ 2 3 ⎢ sθ −cθsθ ⎢ ⎛ ⎞ EA L ⎟ 0 ⎜ ˆ = ⎢ ⎟ K ⎜ D ⎟ 2 ⎟ ⎢ L ⎜ Sim cθ ⎝ L0 ⎠ ⎢ ⎢ ⎢ ⎣ −cθsθ ⎤ ⎥ 2 ⎥ −sθ ⎥ ⎥ cθsθ ⎥⎥ ⎥ 2 ⎥ sθ ⎦ (77) Que coincide con la matriz convencional de la barra biarticulada evaluada en su posición deformada. En la expresión (76) puede desarrollarse el valor de B: ˆ = EA L (BT + BT )(B + B ) = EA L (BT B + BT B + BT B + BT B ) K D L N L N L L L N N L N N 0 0 0 0 Y se obtienen cuatro sumandos. pero introduciendo el factor de proporción entre las longitudes.
Ello permite obtener expresiones más generales de las propiedades del elemento. empleando la coordenada normalizada ξ que varía entre -1 en el nudo inicial y +1 en el nudo final: ⎡1 − ξ ⎢ ⎢ N=⎢ 2 ⎢ 0 ⎢ ⎣ 0 1− ξ 2 1+ξ 2 0 ⎤ 0 ⎥ ⎥ 1 + ξ ⎥⎥ 2 ⎥⎦ La deformación unitaria de Green-Lagrange vale: E11 = dxT dx − dXT dX 2 dXT dX Los diferenciales son: dX = dX dξ dξ 50 . son las habituales para la interpolación lineal. 9. que pueden emplearse para elementos más complejos. se puede formular el elemento de celosía empleando funciones de interpolación y la formulación isoparamétrica estándar en el método de los elementos finitos. y su valor es: T T ˆ = SA L ∂B = SA L ∂BN = SA L 1 AT = SA0 AT K 0 0 0 0 0 0 2 σ ∂U ∂U L0 L0 ⎡+1 0 −1 0⎤ ⎢ ⎥ ⎢ 0 +1 ⎥ 0 1 − SA0 ⎢ ⎥ ˆ Kσ = ⎢ ⎥ 0 +1 0⎥ L0 ⎢−1 ⎢ ⎥ ⎢ 0 −1 0 +1⎥ ⎣ ⎦ (78) Obsérvese que esta matriz es independiente de la orientación de la barra y sólo depende de su nivel de esfuerzo y de su longitud.4 Formulación isoparamétrica Aunque no es necesario.Introducción al análisis de estructuras con no linealidad geométrica 2 2 ⎡U 21 U 21V21 U 21 U 21V21 ⎤ − − ⎢ ⎥ ⎢ 2 2 ⎥ V21 U 21V21 V21 ⎥ − − EA0 ⎢⎢ ⎥ = EA0L0 BT = B N N 2 L0 ⎢⎢ sim U 21 U 21V21 ⎥⎥ ⎢ ⎥ 2 ⎢ ⎥ V21 ⎣ ⎦ ˆ K D2  La segunda matriz es la matriz de rigidez geométrica. Las interpolaciones de coordenadas y desplazamientos son: X = N Xe u = NU Las funciones de interpolación para el elemento de dos nudos.
para este caso particular.Introducción al análisis de estructuras con no linealidad geométrica ⎛ dX d u ⎞ dx = d (X + u) = ⎜ + ⎟ ⎟d ξ ⎜ ⎜ ⎟ ⎝ dξ dξ ⎠ Con lo que el valor de la deformación de Green-Lagrange es: E11 = 1 ⎛ dXT du 1 duT du ⎞ ⎟ ⎜ + ⎟ 2 ⎜ ⎟ ⎜ λ ⎝ dξ dξ 2 dξ dξ ⎠ La constante λ2 representa la expresión del denominador de la deformación unitaria: λ2 = dXT dX dXT dX = dξ dξ dξ dξ Sustituyendo las interpolaciones de deformaciones y coordenadas se obtiene: E11 = T 1 ⎛ 1 T dNT dN ⎞ eT dN dN ⎜ X U + U U⎟ ⎟ 2 ⎜ ⎟ ⎜ λ ⎝ 2 dξ dξ dξ dξ ⎠ λ 2 = XeT La variación de la deformación unitaria es: dNT dN e X dξ dξ δE11 = T T ⎞ 1⎛ ⎟ ⎜XeT dN dN δ U + UT dN dN δ U⎟ ⎜ 2 ⎜ ⎟ λ ⎝ dξ dξ dξ dξ ⎠ Esta expresión define las dos matrices B: δE11 = BL δ U + BN δ U = B δ U BL = 1 eT dNT dN X dξ dξ λ2 1 dNT dN U dξ dξ λ2 BN = Para el elemento de dos nudos sus expresiones son: ⎡ 1 1 ⎤ 0 + 0 ⎥ ⎢− dN ⎢ 2 2 ⎥ =⎢ 1 ⎥⎥ dξ ⎢ 0 −1 0 + ⎥ ⎢ 2 2⎦ ⎣ λ 2 = XeT BL = X XT dNT dN e L2 X = 21 21 = 0 dξ dξ 4 4 ⎡−X Y21 X 21 Y21 ⎤⎥ ⎢⎣ 21 − ⎦ 1 eT dNT dN 1 X = 2 2 λ dξ dξ L0 BN = 1 T dNT dN 1 U U 21 − V21 U 21 V21 ⎤⎥ = 2 ⎡⎢− 2 ⎣ ⎦ λ dξ dξ L0 Son las mismas expresiones ya obtenidas antes de forma directa. 51 .
El movimiento total de la barra se puede descomponer en tres fases: en primer lugar una traslación desde la posición inicial hasta hacer coincidir el nudo inicial con su posición deformada. La figura muestra un elemento plano. u1 = L − L0 52 . de tal manera que este sistema de ejes contiene el movimiento de sólido rígido del elemento y se mueve con él. placas y cáscaras. la cual queda definida por la deformación del nudo inicial u1. empleándose para el estudio de barras. y con él de tal manera que el eje x pasa por la posición deformada de ambos nudos. Figura 7. Elemento de celosía plana. tanto en su planteamiento total como actualizado. Esta limitación hace que esta formulación sea de aplicación más limitada. vigas. La definición de la posición del sistema co-rotacional de ejes para el caso de problemas en 3 dimensiones requiere de técnicas adecuadas. caracterizada por el movimiento del sistema de ejes co-rotacional y una parte de deformación del sólido con respecto a dicho sistema de ejes. emplea un único sistema de ejes global al cual se refieren todos los movimientos y deformaciones del sólido. cuya magnitud es u1 . La limitación principal de la formulación CR está en que se supone a priori que la deformación del elemento respecto del sistema co-rotacional es de pequeña magnitud comparada con el movimiento global. FORMULACIÓN CO-ROTACIONAL La formulación lagrangiana. La formulación co-rotacional (CR) por su parte emplea un sistema de ejes asociado a cada elemento de la estructura. habiéndose empleado diversos métodos para ello. El movimiento total se descompone en una parte de sólido rígido. Formulación co-rotacional.Introducción al análisis de estructuras con no linealidad geométrica 10 ELEMENTO BIARTICULADO. que en este caso consta únicamente de un alargamiento axial. seguida a continuación por una rotación de valor α hasta alcanzar la orientación deformada final. En tercer lugar se produce la deformación de la barra. en el que se define un sistema de ejes co-rotacional x .
V2 u1 U2 V1 U1 Figura 8. La variación del alargamiento de la barra a consecuencia de la variación de las deformaciones nodales es: δu1 = δU 2 cos β + δV2 sin β − δU 1 cos β − δV1 sin β Donde β es el ángulo que forma la barra con respecto al eje x en su posición deformada. por lo que en la práctica es mejor emplearla en la forma: u1 = (L − L0 ) u1 = L + L0 L2 − L2 0 = L + L0 L + L0 = 2 ⎛ T 1 ⎞ ⎟ ⎜ X21 U21 + UT 21 U21 ⎟ ⎜ ⎝ ⎠ 2 L + L0 (X21 + U21 )T (X21 + U21 ) − XT 21 X21 L + L0 10. por lo que es más fácil emplear un método más geométrico.Introducción al análisis de estructuras con no linealidad geométrica u1 = ⎡⎣(X21 + U21 )T (X21 + U21 )⎤⎦ 1/ 2 ⎤ − ⎡⎣ XT 21X21 ⎦ 1/ 2 Esta expresión tiene una mala condición numérica. Variación virtual de las deformaciones. La figura siguiente muestra la configuración una vez aplicada una variación virtual cualquiera. Esta expresión se puede poner: 53 . medido en el sistema co-rotacional: ε= Y su variación: L − L0 u = 1 L0 L0 δε = δu1 L0 Para obtener la variación del alargamiento.1 Deformación unitaria En esta formulación emplearemos la deformación unitaria ingenieril cuyo valor es. el empleo directo de su expresión resulta complejo. que consiste en imponer una variación a las deformaciones de los nudos y determinar cuánto varía el alargamiento a consecuencia de ella. Celosía plana.
10.Introducción al análisis de estructuras con no linealidad geométrica δu1 = ⎡⎢− cos β − sin β cos β ⎣ ⎧ δU 1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δ V ⎪ 1⎪ ⎪ ⎤ sin β ⎥ ⎨ ⎪ = rT δ U ⎦ ⎪δU 2 ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δ V 2 ⎪ ⎪ ⎪ ⎩ ⎭ 10.3 Matriz de rigidez tangente La matriz tangente se puede obtener derivando la expresión de las fuerzas interiores: ˆ = ∂Q = r ∂N + N ∂r = K ˆ +K ˆ K D σ ∂U ∂U ∂U • La primera matriz tangente vale: ˆ = r ∂N = r EA ∂u1 K D L0 ∂U ∂U La derivada del alargamiento axial se obtiene fácilmente a través de la expresión de su variación: δu1 = Luego: ∂u1 δ U = rT δ U ∂U ⇒ ∂u1 = rT ∂U ˆ = r ∂N = EA r rT K D L0 ∂U (79) 54 .2 Vector de fuerzas interiores El trabajo virtual de las fuerzas interiores está producido únicamente por la fuerza axial en la barra N (supuesta positiva a tracción): δWI = δu1 N = (δ U) r N T La definición del vector de fuerzas interiores es: δWI = (δ U) Q T Con lo que su valor es: ⎧ −N cos β ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ −N sin β ⎪ ⎪ ⎪ ⎪ ⎪ Q=rN =⎨ ⎬ ⎪ N cos β ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ N sin β ⎪ ⎪ ⎩ ⎭ Que es la expresión obvia de las componentes cartesianas de la fuerza axial.
Variación de la orientación. Esta expresión coincide con la matriz de rigidez habitual de la barra biarticulada plana. pero evaluada en su posición deformada. 55 . • Por su parte.Introducción al análisis de estructuras con no linealidad geométrica ⎡ c2 sc −c 2 −sc ⎤ ⎢ ⎥ ⎢ 2 2⎥ s −sc −s ⎥ ⎢ sc ˆ = EA ⎢ ⎥ K D L0 ⎢⎢−c 2 −sc c 2 sc ⎥⎥ ⎢ ⎥ ⎢−sc −s 2 sc s2 ⎥ ⎣ ⎦ Donde se ha empleado c = cos β s = sin β . Celosía plana. la matriz de rigidez geométrica es: ˆ = N ∂r = N ∂r ∂β K σ ∂U ∂β ∂ U La primera derivada es inmediata: ⎧ sin β ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ cos β − ⎪ ⎪ ∂r ⎪ ⎪=z =⎨ ⎬ ⎪ ∂β ⎪− sin β ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ cos β ⎪ ⎪ ⎩ ⎭ La segunda derivada en la expresión de la matriz geométrica es la derivada del ángulo de orientación β respecto a las deformaciones de los nudos y se obtiene fácilmente estudiando la variación de dicho ángulo: δβ = δv2 1 = (δV2 cos β − δU 2 sin β − (δV1 cos β − δU 1 sin β )) L L δβ = 1 T z δU L ∂β 1 = zT ∂U L v2 V2 U2 Por lo tanto la derivada buscada es: δβ = 1 T ∂β δU z δU = L ∂U ⇒ V1 U1 Figura 9.
Introducción al análisis de estructuras con no linealidad geométrica La matriz de rigidez geométrica del elemento en esta formulación es: ⎡ s 2 −sc −s 2 sc ⎤ ⎢ ⎥ ⎢ 2 2⎥ sc c sc c − − ⎢ ⎥ ˆ = N z zT = N ⎢ ⎥ K σ L L ⎢⎢−s 2 sc s 2 −sc ⎥⎥ ⎢ ⎥ ⎢ sc −c 2 −sc c 2 ⎥ ⎣ ⎦ (80) 56 .
El movimiento total de la barra se puede descomponer en tres fases: en primer lugar una traslación desde la posición inicial hasta hacer coincidir el nudo inicial con su posición deformada. no hay que considerar deformación lateral de la barra. En tercer lugar se produce la deformación de la barra. En segundo lugar una rotación de valor α hasta alcanzar la orientación deformada final del eje co-rotacional x . que en este caso consta de dos efectos: un alargamiento axial en la dirección del eje x y una deformación por flexión. y co-rotacional con él. Formulación co-rotacional. Viga plana. de tal manera que el eje x pasa siempre por la posición deformada de ambos nudos. la cual queda definida por la deformación del nudo inicial u1. Al haberse tomado los ejes co-rotacionales pasando por la posición deformada de los nudos. caracterizada por los giros de los dos extremos θ1 y θ2. FORMULACIÓN CO-ROTACIONAL La formulación de este elemento emplea un sistema de ejes x .Introducción al análisis de estructuras con no linealidad geométrica 11 ELEMENTO VIGA PLANA. x y 2 1 u2 u1 L0 0 Figura 10. pues ésta está tenida en cuenta en la rotación α.1 Deformación axial y esfuerzo axial { θ2 } T La deformación axial en el sistema co-rotacional es u1 y su determinación es exactamente igual que para el elemento de celosía. L + L0 L2 − L2 0 u1 = (L − L0 ) = L + L0 L + L0 u1 = (X21 + U21 )T (X21 + U21 ) − XT 21 X21 L + L0 = 2 ⎛ T 1 ⎞ ⎟ ⎜X21 U21 + UT 21 U21 ⎟ ⎜ ⎠ 2 L + L0 ⎝ 57 . Los grados de libertad del elemento son: U = U 1 V1 θ1 U 2 V2 11.
medida en el sistema corotacional. suponiendo un comportamiento elástico es: N = Aσ = AE ε = AEu1 L0 11. vale ε= L − L0 u = 1 L0 L0 El valor del esfuerzo axial (supuesto positivo a tracción) producido por esta deformación. por lo que la energía de flexión está asociada únicamente a los giros de los nudos relativos a dicho eje co-rotacional θ1 y θ2 . Deformaciones de flexión. medidos respecto de la orientación inicial de la barra.Introducción al análisis de estructuras con no linealidad geométrica La deformación unitaria ingenieril debida al alargamiento axial. los momentos flectores en ambos extremos de la barra se relacionan con los giros correspondientes mediante la ecuación de rigidez: ⎧θ1 ⎪ ⎫ ⎧ M ⎫ EI ⎡ 4 2⎤ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 1⎪ ⎪= ⎪ ⎢ ⎥ ⎨ ⎬ ⎨ ⎬ ⎢ ⎥ ⎪ ⎪ ⎪ ⎪ 2 4 M θ ⎪ 2⎭ ⎪ L0 ⎢⎣ ⎥⎦ ⎪ ⎩ ⎪ 2⎪ ⎪ ⎩ ⎭ 58 . La definición empleada para el eje co-rotacional x hace que no haya deformaciones laterales en los nudos. Estos giros valen: θ1 = θ1 − α θ2 = θ2 − α Figura 11. θ1 y θ2. Viga plana. Empleando la teoría de Euler – Bernouilli de la flexión.2 Deformación y momentos de flexión La deformación producida por la flexión de la viga queda definida por los dos giros en los extremos. Estos giros se suponen de pequeña magnitud.
Introducción al análisis de estructuras con no linealidad geométrica 11. como se efectuó para el elemento biarticulado: δα = δβ = 1 δv2 = (δV2 cos β − δU 2 sin β − (δV1 cos β − δU 1 sin β )) L L δα = 1 T z δU L Siendo ahora: z = sin β − cos β { 0 − sin β cos β 0 } T La variación de los ángulos de rotación relativos en ambos extremos es: ⎧ ⎫ ⎪δθ1 ⎪ ⎪ ⎪ ⎧δθ ⎫ ⎪ ⎧δα⎫ ⎪ ⎪ 1⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ ⎬ = ⎨ ⎬−⎨ ⎬ = ⎪ ⎪ ⎪δα⎪ δθ ⎪ ⎪ ⎪ ⎩ ⎪δθ2 ⎪ ⎪ ⎭ ⎪ ⎩ ⎪ ⎭ ⎪ 2⎭ ⎪ ⎩ Siendo: ⎡ zT ⎤ ⎡ 0 0 1 0 0 0⎤ 1 ⎢ ⎥ ⎢ ⎥ ⎢ 0 0 0 0 0 1⎥ δ U − L ⎢ zT ⎥ δ U = A δ U ⎢⎣ ⎥⎦ ⎢⎣ ⎥⎦ ⎡0 0 1 0 0 0⎤ 1 ⎡ zT ⎤ ⎥− ⎢ ⎥ A = ⎢⎢ 0 0 0 0 0 1⎥⎥ L ⎢⎢ zT ⎥⎥ ⎣⎢ ⎦ ⎣ ⎦ Agrupando las tres deformaciones virtuales en un vector δ p se puede definir la matriz B del elemento: ⎧ ⎪ ⎪ δu1 ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δ p = ⎨ δθ1 ⎪ ⎬= ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δθ ⎪ 2 ⎪ ⎪ ⎪ ⎩ ⎭ ⎡ rT ⎤ ⎢ ⎥ δU = B δU ⎢A⎥ ⎢⎣ ⎥⎦ 59 . esta variación se puede poner como: δu1 = ⎡⎢− cos β − sin β ⎣ ⎧δU 1 ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δV1 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δθ1 ⎪ ⎪ ⎪ ⎪ = rT δ U 0 cos β sin β 0⎤⎥ ⎪ ⎨ ⎦ ⎪δU 2 ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δV2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ δθ2 ⎪ ⎪ ⎩ ⎪ ⎭ Siendo β el ángulo que forma la barra con respecto al eje x en su posición deformada final. La variación del ángulo de orientación β al variar las deformaciones de los nudos se obtiene fácilmente por consideraciones geométricas.3 Deformaciones virtuales La variación del alargamiento axial es: δε = δu1 L0 Efectuando el mismo desarrollo que para el elemento biarticulado.
4 Trabajo virtual El trabajo virtual de las fuerzas interiores está producido por la fuerza axial y los dos momentos en los extremos. la primera matriz tangente vale: ˆ = BT ∂Q = BT C B K D ∂U (81) 60 .5 Matriz de rigidez tangente Se puede obtener derivando la expresión de las fuerzas interiores: ∂Q ∂BT T ∂Q ˆ ˆ +K ˆ =B + K= Q=K D σ ∂U ∂U ∂U • La relación entre los 3 esfuerzos internos y las deformaciones locales se puede poner en forma compacta como: ⎧N ⎪ ⎫ ⎡ ⎪ ⎪ ⎪ ⎢A 0 ⎪ ⎪ ⎪ ⎪ ⎪M ⎬ ⎪ E ⎢⎢ 0 4I ⎨ 1 = ⎪ ⎪ L0 ⎢ ⎪ ⎪ ⎢ 0 2I ⎪ ⎪ M ⎪ ⎢⎣ ⎪ 2⎪ ⎪ ⎩ ⎭ Q = Cp Y por lo tanto su derivada es: ⎧u ⎫ ⎪ 0 ⎤⎥ ⎪ 1⎪ ⎪ ⎪ ⎪ ⎪ ⎥⎪ ⎪ ⎪ 2I ⎥ ⎨θ1 ⎬ ⎪ ⎥⎪ ⎪ ⎪ ⎪ 4I ⎥⎥ ⎪ ⎪θ2 ⎪ ⎪ ⎦⎩ ⎪ ⎭ ∂Q ∂p =C = CB ∂U ∂U El segundo paso es evidente en base a la definición de la δ p . actuando sobre sus correspondientes deformaciones virtuales: δWI = δu1 N + δθ1 M 1 + δθ2 M 2 Las tres fuerzas interiores independientes se agrupan en un vector: ⎧N ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ Q = ⎨M 1 ⎪ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ M ⎪ 2⎪ ⎪ ⎭ ⎪ ⎩ Con lo que el trabajo virtual de dichas fuerzas interiores se puede poner como: δWI = δu1 N + δθ1 M 1 + δθ2 M 2 = δ pT Q = (δ U)T BT Q El vector de fuerzas interiores del elemento es: δWI = (δ U) BT Q = (δ U) Q T T Q = BT Q 11.Introducción al análisis de estructuras con no linealidad geométrica 11. En consecuencia.
la matriz de rigidez geométrica es: T T T T ˆ = ∂B Q = ∂B1 N + ∂B2 M + ∂B3 M K 1 2 σ ∂U ∂U ∂U ∂U Donde Bi es la fila i-sima de la matriz B. la primera derivada es: ∂BT ∂r ∂r ∂β 1 = = ∂U ∂U ∂ β ∂ U La derivada del vector r es inmediata y corresponde al vector z ya conocido: ⎧ sin β ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ − cos β ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 0 ⎪ ⎪ ∂r ⎪ ⎪=z =⎨ ⎬ ⎪ − β sin ∂β ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ β cos ⎪ ⎪ ⎪ ⎪ ⎪ 0 ⎪ ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ La derivada del ángulo de orientación β respecto a las deformaciones de los nudos se puede obtener a partir de su variación que ya ha sido calculada anteriormente: δβ = Por lo tanto la derivada buscada es: 1 T z δU L ∂β 1 = zT ∂U L δβ = ∂β 1 T z δU = δU L ∂U ⇒ En consecuencia. el primer sumando de la matriz geométrica es: ∂BT ∂r 1 1 = = z zT ∂U ∂U L  La derivada de la segunda fila requiere derivar la primera fila de A transpuesta: ∂BT ∂AT ∂ ⎛ 1 ⎞ 1 ∂z ∂(1/ L) 2 1 ⎜ = = − z⎟ =− −z ⎟ ⎜ ∂U ∂U ∂U ⎝ L ⎠ L ∂U ∂U ∂BT 1 ∂z ∂ β 1 ∂L 2 =− +z 2 ∂U L ∂β ∂U L ∂U 61 .Introducción al análisis de estructuras con no linealidad geométrica Esta expresión coincide con la matriz de rigidez habitual de la viga plana. y se requiere obtener sus derivadas respecto de los grados de libertad del elemento.  Por su parte.  Considerando que B1 = rT . pero evaluada en su posición deformada.
Introducción al análisis de estructuras con no linealidad geométrica Las derivadas necesarias son: ⎧ cos β ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ β sin ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 0 ⎪ ⎪ ∂z ⎪ ⎪ = −r =⎨ ⎬ ⎪ ∂β ⎪− cos β ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ − β sin ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ 0 ⎪ ⎪ ⎪ ⎪ ⎩ ⎭ ∂L ∂u = 1 = rT ∂U ∂U Con lo que la derivada de la segunda fila queda: ∂BT 1 1 2 = 2 r zT + 2 z rT ∂U L L  La derivada de la tercera fila es igual que la de la segunda: ∂BT 1 1 3 = 2 r zT + 2 z rT ∂U L L  Agrupando los distintos valores se obtiene la expresión de la matriz de rigidez geométrica del elemento en esta formulación: ˆ = N z zT + (M 1 + M 2 ) r zT + z rT K σ L L2 ( ) (82) 62 .
Introducción al análisis de estructuras con no linealidad geométrica 12 FLEXIÓN DE PLACAS. y por dos rotaciones θx . Campo de deformaciones en una placa. que se suman a las producidas por las fuerzas contenidas en su plano. θy alrededor de los ejes X. v y uno perpendicular a ella w). Estos desplazamientos laterales dan lugar a su vez a deformaciones unitarias en el plano de la placa. Empleando la teoría de MindlinReissner. Se supone que las fuerzas transversales producen una flexión lateral con unos desplazamientos laterales de magnitud suficiente para no ser despreciables. Y.1 Campo de deformaciones Estudiamos la flexión de placas inicialmente planas. estas rotaciones no son las derivadas de la deformación transversal. despreciando en él los términos cuadráticos en las deformaciones contendidas en el 63 . sometidas a fuerzas tanto transversales como contenidas en el plano de la placa.2 Deformaciones unitarias Se considera una versión degenerada del tensor de deformaciones unitarias de Green – Lagrange. Por lo tanto el problema tiene 5 deformaciones incógnitas: uT = u v w θx { θy } Nota: con objeto de simplificar la notación se emplea la nomenclatura clásica para las coordenadas x ≡ x 1 y ≡ x 2 y para las deformaciones u ≡ u1 v ≡ u2 w ≡ u 3 . con lo que en el estado deformado la placa deja de estar contenida en su plano inicial. El campo de deformaciones en el plano medio de la placa está compuesto por tres desplazamientos (dos contenidos en el plano de la placa u. Las deformaciones en un punto P situado a una distancia z del plano medio son: uP = u + z θy vP = v − z θx wP = w w uP Z Y X x vP w v u y Figura 12. 12. FORMULACIÓN LAGRANGIANA TOTAL 12.
y corresponden a la diferencia entre el giro y la derivada de la deformación transversal: ⎧ ∂v P ⎫ ⎧ ⎫ ∂w P ⎪ ∂w ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ + −θx + ⎪ ⎪ ⎪ ⎪ γ ⎧ ⎫ yz ⎪ ⎪ z y y ⎪ ∂ ∂ ⎪ ⎪ ∂ ⎪ ⎪ ⎪ γ=⎪ ⎨γ ⎪ ⎬=⎪ ⎨ ⎬=⎪ ⎨ ⎬ ⎪ ⎪ ∂w P ⎪ ∂u ⎪ ⎪ ∂w ⎪ xz ⎪ ⎪ ⎪ ⎩ ⎭ ⎪ + P⎪ + θy ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂z ⎪ ⎪ ⎪ ⎩ ∂x ⎭ ⎪ ⎩ ∂x ⎭ Efectuando el mismo desarrollo que en régimen lineal las deformaciones unitarias de cortadura se pueden expresar en función de las deformaciones de los nudos: ⎡ ∂ ⎢0 0 ⎢ ∂y γ=⎢ ⎢ ∂ ⎢0 0 ∂x ⎣⎢ ⎧u ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ⎤⎪ v ⎪ ⎪ ⎪ −1 0 ⎥ ⎪ ⎪ ⎥ ⎪w ⎪ ⎥⎨ ⎪ ⎬ = ∂Q u = ∂Q N U = BQ U ⎪ ⎥⎪ ⎪θx ⎪ ⎪ 0 +1⎥ ⎪ ⎪ ⎦⎥ ⎪ ⎪ ⎪ ⎪ θ⎪ ⎪ ⎪ y⎪ ⎪ ⎩ ⎭ −N k 0 ⎤ 0⎥ ⎥ ⎥ ⎥ Nk ⎥ ⎦⎥ ⎡ ∂N k ⎢0 0 ⎢ ∂y k BQ =⎢ ⎢ ∂N k ⎢0 0 ⎢⎣ ∂x Su variación es: δγ = BQ δ U 67 .Introducción al análisis de estructuras con no linealidad geométrica ⎡ ⎢0 0 0 0 ⎢ ⎢ ⎢ ∂N k ⎢ 0 0 − Bk F = 0 ⎢ ∂y ⎢ ⎢ ∂N k ⎢0 0 0 − ⎢⎣ ∂x ∂N k ⎤⎥ ∂x ⎥⎥ ⎥ 0 ⎥ ⎥ ⎥ ∂N k ⎥ ⎥ ∂y ⎥⎦ • La matriz BF es lineal en las derivadas de la deformación vertical. Cada uno de sus bloques tiene la forma: ⎡ ∂w ∂N k ⎢0 0 ⎢ ∂x ∂ x ⎢ ⎢ ∂w ∂N k ⎢ 0 Bk S = 0 ⎢ ∂y ∂y ⎢ ⎢ ∂w ∂N k ∂w ∂N k ⎢0 0 + ⎢⎣ ∂y ∂ x ∂x ∂y ⎤ 0 0⎥⎥ ⎥ ⎥ 0 0⎥ ⎥ ⎥ ⎥ 0 0⎥ ⎥⎦ 12.4 Deformaciones unitarias de cortadura Su expresión es la misma que en el régimen lineal.
Introducción al análisis de estructuras con no linealidad geométrica 12. sustituyendo las variaciones de las deformaciones unitarias por sus valores en función de las deformaciones nodales: Q= 12. El incremento en la variación de la deformación unitaria sólo depende de la variación de la matriz A: 68 .Kirchhoff y τ el vector de tensiones de cortadura.6 Vector de fuerzas interiores τ =G γ Su expresión se obtiene directamente del principio del trabajo virtual. Por lo tanto la matriz de rigidez tangente es: ˆ = K D ∫ h (BM + BS )T C (BM + BS )dA + ∫ h3 T BF C BF dA 12 • La segunda integral proporciona la matriz de rigidez geométrica. cuyo valor es: ⎧ τyz ⎪ Gγ ⎫ ⎪ ⎧ ⎫ ⎪ ⎪ yz ⎪ ⎪ ⎪ ⎪ = ⎨τ ⎬ ⎨ ⎬ ⎪ ⎪ ⎪ ⎪ ⎪ xz ⎭ ⎪ ⎩ ⎩ ⎪G γxz ⎭ ⎪ 12. las dos últimas integrales son nulas.5 Trabajo virtual interior Teniendo en cuenta los dos tipos de deformaciones unitarias existentes.7 Matriz de rigidez tangente ∫B T S dv + ∫B T Q τ dv El incremento del trabajo virtual interior es: Δ (δWI ) = ∫ (δE) ∫ (δE) T ΔS dv + ∫ Δ(δ E)T S dv + ∫ (δγ)T Δτ dv • La primera integral proporciona el primer sumando de la matriz de rigidez tangente: ∫ (δE) T ΔS dv = T ˆ = (δ U)T K ˆ U ˆ C ΔE dv = (δ U)T ∫ BT C B dv U D El valor de esta matriz es: ˆ = K D ∫ B C B dv = ∫ (B T T M S T M S M + z BF + BS )T C (BM + z BF + BS )dv M ˆ = K D ∫ (B + B ) C (B ∫ z (B + B ) C B + BS )dv + ∫ z 2 BFT C BF dv + dv + ∫ z BFT C (BM + BS )dv F Efectuando la integración en la coordenada z. su valor es: δWI = ∫ δE T S dv + ∫ δγT τ dv Siendo S el vector de tensiones de Piola . Para ello se desarrolla su integrando. pues z está medida desde el centro de gravedad y los demás términos del integrando son independientes de z.
Introducción al análisis de estructuras con no linealidad geométrica Δ(δ E) = 0 + 0 + Δ(A) δ H ˆ: El incremento de A. el integrando vale: ⎧ ⎪ ∂δw Δ(δ E)T S = ⎪ ⎨ ⎪ ⎪ ∂x ⎩ ∂δw ⎫ ⎪ ⎡⎢S11 ⎪ ⎬⎢ ∂y ⎪ ⎪ ⎣⎢S12 ⎭ ⎧ ˆ⎫ ∂w ⎪ ⎪ ⎪ ⎪ ⎪ S12 ⎤ ⎪ ⎪ ∂ x ⎪ ⎪ ⎪ ⎥⎨ ⎬ ⎥ ˆ ∂ w ⎪ ⎪ S 22 ⎥ ⎪ ⎪ ⎦⎪ ⎪ ⎪ ⎪ ⎪ ∂y ⎭ ⎪ ⎩ Las derivadas del incremento de la flecha lateral se pueden poner. cada uno de los cuales es de la forma: ⎡ ∂N k ⎢0 0 ⎢ ∂x k BW =⎢ ⎢ 0 0 ∂N k ⎢ ∂y ⎢⎣ De forma análoga. sólo depende del incremento de la deformación lateral w ⎡ ˆ ⎢ 0 0 0 0 ∂w ⎢ ∂x ⎢ ⎢ ΔA = ⎢ 0 0 0 0 0 ⎢ ⎢ ⎢ ˆ ∂w ⎢0 0 0 0 ⎢⎣ ∂y 0 ˆ ∂w ∂y ˆ ∂w ∂x ⎤ 0 0 0 0⎥⎥ ⎥ ⎥ 0 0 0 0⎥ ⎥ ⎥ ⎥ 0 0 0 0⎥ ⎥⎦ El incremento en la variación de la deformación unitaria es por lo tanto: ⎧ ⎫ ⎪ ⎪ ˆ ∂δw ∂w ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x x ∂ ∂ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ˆ ∂δw ∂w ⎪ ⎪ Δ(δ E) = ⎨ ⎬ ⎪ ⎪ ∂y ∂ y ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ˆ ˆ δ δ w w w w ∂ ∂ ∂ ∂ ⎪ ⎪ ⎪ ⎪ + ⎪ ⎪ y x x y ∂ ∂ ∂ ∂ ⎪ ⎪ ⎩ ⎭ Finalmente. y está compuesta por n bloques. se cumple que: ⎤ 0 0⎥ ⎥ ⎥ 0 0⎥⎥ ⎥⎦ ⎧ ∂δw ⎪ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂x ⎪ ⎪ ⎪ ⎨ ⎬ = BW δ U ∂ δ w ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ y ∂ ⎪ ⎪ ⎩ ⎭ 69 . empleando la interpolación de deformaciones en la forma: ∂N k ˆ ⎫ ⎧ ∂w ⎫ ⎧ ˆ⎪ ⎪ ⎪ ⎪ ⎪ Wk ⎪ ⎪ ⎪ ∑ ⎪ ⎪ ⎪ ⎪ ⎪ ∂ x ⎪ ⎪ ∂ x ⎪ k ⎪ ⎪ ⎪ ⎪ ˆ ⎨ ˆ⎬=⎨ ⎬ = BW U ∂ w ⎪ ⎪ ⎪ ⎪ N ∂ k ˆ ⎪ ⎪ ⎪ ⎪ Wk ⎪ ⎪ ⎪ ⎪ ∑ ⎪ ⎪ ⎪ ⎪ ∂ y ∂ y ⎪ ⎪ ⎩ ⎭ ⎩ ⎪ k ⎪ ⎭ La matriz BW es de tamaño 2 x 5n.
que tiene la misma expresión que en caso lineal. se obtiene: ˆ = K σ ∫ ⎡N 11 (BW )T ⎢⎢ ⎢⎣N 12 N 12 ⎤ ⎥ B dA N 22 ⎥⎥ W ⎦ Donde N11 y N22 son los esfuerzos axiales y N12 el esfuerzo cortante en el plano de la misma. ∫ (δγ) T ˆ Δτ dv = (δ U)T ∫ BT Q G BQ dv U ˆ = K Q ∫B T Q G BQ dv La matriz de rigidez tangente del elemento es: ˆ =K ˆ +K ˆ +K ˆ K D σ Q 70 . Esfuerzos en el plano en una placa. • La tercera integral proporciona la matriz de rigidez asociada al esfuerzo cortante. por unidad de anchura: N ij = ∫ S dz ij S22 S11 Z Y X S12 N22 N12 N11 Figura 13.Introducción al análisis de estructuras con no linealidad geométrica Con lo que el incremento del trabajo virtual interior se puede poner como: ∫ ⎡S11 S12 ⎤ ⎥ B dv U ˆ Δ(δ E)T S dv = (δ U)T ∫ (BW )T ⎢⎢ ⎥ W S S 12 22 ⎥⎦ ⎣⎢ ⎡S11 S12 ⎤ ⎥ B dv (BW )T ⎢⎢ ⎥ W S S 12 22 ⎥⎦ ⎣⎢ Esta expresión define la matriz de rigidez geométrica: ˆ = K σ ∫ Efectuando la integración en z. dado que BW no depende de z.
Si se desea aplicar una cantidad variable de carga en cada incremento.2. de tal manera que al final del paso de carga n se obtienen las deformaciones en el instante n+1: Un +1 = Un + ΔUn (84) siendo ΔUn el incremento de la deformación que se produce en el paso de carga n. Cada paso de carga de la secuencia se identifica mediante un subíndice n=0. Las cargas aplicadas en el paso n se denominan Pn y puede considerarse que el incremento de carga aplicado en cada caso es constante o variable. Pet el vector de fuerzas exteriores y siendo K et Q el vector de fuerzas interiores (se ha añadido el superíndice e para indicar que se trata de magnitudes propias de un elemento). De esta manera se obtiene toda la respuesta de la estructura ante un sistema de cargas creciente. que define los valores relativos entre las distintas componentes de la fuerza. se ha obtenido en la forma: ˆe U ˆ e = Pe (t +Δt ) − Qe K ˆ et la matriz de rigidez tangente del elemento... por incrementos y en cada uno de dichos incrementos se busca el estado de equilibrio. y P es un vector de fuerzas de referencia. en un solo paso entre el instante inicial t =0 y el instante final no es posible en general. La resolución de la ecuación anterior para la carga total aplicada. 71 . linealizada en un instante cualquiera t. El vector U ˆ es la el incremento de deformación entre t y t + Δt en todos los nudos de la estructura y K matriz de rigidez tangente en el instante t.Introducción al análisis de estructuras con no linealidad geométrica 13 RESOLUCIÓN DE LAS ECUACIONES INCREMENTALES La ecuación incremental de equilibrio de un elemento finito. ésta se representa en la forma: Pn = λn P (85) En este caso λ es un parámetro sin dimensiones que define el valor real de la fuerza en el paso n.1. en el que las cargas se van aplicando de forma paso a paso. Ensamblando los términos correspondiente a los distintos elementos finitos se obtiene la ecuación incremental de equilibrio de toda la estructura: ˆ U ˆ = Pt +Δt − Q K (83) El término independiente contiene las fuerzas exteriores aplicadas P en el instante t + Δt y ˆ contiene las fuerzas interiores Q en los elementos de la estructura en el instante t. Esta expresión es válida tanto para el planteamiento total como para el actualizado. En el primer caso. y sólo cambian en ellos los valores concretos de la matriz tangente y del vector de fuerzas interiores. La obtención de la respuesta de un sistema no lineal se efectúa en la práctica empleando un proceso incremental. la carga en un paso cualquiera será: Pn = n PP Siendo PP la carga aplicada en cada paso.
1 Método incremental puro En este método se efectúa la fase de predicción de las deformaciones en el incremento de carga. y en cualquier caso es más ventajoso emplear los métodos que se explican a continuación. Kn Pn+1-Qn Pn Un pred Un Un Un+1 U Figura 14. 13. pues no se satisface el equilibrio en los distintos puntos obtenidos. Puede mejorarse fácilmente a base de calcular el residuo no equilibrado en cada iteración y añadirlo a las fuerzas a aplicar al siguiente incremento. Esta predicción va seguida seguido de un proceso de corrección de las deformaciones. que se va acumulando a medida que se aplican nuevos incrementos de carga. pero no se efectúa ningún proceso de corrección del error cometido. Q(U) Pn+1 corr. al final del cual lógicamente no habrá equilibrio. que se encarga de buscar el equilibrio al final de dicho paso de cargas. En función de cómo se haga la fase de corrección se plantean diversos métodos. los cuales permiten garantizar el equilibrio.Introducción al análisis de estructuras con no linealidad geométrica Todos los métodos de resolución combinan el proceso incremental de aplicación de cargas con un proceso iterativo dentro de cada paso de cargas. Se trata de un método no exacto. Al ser el error acumulativo. 72 . Fase de predicción y corrección. El proceso iterativo consta de un primer paso de predicción del incremento de deformaciones producido por el incremento de cargas aplicado. Método incremental. es decir para t ≡ tn : ˆ ΔU = P − Q K n n n +1 n (86) El incremento de deformación así obtenido tiene un error. hasta satisfacer el equilibrio en la nueva posición. sólo puede emplearse con incrementos de carga muy pequeños. El incremento de deformación ΔUn producido en un incremento de carga se calcula apoyándose en la ecuación incremental al comienzo de dicho paso de carga. lo cual permite estimar el error producido.
al final de la iteración k. apoyándose en la solución conocida en el instante anterior n. la derivada del residuo sólo corresponde a la −1 derivada de las fuerzas internas Qk n +1 : −1 ⎞ ⎛ ∂Qk k −1 n +1 ⎟ ˆ k ⎜ ⎟U = 0 Rk ≈ R − ⎜ n +1 n +1 ⎟ ⎟ ⎜ ⎝ ∂U ⎠ La derivada de las fuerzas internas respecto de las deformaciones proporciona la matriz de ˆ k −1 : rigidez tangente K n +1 −1 ∂Qk n +1 ˆ k −1 =K n +1 ∂U 73 .Introducción al análisis de estructuras con no linealidad geométrica Pn+2 Pn+1 Error Error Pn Un U Un Un+1 Un+2 Figura 15. Si las fuerzas no dependen de la deformación.…).2. 13. mediante una secuencia de iteraciones (k=1.2 Método de Newton-Raphson En este método se emplea un proceso iterativo completo de predicción – corrección hasta alcanzar el equilibrio en el instante n+1. Método incremental puro. Las deformaciones en el instante n+1 se calculan por aproximaciones sucesivas. Para un instante cualquiera (k) de la iteración las ecuaciones de equilibrio se pueden poner en forma de residuo: k Rk n +1 = P n +1 − Qn +1 = 0 Desarrollando este residuo en serie alrededor del punto de iteración anterior (k-1) se obtiene: −1 ⎞ ⎛ ∂R k k −1 n +1 ⎟ ˆ k ⎜ ⎟ Rk ≈ R + U =0 ⎜ n +1 n +1 ⎟ ⎟ ⎜ ⎝ ∂U ⎠ ˆ k es el incremento de deformación producido en la iteración k: En esta expresión U ˆ k = Uk − Uk −1 U n +1 n +1 (87) y Uk n +1 es la estimación de las deformaciones en el instante n+1.
que debe efectuarse en cada paso de la iteración. Figura 16.3 Método de Newton modificado En el método de Newton-Raphson. Método de Newton-Raphson. pero no en el estado que se toma como inicio en el incremento. 13. la parte más costosa es la factorización de la matriz de rigidez tangente. pues la diferencia entre ambas formulaciones está en la situación que se toma como referencia para las distintas magnitudes. que siempre es el último conocido. se plantea el método de Newton-Raphson modificado. en el cual la matriz de rigidez tangente en la primera iteración se utiliza en todas las iteraciones posteriores. tanto para la formulación total como para la actualizada. Como alternativa e dicho método. Nótese que ambas magnitudes se evalúan para los últimos valores actualizados de las deformaciones calculados a medida que progresa la iteración (al final de la iteración anterior). Por lo tanto la ecuación de la iteración es: 74 .Introducción al análisis de estructuras con no linealidad geométrica Por lo tanto la ecuación a resolver en la iteración k es: k −1 ˆ k −1 ˆ k Rk n +1 ≈ P n +1 − Qn +1 − Kn +1 U = 0 ˆ k −1 U ˆ k = P − Qk −1 K n +1 n +1 n +1 (88) ˆ k −1 y el vector de fuerzas interiores Qk −1 están En esta ecuación la matriz tangente K n +1 n +1 evaluados para la última estimación conocida (k-1) de las deformaciones en el instante n+1 −1 que son las de la iteración anterior Uk n +1 . Como condiciones para el comienzo de la iteración se emplean las del último estado de equilibrio conocido: 0 Un +1 ≡ Un 0 Qn +1 = Qn Nótese que la iteración se inicia apoyándose en el último estado de equilibrio conocido. no para los valores al inicio de la misma.
basados en dos ideas: el control de fuerza y el control de desplazamiento. 13. Los principales problemas que plantean los algoritmos de control de desplazamiento son la 75 . Los algoritmos de control de fuerza corresponden a lo ya explicado anteriormente.Introducción al análisis de estructuras con no linealidad geométrica ˆ0 U ˆ k = P − Qk −1 K n +1 n +1 n +1 Este método tiene un menor costo computacional en cada iteración. a base de aumentar la fuerza exterior paulatinamente. La determinación de la curva fuerza deformación completa en estos casos requiere por del empleo de técnicas que permitan identificar la presencia de un punto límite y pasarlo eficazmente. pero se aumenta el número de iteraciones necesarias para alcanzar la convergencia. pero al llegar al punto B este método fallará también. pues la iteración de Newton falla en las proximidades de los puntos límite. el sistema no tiene ese tipo de respuesta monótona. En estos casos no es posible aplicar una estrategia simple de aumentar de forma continua la carga. es decir que se puede aumentar la carga aplicada y se obtiene un aumento de deformación. este método puede ser más ventajoso que el inicial o no. incluso simples. En una estructura cuya respuesta sea como la de la figura estos métodos fallan en las proximidades del punto A. A partir del punto A puede emplearse un método de control de desplazamiento. Se han desarrollado muchos algoritmos que permiten pasar puntos límites. Método de Newton-Raphson modificado.4 Métodos restringidos La combinación del proceso incremental de carga y de la iteración de Newton es muy eficiente para obtener la respuesta de sistemas no lineales cuando ésta es creciente. sino que existen puntos límites en los que la respuesta pasa de ser creciente a decreciente o viceversa y la estructura muestra fenómenos de snap-through o snapback. Figura 17. Dependiendo del costo de la factorización y de las restantes operaciones. Sin embargo en muchísimas aplicaciones.
Existen varios métodos restringidos. dentro de la curva de respuesta del sistema. Para resolver estos problemas se han desarrollado los denominados métodos restringidos. Control de fuerza y de desplazamiento. φ (ΔUk n . es decir que los nuevos incrementos de deformación se buscan en la intersección con dicha perpendicular a la tangente. que se diferencian en la ecuación de restricción que añaden al sistema. es decir el nivel de carga. se considera dicho nivel de carga. 13. y su dificultad para tratar fenómenos de snap-back. se limita el máximo incremento a efectuar por medio de Δs y de la ecuación de restricción se determina la λ (es decir la carga) a aplicar. aunque ésta muestre cambios de dirección.4. La idea fundamental en que se apoyan es modificar el nivel de carga aplicada en cada paso del proceso incremental de carga en vez de mantenerlo constante. Δs ) = 0 Así pues en estos métodos restringidos. La ecuación de restricción φ relaciona el incremento de desplazamiento que es posible alcanzar en cada iteración ΔUk n con alguna distancia máxima en la curva de respuesta de la estructura Δs.1 Método del plano normal En este método la iteración para obtener el nuevo equilibrio en el instante n+1 (es decir la fase de corrección) se efectúa sobre la perpendicular a la tangente al último equilibrio alcanzado n. Figura 18. que está definido por el parámetro λ. 76 .Introducción al análisis de estructuras con no linealidad geométrica elección del desplazamiento usado para controlar el proceso. la cual se introduce como dato en el método. Para poder modificar el nivel de carga aplicada. El valor de esta variable. se determina añadiendo una ecuación de restricción que obligue al método iterativo a moverse hacia la posición de equilibrio. como una variable más del problema. En todos ellos la variable que define el nivel de carga se actualiza en cada iteración del proceso en la forma: k k −1 ˆk λn = λn +λ ˆk el incremento del parámetro que define la carga en la iteración k del paso de carga Siendo λ n.
Método del plano normal actualizado. 13. En la práctica pueden emplearse varios de ellos. con lo que se consigue localizar mejor los puntos límites. P s U Figura 20. 13. que se basan en comparar la norma de alguna magnitud con algún valor de referencia considerado despreciable.2 Método del plano normal actualizado Este método es una variante del anterior.5 Criterios de convergencia Para terminar la iteración de búsqueda del equilibrio es necesario emplear un criterio adecuado. que indique que se ha llegado a la convergencia de la solución. Es decir: ˆk U 2 ≤ εD Un +1 2 77 . • El método más simple consiste en imponer que la norma del incremento de desplazamiento producido en una iteración sea muy inferior a la norma de la deformación total al final del caso de carga.4. y en él la iteración para obtener el nuevo equilibrio se efectúa sobre la perpendicular a la tangente en la última iteración efectuada k-1. Método del plano normal.Introducción al análisis de estructuras con no linealidad geométrica P s U Figura 19.
Introducción al análisis de estructuras con no linealidad geométrica Donde εD es la tolerancia. Por eso resulta interesante introducir un criterio basado en la fuerza no equilibrada durante la iteración. pero no garantiza el equilibrio de fuerzas. se aproxima por el último valor de ella que se haya calculado: ˆk U 2 ≤ εD Uk n +1 2 En algunos casos la solución obtenida con este método puede estar lejos de la convergencia. como ocurre cuando la deformación cambia muy poco en cada paso de carga. se puede usar un criterio en el que se evalúa el incremento de energía interna en cada iteración (es decir el trabajo hecho el incremento de deformación y por las fuerzas no equilibradas). • Para evitar los problemas de los métodos anteriores. Como la deformación al final del paso de carga no es conocida. y se compara con el incremento de energía interna inicial en el paso de carga: ˆ ) (P (U k T n +1 ˆ1 − Qk n +1 ) ≤ εE (U ) (P n +1 − Qn ) T 78 . • El cumplimiento del criterio anterior garantiza en todo caso que las deformaciones cambian poco. pero continúa cambiando durante muchos pasos. Por ejemplo se puede imponer que la norma del residuo al final de la iteración sea despreciable frente al residuo con el que se comenzó la iteración: Rk n +1 k Pn +1 − Qn +1 2 0 ≤ εR Rn +1 2 2 ≤ εR Pn +1 − Qn 2 El principal problema de este método es que no considera la contribución de la deformación al criterio de terminación. pero las deformaciones sigan aumentado en cada paso de carga. como ocurre en el caso de materiales con un módulo de endurecimiento por deformación muy bajo en los que las fuerzas cambien muy poco.
Las fuerzas exteriores. y buscar los nuevos incrementos de deformación en la intersección con dicho círculo. que son desconocidas. pueden expresarse en la forma: Pn +1 = λn +1 P Donde λ es un parámetro sin dimensiones que define el valor real actual de la fuerza. El residuo queda: Rn +1 = λn +1P − Qn +1 = 0 Esta ecuación se debe satisfacer en cualquier instante. Este método fue propuesto inicialmente por Riks y Wempner y posteriormente modificado por Crisfield. la ecuación de equilibrio es: k k k Rk n +1 = λn +1 P − Qn +1 (Un +1 ) = 0 k k Donde λn +1 es el valor de λ en la iteración k del caso de carga n+1 y Un +1 son las deformaciones totales tras la iteración k del caso de carga n+1.2. y P es un vector de fuerzas de referencia. 79 .Introducción al análisis de estructuras con no linealidad geométrica 14 MÉTODO DE LA LONGITUD DEL ARCO Entre los métodos restringidos. La ecuación de equilibrio en un instante cualquiera n+1 del proceso de carga puede ponerse en forma de residuo Rn+1. que define los valores relativos entre las distintas componentes de la fuerza. Sea k una iteración cualquiera (k=1. y al ser no lineal. P s u Figura 21. suponiendo que son independientes de la deformación. En ambos casos se puede emplear el método de Newton normal o el modificado. El planteamiento de Crisfield se desarrolla a continuación. como: Rn +1 = Pn +1 − Qn +1 = 0 Siendo Pn +1 las fuerzas exteriores aplicadas y Qn +1 las fuerzas interiores producidas por las tensiones en los elementos. se resuelve por iteraciones sucesivas. Método de la longitud del arco.…) en la búsqueda del equilibrio para el estado de carga n+1. Consiste en utilizar un círculo de radio Δs con centro en el último equilibrio obtenido. que son funciones no lineales de las deformaciones Un+1. el conocido como método de la longitud del arco es uno de los más usados en la práctica.
se puede desarrollar en serie de Taylor alrededor de su valor en la iteración anterior: −1 ⎞ −1 ⎞ ⎛ ∂R k ⎛ ∂R k k −1 n +1 ⎟ ˆk n +1 ⎟ ˆ k ⎜ ⎜ ⎟ ⎟ Rk ≈ R + + λ ⎜ ⎜ n +1 n +1 ⎟ ⎟U = 0 ⎟ ⎟ ⎜ ⎜ ∂U ⎠ ⎝ ∂λ ⎠ ⎝ ˆk el incremento de la variable λ y U ˆ k el incremento de las deformaciones al Siendo λ efectuarse la iteración k. es decir que se limita el incremento de 80 . De forma similar ΔUk n es el incremento de deformación acumulado tras efectuarse la iteración k. sino que puede mantenerse el del primero. Las derivadas necesarias son: −1 ⎞ ⎛ ∂R k n +1 ⎟ ⎜ ⎟ ⎜ ⎟= P ⎟ ⎜ ⎝ ∂λ ⎠ −1 ⎞ −1 ⎞ ⎛ ∂R k ⎛ ∂Qk n +1 ⎟ n +1 ⎟ ⎜ ⎜ ˆ k −1 ⎟ ⎟ = − = −K ⎜ ⎜ n +1 ⎟ ⎟ ⎟ ⎟ ⎜ ⎜ ⎝ ∂U ⎠ ⎝ ∂U ⎠ Luego la ecuación a resolver en cada iteración es: ˆ k −1 U ˆ k = Rk −1 + P λ ˆk K n +1 n +1 Despejando el incremento de deformación se obtiene: ˆ k = (K ˆ k −1 ) Rk −1 + (K ˆ k −1 ) P λ ˆk = Uk + Uk λ ˆk U n +1 n +1 n +1 n +1 n +1 −1 −1 (89) El primer sumando del incremento de deformación se puede calcular fácilmente y representa la deformación producida por la parte del residuo no equilibrado en la iteración anterior: ˆ k −1 Rk −1 Uk n +1 = (Kn +1 ) n +1 −1 (90) El segundo sumando no puede evaluarse hasta no conocer el valor de λ pero su coeficiente puede evaluarse con sencillez. las deformaciones U y el parámetro λ. y por lo tanto el valor del incremento • Suponiendo por el momento conocido el valor de λ k ˆ en esta iteración. de deformación U Para las deformaciones la actualización es: k −1 ˆk ΔUk +U n = ΔUn (92) −1 Donde ΔUk es el incremento de deformación acumulado a lo largo de las (k-1) iteraciones n anteriores. no es necesario recalcular este término a cada paso de la iteración. y representa la deformación producida por las fuerzas básicas: ˆ k −1 P Uk n +1 = (Kn +1 ) −1 (91) Obsérvese que si se emplea el método de Newton modificado.Introducción al análisis de estructuras con no linealidad geométrica • Considerando que el residuo es una función de dos variables. ˆ . De manera análoga se actualiza el parámetro λ: k k −1 ˆk λn +1 = λn +1 + λ (93) ˆk se efectúa introduciendo una ecuación que imponga la condición de • El cálculo de λ distancia máxima recorrida en este paso de carga. se procede a actualizar los valores de las incógnitas.
ˆk = λ ˆk λ (i ) i = indice(max(cos ϕ1. cos ϕ2 )) 81 . Si se denomina Δs a la distancia máxima a recorrer. Para ello. la condición es: T k (Δs )2 = (ΔUk n ) (ΔUn ) (94) Sustituyendo los incrementos por sus valores y operando: −1 ˆ k )T (ΔUk −1 + U ˆk) (Δs )2 = (ΔUk +U n n −1 k k −1 k ˆk T ˆk + Uk + Uk (Δs )2 = (ΔUk n n +1 + Un +1 λ ) (ΔUn n +1 + Un +1 λ ) −1 k T k −1 k −1 k T k k T k ˆk ˆk 2 + Un + Uk + Un (Δs )2 = (ΔUk +1 ) (ΔUn +1 ) Un +1 λ + (Un +1 ) Un +1 (λ ) n n +1 ) + 2(ΔUn ˆk : La ecuación anterior es una ecuación de segundo grado en λ ˆk )2 + a λ ˆk + a = 0 a1(λ 2 3 T k a1 = (Uk n +1 ) Un +1 T k −1 a2 = 2(ΔUk + Uk n n +1 ) Un +1 T k −1 −1 2 a 3 = (ΔUk + Uk + Uk n n +1 ) (ΔUn n +1 ) − (Δs ) ˆk y λ ˆk . en primer lugar se determina cuál sería el incremento de deformación producido por cada una de las soluciones: k −1 k ˆk ΔUk + Uk n (1) = ΔUn n +1 + Un +1 λ(1) k −1 k ˆk ΔUk + Uk n (2) = ΔUn n +1 + Un +1 λ(2) A continuación se calcula la proyección de dichos incrementos de deformación sobre el incremento de la iteración anterior.Introducción al análisis de estructuras con no linealidad geométrica deformación acumulado en todas las iteraciones efectuadas en este caso de carga. que de alguna manera estima el ‘ángulo’ entre ambos vectores: −1 T k k −1 T k −1 k ˆk ΔUn + Uk (Δs ) cos ϕ1 = (ΔUk n ) ΔUn (1) = (ΔUn ) n +1 + Un +1 λ(1) 2 2 ( ) −1 T k k −1 T k −1 k k ˆk ΔUn + Un (Δs ) cos ϕ2 = (ΔUk n ) ΔUn (2) = (ΔUn ) +1 + Un +1 λ(2) ( ) Se elige aquélla solución que produzca el menor ángulo. De entre ellas se elige aquélla Resolviendo esta ecuación se obtienes dos raíces λ (1) (2) que producirá un incremento de deformación acumulado más próximo al incremento de deformación acumulado en la iteración anterior. es decir el mayor valor del coseno de ϕ. que será un escalar.
Iteración en el método de la longitud del arco. suponiendo que las deformaciones iniciales U0 son nulas y por lo tanto las fuerzas interiores también son nulas. En su lugar es más sencillo definir un valor de λ al comienzo de la iteración. Muchas veces se supone λ1 = 1 . a partir del valor de 1 λ1 supuesto: ⎡ Δs = ⎢ U1 ⎢⎣ 1 ( ) T ⎤ U1 1 ⎥ ⎥⎦ 1/ 2 λ11 (95) 82 .1. que no resulta fácil pues depende del problema estudiado. En la primera iteración del primer caso de carga (n=0. y en base a él determinar el Δs. 14. con lo cual en esta definiendo para ello el valor de λ11 = λ 1 primera iteración se aplica toda la carga básica. ˆ1 .1 Comienzo de la iteración en el primer caso de carga Una pequeña dificultad del método está en la definición del valor de Δs. es: ˆ0 U ˆ 1 = P λ1 K 1 1 1 Al ser conocido λ1 se obtienen los incrementos de desplazamiento iniciales: ˆ 1 = (K ˆ 0 )−1 P λ1 = U1 λ1 U 1 1 1 1 La condición de longitud de arco máximo es: ˆ 1) U ˆ 1 = U1 (Δs )2 = (U 1 T ( ) T 1 1 2 U1 (λ1 ) De esta expresión se puede obtener el Δs a emplear en este caso de carga. La ecuación de equilibrio incremental en esta primera iteración del primer caso. k=1) se toma λ10 = 0 como punto de partida y se define como dato el valor de la carga aplicada en esta primera iteración.Introducción al análisis de estructuras con no linealidad geométrica k n+1 k k-1 n+1 Rn+1 Un+1 k k k k-1 Un+1 k r= s Un+1 n Un Un Un k-1 k Un+1 k-1 Un+1 k Figura 22.
5 F EA0 L0 Y2 3 B C lineal A V/Y2 -2 -1. pero es más flexible a compresión.m.5 -1 -0. El análisis lineal de esta estructura. arroja el siguiente resultado para la relación fuerza – deformación: F lin ⎛Y2 ⎞ ⎟ ⎟ = EA0 ⎜ ⎜ ⎟ ⎜ ⎟ ⎝L ⎠ 0 3 ⎛V lin ⎞ ⎟ ⎜ ⎟ ⎜ ⎟ ⎜ ⎟ ⎝Y ⎠ 2 Se observa que esta respuesta lineal corresponde al primer término de la solución no lineal.Introducción al análisis de estructuras con no linealidad geométrica 1 0.4 0.5 Figura 24. 86 .6 0. El modelo numérico de esta estructura para su simulación mediante el procedimiento nolin está en el archivo modelo1.8 -1 -2.2 -0. además de presentar el fenómeno de la inestabilidad. suponiendo que el estado deformado coincide con el inicial.6 -0. Respuesta de la barra apoyada – deslizante.5 0 0. El modelo no lineal muestra que la estructura es más rígida a tracción que en el modelo lineal.4 -0.8 0.2 0 -0.
87 . El desarrollo es el mismo. disminuye la zona descendente de la curva de respuesta. de tal forma que a medida que se aumenta su rigidez. de constante KM. Se observa que se mantiene la posibilidad del snap-through. salvo que a la fuerza exterior se le debe sumar la fuerza necesaria para deformar el muelle: ⎛Y2 ⎞ ⎟ ⎟ F = EA0 ⎜ ⎜ ⎟ ⎜ ⎟ ⎝L ⎠ 0 3 2 3⎞ ⎛V ⎞ ⎛V ⎞ 3⎛ V 1 ⎟ ⎜ ⎟ ⎟ ⎜ ⎟ + ⎜ ⎟ ⎜ ⎟ + K MV ⎜ ⎟ ⎜Y + 2 ⎜ ⎟ ⎟ ⎟ ⎜ ⎜ ⎟ ⎟ Y 2 Y ⎟ ⎝ ⎠ ⎝ ⎠ ⎜ 2 2 ⎝ 2 ⎠ La rigidez del muelle se puede representar en forma relativa a la rigidez de la barra mediante la constante K : KM EA0 (Y2 )2 =K (L0 )3 Con esto la respuesta del sistema es: 3 2 3⎞ ⎛V ⎛ L0 ⎞ F ⎜ 3⎛ V⎞ 1⎛ V⎞ V ⎟ ⎜ ⎟ ⎟ ⎟ ⎜ ⎜ ⎟ K =⎜ + ⎜ ⎟ + ⎜ ⎟ + ⎜ ⎟ ⎟ ⎜ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠ ⎜Y2 ⎠ ⎜Y2 ⎠ ⎜Y2 2 ⎝ EA0 ⎜ 2⎝ Y2 ⎟ ⎝Y2 ⎠ ⎝ La figura siguiente muestra la respuesta. aunque ahora la respuesta es mucho más suave a consecuencia de la presencia del muelle. pero apoyada en un muelle lineal.Introducción al análisis de estructuras con no linealidad geométrica 15. Figura 25. en función de la rigidez del muelle.2 Ejemplo 2. Barra apoyada elásticamente. Barra deslizante apoyada elásticamente Se estudia la misma barra que en el ejemplo anterior.
15. 88 . Longitud total: 500 cm Módulo de elasticidad: 800. sin la técnica del seguimiento de path y se han encontrado problemas de convergencia.m. horizontal +20 kg. ancho 2 Fuerzas máximas en el extremo: vertical -50 kg.000 kg/cm2.Introducción al análisis de estructuras con no linealidad geométrica Figura 26. y sometida a dos fuerzas puntuales en su extremo. Voladizo muy flexible Se estudia una viga en voladizo vertical. modelizada con elementos viga de dos nudos en formulación co-rotacional. con incrementos fijos de la carga. Respuesta de la barra apoyada elásticamente. Se ha utilizado un método de Newton puro. Número de elementos viga: 15 Dimensiones transversales: canto 2. La figura siguiente muestra el proceso de deformación de la viga. Número de incrementos de carga: 50 El modelo numérico de esta estructura para su simulación mediante el procedimiento nolin está en el archivo modelo5.3 Ejemplo 3.
Las deformaciones finales de este punto son DX= 356 cm y DY=-643 cm.Introducción al análisis de estructuras con no linealidad geométrica 400 300 200 100 0 0 200 400 500 -100 20 50 Figura 27. 89 . La figura siguiente muestra la relación fuerza / desplazamiento para el punto extremo de la viga. Deformación de un voladizo flexible. 50 40 FY 30 20 FX 10 UX . Curvas fuerza – deformación de un voladizo flexible. UY 0 0 100 200 300 400 500 600 700 Figura 28.
con λ1 = 1 . en la que puede verse que presenta un fenómeno de snap-back. 1 Se ha empleado el método del seguimiento del path. En el nudo 1 se aplica una fuerza horizontal.4 Ejemplo 4.5. seguidas por una barra inclinada. El modelo numérico de esta estructura para su simulación mediante el procedimiento nolin está en el archivo modelo6. exponente γ = 5. Figura 29. uno vertical debido a la barra 8-9 y otro horizontal más flexible formado por las 6 barras horizontales. Celosía Este ejemplo corresponde a una celosía muy simple. la estructura se puede considerar formada por una barra. Celosía simple. y el nudo 9 está fijo. 90 3 5 . La carga de referencia aplicada en cada paso de carga es de 40. apoyada en dos muelles. con EA = 3 ⋅ 106 N . número de iteraciones deseado J des = 5 . La figura siguiente muestra la relación entre la fuerza aplicada y la deformación horizontal del nudo 1. y finalmente una barra vertical.m. En cada iteración se ha limitado el incremento de λ a 0. En total hay 8 barras y 9 nudos.000 N. estudiada por varios autores. Con esta disposición. Tiene forma de L. Todas las barras son de las mismas propiedades. la barra inclinada. con 6 barras horizontales una a continuación de la otra.Introducción al análisis de estructuras con no linealidad geométrica 15.
de 120 cm de lado. Tanto el poste como la viga se han modelizado con 10 elementos viga iguales. El valor de referencia aplicado en cada paso de carga es de 100 N. Todas las barras son de las mismas propiedades. La estructura está sometida a una carga puntual vertical situada sobre la viga. estudiado por Lee. Pórtico biarticulado Este ejemplo corresponde a un pórtico biarticulado en L. de 12 cm de longitud. A = 6 cm 2 . dando un total de 20 vigas y 21 nudos. 15. El modelo de esta estructura está en el archivo lee_frame. 91 .5 Ejemplo 5. con E = 720 kN/cm2 . Respuesta con snap-through de una celosía flexible. a 24 cm del poste.Introducción al análisis de estructuras con no linealidad geométrica Nudo: 1 4 x 10 3 2 Fuerza X 1 0 -1 -2 Desplazamiento X 5 -3 0 2 4 6 8 10 12 14 16 Figura 30.m. I = 2 cm 4 .
en la que puede verse una respuesta muy no lineal y un claro fenómeno de snap-back. exponente γ = 5. Pórtico biarticulado. 92 . número de iteraciones deseado J des = 5 . Deformación del pórtico biarticulado flexible. con λ1 = 1 . La figura siguiente muestra la evolución de la estructura en los primeros incrementos de carga. La figura siguiente muestra la relación entre la fuerza aplicada y la deformación vertical del nudo sobre el que se aplica la fuerza (nudo 13). 1 Se ha empleado el método del seguimiento del path.Introducción al análisis de estructuras con no linealidad geométrica P=100 24 10 vigas de 12 cm 120 cm 10 vigas de 12 cm 120 cm Figura 31. Figura 32. En cada iteración se ha limitado el incremento de λ a 2.
Respuesta del pórtico biarticulado flexible. 93 .Introducción al análisis de estructuras con no linealidad geométrica Fuerza Y Figura 33.
apareciendo en el sólido un campo de velocidades u y uno de aceleraciones u que dan lugar a las correspondientes fuerzas de inercia.2 Ecuaciones de equilibrio. Formulación lagrangiana total Considerando un elemento finito.Introducción al análisis de estructuras con no linealidad geométrica 16 DINÁMICA Se considera en este apartado el caso de que las cargas aplicadas sean variables con el tiempo. la hipótesis de discretización permite establecer las relaciones entre los campos de deformaciones y aceleraciones dentro del elemento y los valores nodales de dichas deformaciones y aceleraciones: u = NU u = NU Además se cumple que: δu = N δU Por lo tanto el trabajo virtual de las fuerzas de inercia es: 94 . (tensiones de Cauchy y deformaciones unitarias infinitesimales) o de las magnitudes en el estado inicial (deformaciones de Green-Lagrange y tensiones de Piola-Kirchhoff) δWI = Vt ∫ σ : δε dv t = V0 ∫ S : δE dv 0 16.1 Principio del trabajo virtual en dinámica El principio del trabajo virtual en régimen dinámico indica que la condición necesaria y suficiente para que exista equilibrio es que la suma del trabajo virtual de las fuerzas interiores δWI y el trabajo virtual de las fuerzas de inercia δWIN sea igual al trabajo virtual de las fuerzas exteriores δWE para cualquier variación virtual de las deformaciones δ u . compatible con las condiciones de ligadura: δWI + δWIN = δWE • El trabajo virtual de las fuerzas de inercia en el instante t es: δWIN = −∫ δ uT ρ u dv t vt (97) Esta expresión se puede referir al estado inicial considerando la conservación de la masa: ρ dv t = ρ 0 dv 0 Con lo que queda: δWIN = −∫ δ uT ρ 0 u dv 0 v0 • El trabajo virtual de las fuerzas interiores puede ponerse en función de las magnitudes en el estado t. teniendo sentido la derivada respecto a él. El parámetro t. 16. juega por lo tanto ahora el papel de tiempo. dando lugar a una respuesta dinámica en la que las deformaciones del sólido varíen con el tiempo.
las ecuaciones diferenciales de equilibrio de la estructura completa se obtienen por ensamblado de las ecuaciones de los distintos elementos finitos. cuya contribución a las ecuaciones de equilibrio se representa mediante una matriz de amortiguamiento C. calcular las deformaciones en el instante tn + h aproximando la aceleración y velocidad en tn mediante un operador de diferencias centrales. son de la misma forma que las de un elemento. uno de los métodos más habituales es el método de las diferencias centrales. Para mayor generalidad.3 Método explícito basado en diferencias centrales Para la integración numérica de las ecuaciones de equilibrio. Sustituyendo (101) y (102) en la ecuación de equilibrio en tn y reordenando se obtiene: 95 . se considera la posibilidad de que sobre el sistema existan también efectos de amortiguamiento.Introducción al análisis de estructuras con no linealidad geométrica t δWIN = −∫ (δ ut )T ρ 0 ut dv 0 = −(δ Ut )T ∫ ρ 0 NT Ndv 0 Ut v0 v0 En esta expresión se ha definido la matriz de masas del elemento. pero corresponden a toda la estructura. que es constante y se evalúa en el estado inicial: M= ∫ρ v0 0 NT Ndv 0 El trabajo virtual interior se expresa en función del vector de fuerzas nodales equivalentes a los esfuerzos interiores en el elemento Q: δWI = δ UT Q El trabajo virtual de las fuerzas exteriores define el vector de fuerzas exteriores nodales equivalentes P: δWE = δ UT P (98) Sustituyendo en la expresión del principio del trabajo virtual y considerando que la variación de las deformaciones nodales es arbitraria se obtiene la ecuación de equilibrio del elemento: MU + Q = P (99) Finalmente. La idea es considerar conocido el equilibrio en el instante tn. en la forma: Un = 1 Un −1 − 2 Un + Un +1 h2 ( ) (101) (102) Un = 1 (Un +1 − Un−1 ) 2h Un −1 ≡ Ut −h En estas expresiones se ha introducido la notación: Un ≡ Ut Un +1 ≡ Ut +h El mismo criterio de notación se aplica a las demás magnitudes. y a partir de él. Para un instante de tiempo cualquiera t. que supondremos constante: MU + C U + Q = P (100) 16.
por lo que el método tiene un carácter explícito. pero en dichos pasos se acumula un gran error en la solución. Además esta condición debe satisfacerse en todos los instantes de tiempo durante la simulación. Tiene pues innumerables ventajas que explican su amplia utilización. El proceso de integración es una secuencia de pasos iguales en el tiempo. De hecho puede plantearse el método empleando la ecuación (104) en lugar de la (103) para calcular la deformación en el paso siguiente. sino únicamente el vector de fuerzas interiores. En el caso de que M y C sean diagonales.4 Estabilidad del método de diferencias centrales El principal inconveniente del método explícito de diferencias centrales es que es condicionalmente estable. es necesario a continuación utilizar la ecuación de equilibrio en n+1 para hallar la aceleración en el nuevo estado. no es necesario ensamblar la matriz de rigidez. aplicando de forma repetitiva las ecuaciones anteriores. que distorsiona la solución total obtenida. por lo que el método es explícito. es decir que se debe emplear un tamaño de paso inferior a un paso crítico para que el método sea estable. Esta naturaleza se pone de manifiesto si de las dos ecuaciones (101) y (102) se despeja la deformación en n+1. ni es necesario emplear la matriz de rigidez tangente. El método no requiere ninguna iteración para alcanzar el equilibrio. ni siquiera es necesario resolver ningún sistema de ecuaciones. los requerimientos de almacenamiento de datos son muy pequeños. Desde el punto de vista de la implementación. y si las matrices de M y C son diagonales. obteniéndose: Un +1 = Un + hUn + h2 Un 2 (104) Esta ecuación indica que la deformación en n+1 se puede determinar directamente a partir del estado en n. En esto la respuesta es diferente al análisis estático. no se observa un fenómeno de inestabilidad obvio en la solución total. 16. lo cual suele ser habitual en formulaciones de masas consistentes. 96 . Nótese que la respuesta en tn+h. sin necesidad de aplicar ninguna ecuación de equilibrio. pues casi todas las operaciones se pueden efectuar a nivel de elemento. en el que la respuesta muestra claramente un crecimiento incontrolado si no se satisface el criterio de estabilidad.Introducción al análisis de estructuras con no linealidad geométrica ⎛1 1 ⎞ 2 ⎛1 1 ⎞ ⎜ 2 M+ C⎟ Un +1 = Pn − Qn + 2 M Un − ⎜ M − C⎟ ⎜ 2 ⎟ ⎟ Un −1 ⎜ ⎝h ⎝h 2h ⎠ 2h ⎠ h (103) De esta ecuación se obtiene la deformación en tn+h. que se corresponde con el menor periodo de oscilación Tmin. y a continuación la velocidad y aceleración se obtienen de (102) y (101). Este método es un caso particular de la familia de métodos de Newmark. Dicho paso crítico vale: hCR = 2 ωmax = Tmin π Siendo ωmax la máxima frecuencia propia existente en la malla de elementos finitos. Si la condición de estabilidad no se cumple durante unos pocos de pasos del proceso total. se obtiene apoyándose en el equilibrio en t. aunque si esto se hace así.
y sean ωie las frecuencias propias del elemento e. Por esta razón se trata de obtener estimaciones o límites superiores de dicha frecuencia máxima que sean fáciles de calcular. el problema de autovalores que proporciona las frecuencias propias de un elemento finito de este tipo es: ˆ (K D ˆ ) φ = ω2 M φ +K σ 97 . K φi = ωi2 M φi 2 ωmax = max (ωi2 ) i Sean Ke y Me las matrices de masas y rigidez de los distintos elementos de la malla.i 2 Por aplicación del cociente de Rayleigh se puede demostrar que se cumple siempre que: 2 E 2 ωm ax ≤ (ωmax ) Es decir que la máxima frecuencia individual que presentan los distintos elementos finitos desacoplados unos de otros es mayor que la máxima frecuencia del sistema ensamblado. y sean ω2 las frecuencias propias de dicho sistema. Esto implicar resolver un problema de valores y vectores propios de tamaño igual al número de gados de libertad del sistema. y no la lateral.1 Paso de integración crítico en problemas unidimensionales Consideramos un problema unidimensional. los cuales corresponden al elemento biarticulado ya estudiado. modelizado con elementos de dos nudos. solución de los problemas de autovalores individuales de los distintos elementos: Ke φie = (ωie ) Me φie 2 La mayor de todas estas frecuencias de los distintos elementos desacoplados es: E (ωmax )2 = max (ωie ) e . y para la cual existen de hecho soluciones analíticas.4. lo cual resulta prohibitivo en las aplicaciones reales. cuyo ensamblaje da lugar a las K y M anteriores. cuyo cálculo tiene un costo prohibitivo y sea ωmax la mayor de todas estas frecuencias. resulta del máximo interés determinar el valor de la frecuencia máxima presente en la malla de elementos. Esto proporciona un límite superior de la frecuencia máxima del sistema ωmax que es muy fácil de evaluar. Figura 34. pero considerando sólo la deformación axial.Introducción al análisis de estructuras con no linealidad geométrica A la vista del paso crítico para garantizar la estabilidad de los métodos. Sean K y M las matrices de rigidez y masas del sistema estructural estudiado. 16. Empleando una formulación lagrangiana total. Elemento unidimensional.
La tabla siguiente muestra los valores más habituales. Elemento Unidimensional de 2 nudos Unidimensional de 2 nudos Matriz M Diagonal Consistente ωe max 2c L 2 3c L hCR L c L 3c 98 . quien la formuló para modelos de diferencias finitas. despreciando por una parte la variación de la longitud y por otra la tensión S frente al módulo de elasticidad E: ⎛ E (L / L0 )2 + S ⎞ ⎟ ⎟ c0 = ⎜ ⎜ ⎟ ⎟ ⎜ ρ0 ⎝ ⎠ 2 ωmax 1/ 2 ≈ E ρ0 Por lo tanto el paso crítico para un modelo unidimensional con este tipo de elemento es: hCR = = L0 c0 Esta expresión se suele denominar condición de Courant. En cada caso se obtiene el valor de la frecuencia máxima del elemento.2 Pasos críticos de integración para diversos elementos finitos Se puede efectuar un análisis similar al efectuado para el elemento unidimensional.4. aunque habitualmente se simplifica.Introducción al análisis de estructuras con no linealidad geométrica Considerando sólo los términos correspondientes a la deformación en la dirección del elemento y empleando la matriz de masas diagonal. para otros tipos de elementos finitos sencillos en los que se conozca la expresión analítica de sus matrices de rigidez y masas. la ecuación anterior es: ⎛ ⎛ L ⎞2 ⎞ ⎧ ⎫ ⎪φ1X ⎪ ⎪ A0 ⎜ ⎟ ⎢⎡ 1 −1⎥⎤ ⎪ ⎜ ⎜ 2 ρ0A0L0 ⎟ ⎟ ⎟ E S + = ω ⎜ ⎟ ⎟ ⎨ ⎬ ⎜ ⎢ ⎥ ⎟ ⎟ φ ⎪ L0 ⎜ 2 ⎟ ⎣⎢−1 1 ⎥⎦ ⎪ ⎝ L0 ⎠ ⎝ ⎜ ⎠ ⎪ ⎪ ⎩ 2X ⎭ ⎫ ⎡ 1 0⎤ ⎧ ⎪φ1X ⎪ ⎪ ⎪ ⎢ ⎥⎨ ⎢ 0 1⎥ ⎪φ ⎬ ⎪ 2X ⎪ ⎪ ⎣⎢ ⎦⎥ ⎩ ⎭ Esta ecuación tiene dos soluciones. 16. La segunda corresponde a la frecuencia máxima del elemento y su valor resulta ser: ωmax = 2c0 L0 1/ 2 Siendo c0 la velocidad de propagación de las ondas elásticas en el material: ⎛ E (L / L0 )2 + S ⎞ ⎟ ⎟ c0 = ⎜ ⎜ ⎟ ⎟ ⎜ ρ0 ⎝ ⎠ En principio esta velocidad depende del nivel de tensión y de la longitud deformada. la primera es ω=0. Esta condición lo que especifica es que el paso de integración debe ser como mínimo aquel tiempo que permita la propagación de una onda elástica de velocidad c0 dentro del elemento de longitud L0. que condiciona el paso crítico del método de las diferencias centrales. que no interesa.
Se debe comprobar además el valor correspondiente a la deformación axial. 16. que es condicionalmente estable. lo cual implica la realización de un proceso iterativo para hallar la solución.Introducción al análisis de estructuras con no linealidad geométrica Unidimensional de 3 nudos Viga plana de 2 nudos (1) Diagonal Diagonal Diagonal 2 6c L L 6c A L2 48I c L 1−ν c Cuadrado plano. como en un elemento de 2 nudos. se caracteriza por calcular los desplazamientos y velocidades en el instante tn+h apoyándose en el estado conocido anterior. La familia de Newmark. β=1/4 y γ=1/2 corresponde a una aceleración media constante en el intervalo y es el método originalmente propuesto por Newmark.5 Métodos implícitos de integración de paso simple Todos los métodos empleados para la integración numérica de ecuaciones diferenciales de segundo orden pueden emplearse para la resolución de problemas no lineales. De entre ellos. Si se emplea β=1/6 y γ=1/2 se obtiene un método con interpolación lineal de las aceleraciones. En todos ellos se plantea el equilibrio en el instante tn+h. que es incondicionalmente estable. mediante un desarrollo en serie de los mismos hasta términos de orden 2. la familia de métodos de Newmark o el método de Wilson son unos de las más populares. en la forma: Un +1 = Un + ∫ U(τ )d τ 0 h Un +1 = Un + h Un + ∫ U(τ )(h − τ )d τ 0 h Las integrales se evalúan mediante una regla de cuadratura. Así. se obtienen diferentes métodos. con lo que las aproximaciones de posición y velocidad son Un +1 = Un + (1 − γ ) h Un + γ h Un +1 Un +1 = Un + h Un + h2 2 ⎡(1 − 2β ) U + 2β U ⎤ n n +1 ⎥ ⎢⎣ ⎦ (105) (106) Adoptando diferentes valores de los parámetros β y γ. (1) Corresponde sólo al efecto de flexión. Despejando la aceleración de (106) se obtiene su valor en función de las deformaciones: Un +1 = ⎛1 ⎞ 1 1 ⎟U ⎜ − 1⎟ U − Un ) − Un − ⎜ n 2 ( n +1 ⎟ ⎜ 2 βh βh β ⎝ ⎠ (107) Sustituyendo este valor en (105) se obtiene la velocidad en función de las deformaciones: 99 . Tensión plana.
Sea k=1.. apoyándose en la solución conocida en el instante anterior n. por lo que debe emplearse un proceso iterativo para obtener la respuesta en el instante n+1.. Así. al final de la iteración k. • El valor de la matriz de rigidez efectiva se obtiene derivando (109) (particulariza para k-1 en vez de k) y es: Kn +1 = k −1 −1 −1 k −1 k −1 ∂R k ∂Uk ∂Un ∂Qn n +1 n +1 +1 +1 M C = + + k −1 k −1 k −1 k −1 ∂Un +1 ∂Un +1 ∂Un +1 ∂Un +1 Se ha supuesto que las fuerzas exteriores P no dependen de las deformaciones. derivando (107) respecto de las deformaciones Un+1 se obtiene: 100 . la secuencia de iteraciones. en una cualquiera de ellas las ecuaciones de equilibrio se pueden poner en forma de residuo: k k k Rk n +1 ≡ M Un +1 + C Un +1 + Qn +1 − P n +1 = 0 (109) • Desarrollando este residuo en serie alrededor del punto de iteración anterior (k-1) se obtiene: −1 ⎞ ⎛ ∂R k −1 k −1 n +1 ⎟ ⎜ ⎟ Uk − Uk Rk ≈ R + ⎜ n +1 n +1 n +1 ) = 0 k −1 ⎟ ( n +1 ⎜ ⎟ ⎜ ⎝ ∂Un +1 ⎠ La derivada del residuo respecto de las deformaciones define la matriz de rigidez efectiva del sistema: −1 ∂Rk k −1 n +1 = Kn +1 k −1 ∂Un +1 Por lo tanto la ecuación a resolver en cada iteración es: k −1 ˆ k = −Rk −1 Kn +1 U n +1 (110) ˆ k es el incremento de deformación producido en la iteración k: En esta expresión U ˆ k = Uk − Uk −1 U n +1 n +1 (111) y Uk n +1 es la estimación de las deformaciones en el instante n+1.2.Introducción al análisis de estructuras con no linealidad geométrica Un +1 = Un + (1 − γ ) h Un + Un +1 = γ γ βh (Un +1 − Un )− Un − γh ⎜ ⎜ ⎜ ⎛γ ⎞ ⎛γ ⎞ γ β ⎛1 ⎞ − 1⎟ ⎟ ⎟ Un ⎝ 2β ⎠ (108) ⎜ − 1⎟ − 1⎟ ⎟ ⎟ (Un +1 − Un ) − ⎜ ⎜ ⎜ ⎟ Un − h ⎝ ⎟ Un ⎜β ⎜ 2β βh ⎝ ⎠ ⎠ Las ecuaciones anteriores deben combinarse con la ecuación de equilibrio dinámico del sistema. que es no lineal. El último sumando introduce la matriz de rigidez tangente: Kn +1 = M k −1 −1 −1 ∂Uk ∂Uk n +1 n +1 ˆ k −1 C + +K n +1 k −1 k −1 ∂Un +1 ∂Un +1 Las derivadas de las deformaciones que aparecen en esta expresión se pueden obtener derivando las aproximaciones del método de integración.
Un . Un +1 ) (117) Este sistema de ecuaciones lineales se emplea en un proceso iterativo k=1.Introducción al análisis de estructuras con no linealidad geométrica −1 ∂Uk 1 n +1 = k −1 ∂Un +1 βh 2 ⎡ ∂Uk −1 ⎢ n +1 − ∂Un k −1 ⎢ ∂Uk −1 ∂Un +1 ⎢⎣ n +1 ⎤ ⎞ ∂Un 1 1 ⎥ − 1 ∂Un − ⎛ ⎜ I = − 1⎟ ⎟ ⎜ − 1 k ⎥ βh ∂Uk −1 ⎟ ⎜ βh 2 ⎝ 2β ⎠ ∂Un +1 n +1 ⎥⎦ Derivando (108) se obtiene: −1 −1 ⎛γ ⎞ ∂Un ⎛γ ⎞ ∂Uk γ ∂Uk γ ∂Un n +1 n +1 ⎟ ∂Un = γ I ⎜ − 1⎟ = − −⎜ − 1⎟ −h⎜ ⎟ ⎜ − 1 −1 k −1 k −1 k −1 k ⎟ ⎟ ⎜ ∂Un +1 βh ∂Un +1 βh ∂Un +1 ⎜ βh ⎝β ⎠ ∂Un +1 ⎝ 2β ⎠ ∂Uk n +1 Por lo tanto el valor final de la matriz de rigidez efectiva es: Kn +1 = k −1 1 γ ˆ k −1 M+ C+K n +1 2 βh βh (112) • El término independiente de la ecuación (110) es el residuo al final de la iteración. y éste a continuación en la ecuación de la iteración (110) se obtiene: ⎛1 ⎞ k −1 ⎟M U ˆ k = P − Qk −1 − 1 M (Uk −1 − U ) + 1 M U + ⎜ ⎜ − 1⎟ Kn +1 U n +1 n +1 n +1 n n n ⎟ ⎜ 2β βh 2 βh ⎝ ⎠ ⎛γ ⎞ ⎛γ ⎞ γ −1 ⎜ − 1⎟ C (Uk C Un + h ⎜ − − 1⎟ ⎟ ⎟ ⎜ n +1 − Un ) + ⎜ ⎟ ⎟ C Un ⎜β ⎜ 2β βh ⎝ ⎠ ⎝ ⎠ (116) Todos los sumandos del término independiente son conocidos. para k=1. Un .… hasta alcanzar la convergencia. agrupando todos los sumandos del término independiente en un vector de fuerzas efectivas. por lo que ese sistema de ecuaciones lineales puede resolverse y proporciona el valor del incremento de deformación a aplicar en la iteración k. que depende del estado anterior y de última estimación de las deformaciones en el estado actual: k −1 −1 k −1 ˆ k = Pk Kn +1 U n +1(Un . cambiado de signo: −1 k −1 k −1 k −1 Rk n +1 = M Un +1 + C Un +1 + Qn +1 − P n +1 (113) La expresión de la aceleración en la iteración k-1 está dada por (107): −1 Uk n +1 = ⎛1 ⎞ 1 1 ⎟U ⎜ − 1⎟ Uk −1 − Un ) − Un − ⎜ n 2 ( n +1 ⎟ ⎜ βh βh ⎝ 2β ⎠ (114) El valor de la velocidad en la iteración k-1 viene dado por (108): −1 Uk n +1 = γ βh k −1 ⎜ (Un +1 − Un ) − ⎜ ⎜ ⎛γ ⎞ ⎞ γ ⎟U − h ⎛ ⎟U ⎜ − 1⎟ − 1⎟ ⎜ n n ⎟ ⎟ ⎜ ⎝β ⎠ ⎝ 2β ⎠ (115) Sustituyendo los valores de la aceleración (114) y la velocidad (115) en la expresión del residuo (113). Se puede poner en forma compacta.2. Como punto de partida del mismo. se emplea: 0 Un +1 = Un Se observa que la ecuación a resolver para efectuar la integración de las ecuaciones diferenciales de equilibrio dinámico es del mismo tipo que la ecuación obtenida en el análisis 101 . bien del último paso de integración n o de la iteración anterior.
En el caso de emplear un criterio basado en el residuo. Por lo tanto todos los métodos y estrategias de iteración empleados en el análisis estático para este tipo de sistemas de ecuaciones son aplicables en este caso dinámico. 16. en el caso de emplear un criterio basado en el incremento de energía interna en cada iteración se deben añadir los términos correspondientes a dichas fuerzas.Introducción al análisis de estructuras con no linealidad geométrica estático. debe considerarse en él los términos correspondientes a las fuerzas de inercia y amortiguamiento: Rk n +1 2 0 ≤ εR Rn +1 2 k k Pn +1 − Qk n +1 − M Un +1 − C Un +1 2 ≤ εR Pn +1 − Qn − M Un − C Un 2 De la misma forma. ˆ ) (P (U k T n +1 k k k ˆ 1 P − Q − MU − C U − Qn n +1 n n n +1 − MUn +1 − C Un +1 ≤ εE (U ) T ) ( ) 102 . excepto por los valores de la matriz y el vector de fuerzas efectivas. que ahora incluyen términos debidos a la inercia y al amortiguamiento.6 Criterios de convergencia Para finalizar la iteración de búsqueda del equilibrio se pueden emplear los mismos tipos de criterios de convergencia empleados en el caso estático.
hasta alcanzar el equilibrio con una deformación final de 43. la cual coincide con el valor hallado estáticamente.1 106 kg/cm2 ρ=7860 kg/cm3 La matriz de masa se ha supuesto diagonal.1. aplicada en forma escalón en t=0.Introducción al análisis de estructuras con no linealidad geométrica 17 EJEMPLOS DINÁMICOS 17. El paso mínimo para garantizar la estabilidad es hcr ≥ L L = = 0. Al ser la carga aplicada superior al valor que provoca el snap-through.32 cm. La figura siguiente muestra la evolución de la deformación vertical del punto de aplicación de la carga. con un factor de proporcionalidad de valor 5. Barra apoyada . que corresponde a una frecuencia máxima en la estructura de valor 483 rad/s.013 s. Barra apoyada – deslizante Se trata de la misma barra estudiada en el apartado 15. En la respuesta se observa un periodo de oscilación del orden de 0. la cual produce el paso de integración crítico de valor 2/481=0. Se ha efectuado una simulación dinámica empleando un integrador explícito basado en diferencias centrales.004 s c (E / ρ) El modelo numérico para su simulación mediante el procedimiento dynex se encuentra en el archivo modelo1D.0041 s 103 . El amortiguamiento se supone proporcional a la matriz de masas.1 Ejemplo 1. Las propiedades de la barra son: A0: 2 cm2 E: 2. pero sometida a una carga dinámica. La configuración geométrica particular estudiada se muestra en la figura siguiente: Figura 35. se observa un salto brusco en la fase inicial de la respuesta.m. con paso de integración h=0. es decir: C = 5M Se aplica una carga exterior vertical de valor -800 kg.deslizante.002 s.
El paso mínimo para garantizar la estabilidad es hCR=6. Se han efectuado simulaciones dinámicas empleando: • Un integrador explícito basado en diferencias centrales (procedimiento dynex). La matriz de masa se ha supuesto diagonal. 17. β=1/4 (procedimiento dynim). El modelo numérico de esta estructura para su simulación mediante los procedimientos dynex y dynim está en el archivo modelo5D.1 0.9 1 Figura 36.8 0. La viga está modelizada mediante elementos viga de dos nudos en formulación co-rotacional.deslizante. Deformación vertical de la barra apoyada .2 0. con un factor de proporcionalidad de valor 5.m. Voladizo muy flexible.4 0. La densidad empleada es ρ=2700 kg/cm3 . Los resultados obtenidos en ambos casos son prácticamente iguales.18 10-5 s. Se estudia la respuesta dinámica de la viga en voladizo vertical ya analizada en un ejemplo anterior en régimen estático (apartado 15.3 0.Nudo 2 -10 -20 -30 -40 -50 Tiempo -60 0 0. • Un integrador implícito de Newmark con γ=1/2. 104 . El amortiguamiento se supone proporcional a la matriz de masas. es decir: C = 5M Se aplican las mismas fuerzas que en caso estático (FY=-50 kg y FX=20 kg. con paso de integración h=1 10-3 s.3).7 0. en el extremo superior de la viga) en forma escalón en t=0.5 0.2 Ejemplo 2.6 0. con paso de integración h=6 10-5 s. y sus propiedades son las mismas que en el análisis estático.Introducción al análisis de estructuras con no linealidad geométrica 0 Deformación Y .
105 . Evolución dinámica del voladizo vertical. que lógicamente coincide con la obtenida en el análisis estático. y la estructura adopte una configuración deformada final estática. Figura 37.Introducción al análisis de estructuras con no linealidad geométrica La figura siguiente muestra el proceso de deformación de la viga. La figura siguiente muestra la evolución en el tiempo de la deformación horizontal del extremo superior de la viga. La presencia de amortiguamiento hace que con el paso del tiempo la velocidad y aceleración se anulen.
Se emplea la matriz de masa diagonal y no se considera el amortiguamiento. 106 . según la ley: Ftot = 10000 t (kg / s ) . Se aplica una carga exterior transversal al cable.3 Ejemplo 3. y se modeliza mediante un total de 20 barras biarticuladas. 17. Figura 39. Cable pretensado Se estudia un cable pretensado y sometido a una carga transversal distribuida que varía linealmente con el tiempo. distribuida en toda su longitud. con una resultante total que varía linealmente con el tiempo. Las propiedades del cable son: A0: 2 cm2 E: 2 ⋅ 106 kg/cm2 ρ=7860 kg/cm3 La fuerza de pretensión es N0 = 2000 kg. Cable pretensado.Introducción al análisis de estructuras con no linealidad geométrica Figura 38. Deformación dinámica del extremo superior de la viga. El cable tiene una luz de 20 m.
Deformación vertical del punto central del cable pretensado. 100 90 80 70 60 50 40 30 20 10 Tiempo (s) 0. requiriendo un número medio 2 de iteraciones por cada paso. que muestra un comportamiento muy no lineal desde los primeros instantes del movimiento (la respuesta lineal es una cúbica. tras 0.7 s de integración. con paso de integración h=1 10-3 s. El paso mínimo para garantizar la estabilidad es hCR=2. β=1/4 (procedimiento dynim).4 0. 107 . con paso de integración h=1. El modelo numérico para la simulación se encuentra en el archivo cableD.m.1 0. mostrada a efectos comparativos). 10-4 s. Con ambos integradores los resultados son coincidentes. observándose una diferencia en la posición del orden del 0.6 0.15%.5 0.5 10-4 s. Se han efectuado dos simulaciones dinámicas empleando: • Un integrador explícito basado en diferencias centrales (procedimiento dynex). La Figura 40 muestra la evolución de la deformación vertical del punto central del cable.7 Deformación vertical (cm) Nudo central Lineal Figura 40.3 0. concentrado en cada uno de ellos la parte de cable que le corresponde. La resolución del sistema de ecuaciones no lineales en cada paso de integración se efectúa por el método de Newton. esta carga total se ha aplicado sobre los nudos.Introducción al análisis de estructuras con no linealidad geométrica En la implementación del modelo. • Un integrador implícito de Newmark con γ=1/2.2 0.
108 .. 1991. Prentice-Hall.. Oxford University Press. Glasgow. 2004. R. Finite Element Procedures.. 2004. Vol. (ed. A.. M.. Elsevier Science Publishers. Colorado. 1991. J. 1983. Non linear Finite Element Methods. Taylor R. NAFEMS. 1... An Introduction to Nonlinear Finite Element Analysis. Gruttmann F. Introduction to Nonlinear Finite Element Analysis. Prentice-Hall.. 20.... Idelshon S. 2001. Newmark N.. A new arc-length method for handling sharp snap-backs. C. USA. 2. A.. Schoop H. Wagner W. A robust non-linear mixed hybrid quadrilateral shell element. Hinton E. Zienkiewicz O. J. Computers and Structures.Introducción al análisis de estructuras con no linealidad geométrica 18 BIBLIOGRAFÍA Asghar Bhatti M. Bathe K. 2006. Communications in Numerical Methods in Engineering. Moran B. Vol. 705-709. Cambridge University Press. Non-linear Finite Element Analysis of Solids and Structures. Geradin M. 1992. 67-94.. 2ª Ed. A non-linear triangular curved shell element. September 2005. 2000. 2008. 1959. C. 1991.. Crisfield M. N. Crisfield M. Wood R.. University of Liège. The Finite Element Method.. Reddy J. Wiley. 417-471. Bonet J. Belytschko y T. J. T. 251-264. Nonlinear Continuum Mechanics for Finite Element Analysis. McGraw-Hill. 66. Int. Ed. Non-linear Finite Element Analysis of Solids and Structures – Vol. 2 Advanced Topics. 1965. Foundations of Solid Mechanics. Advanced Topics in Finite Element Analysis of Structures. Crisfield M. J. B. 1996. Nº 5. 85. Felippa C. J.. Wiley.).. Belytschko T. UK. Vol. University of Colorado.. Nonlinear Finite Elements for Continua and Structures. D. K. Fung Y. 1994.. Rixen D. J. 1998. A.. Geradin M.. Wenzel T. Vol. Wiley. A Method of Computation for Structural Dynamics. Kotronis P. Collin F. Wiley. Implicit Finite Element Methods. Vol. ASCE Journal of Engineering Mechanics Division.. Computational Methods for Transient Analysis.. Vol. Hughes. 2005. Boulder. Hogge M. Wiley.. Liu W. 635-666. J. 64.. Journal on Numerical Methods in Engineering.. Internal report Géomac/3S. L. Mechanical Vibrations. Department of Aerospace Engineering Sciences. Hellweg H. A.. Implementation of path following techniques into the finite element code Lagamine.
Por otra parte cada una de de ellas tiene ventajas e inconvenientes respecto a ser más o menos compactas. que en ocasiones son diferentes y otras veces coincidentes.Introducción al análisis de estructuras con no linealidad geométrica 19 ANEJO 1. Sin embargo existen muchísimas excepciones. como matricial: el tipo de representación quedará definido por el contexto y por los operadores empleados. En ambos casos se emplea la negrilla para vectores y matrices. y en muchas ocasiones las expresiones obtenidas son casi iguales. y los tensores de orden 2 como una matriz de 2 dimensiones. intuitivas.1. el ‫ ׃‬para el producto contracto (contracción de dos índices) y el ⊗ para el producto tensorial.3 Notación de matrices Es la más habitual en textos de ingeniería mecánica y de estructuras por su equilibrio entre claridad. En esta notación suele ser habitual denominar a los tensores de orden 1 con letras negrillas minúsculas. se empleará la misma letra o símbolo para denominar a una misma magnitud tanto en su representación tensorial. 109 . 19. 19. en cuyo caso se añade a la representación vectorial una barra sobre el símbolo. para no complicar la notación. Se empleará sólo cuando sea necesario. las expresiones obtenidas son válidas en sistemas de coordenadas no cartesianos. y es la que se empleará preferentemente. NOTACIÓN La mayor parte de las magnitudes empleadas en mecánica de sólidos tienen carácter tensorial.1 Notación de subíndices Es muy empleada en los textos de mecánica de los medios continuos. Corresponde a una representación directa de la notación de tensores. con lo que las expresiones son mucho más compactas y fáciles de recordar. fácil implementación y compacidad similar a la notación de tensores estricta. En algunos casos puede producirse alguna confusión entre la representación tensorial y matricial.1. En este caso los subíndices no se muestran explícitamente. que normalmente son aplicables sólo en coordenadas cartesianas. Su principal inconveniente es que da lugar a expresiones muy farragosas.1. Tiene las ventajas de su generalidad y la facilidad de transformarse en algoritmos implementables en lenguajes de programación. y a los de orden 2 o superior con negrilla en mayúsculas. y para su manejo existen distintas notaciones. Además. En esta notación se introducen símbolos específicos para las operaciones entre tensores: el ⋅ para el producto escalar (contracción de un índice). o fáciles de transformar en algoritmos. al ser las magnitudes tensoriales independientes del sistema de referencia. Como es habitual se supone que los tensores de orden 1 (vectores) se representan como una matriz de una columna. unido al hecho de la escasa formación en su utilización más allá de los casos simples. 19.2 Notación de tensores Es muy utilizada asimismo en textos de mecánica de los medios continuos. donde el primer índice corresponde a la fila. en los cuales la nomenclatura de subíndices es imprescindible. Además.
Notación de subíndices: Dij = ai bj Notación de tensores: D = a ⊗ b Notación de matrices: D = a bT 20. Se emplea la misma expresión en notación de tensores y notación de matrices: k j . 20. -1 si la permutación es impar y 0 si hay índices repetidos.2 Operaciones entre tensores de orden 2 • Producto ordinario. • Producto contracto o producto escalar de dos tensores de orden 2.1 Resumen de álgebra de vectores y tensores A continuación se resume la notación empleada para las operaciones más importantes.Introducción al análisis de estructuras con no linealidad geométrica 20 ANEJO 2. pero es más claro ponerlo para indicar que se contrae un índice entre ambos tensores.1 Operaciones entre vectores • Producto escalar de vectores.j. Notación de subíndices: ci = ∑ eijk ai bj El símbolo e representa el tensor alternador de orden 3.1. Notación de subíndices: s = ∑ ai bi Notación de tensores: s = a ⋅ b = b ⋅ a Notación de matrices: s = aT b = bT a • Producto vectorial de vectores. Produce un tensor D de orden 2. PRELIMINARES MATEMÁTICOS 20. j 110 . para dar lugar a un escalar.1.k i s = A : B = ∑ Aij Bij = B : A i. o producto diádico. Notación de tensores: c = a × b Notación de matrices: c = a b La notación a corresponde a la matriz hemisimétrica asociada el vector a. muchas veces se omite el símbolo ⋅ entre los tensores. o composición de tensores. Contrae dos índices. Notación de tensores: D = A ⋅ B Notación de matrices: D = A B Notación de subíndices: Dij = ∑ Aik Bkj En la notación de tensores. • Producto tensorial de vectores. Da lugar a otro tensor del mismo orden. Se emplean minúsculas para los tensores de orden 1 y mayúsculas para los de orden 2.k} es par. definido como eijk = 1 si la permutación {i.
20. el producto contracto se puede poner como: A : B = ∑ Aij Bij = tr(AT B) = tr(B AT ) = tr(A BT ) = tr(BT A) i. es un vector definido por las tres derivadas parciales del escalar respecto a las tres coordenadas del espacio. muchas veces se omite el símbolo ⋅ entre los tensores.3 Gradiente • El gradiente de un campo escalar f.2 Traza Para un tensor de orden 2. pero es más claro ponerlo para indicar que se contrae un índice. Notación de tensores: c = A ⋅ b Notación de matrices: c = A b En la notación de tensores. se define como el escalar: s = tr(A) = A : I = ∑ Aii i Sus propiedades más importantes son: tr(A) = tr(AT ) tr(A B) = tr(B A) tr(a ⊗ b) = tr(a bT ) = ∑ (a ⊗ b) i ii = ∑ aibi = a ⋅ b = aT b i Empleando la traza. Pueden emplearse las notaciones siguientes: ⎧ ⎪ ⎪ ∂f ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ∂ ⎪ 1⎪ ⎪ ⎪ ∂f ⎪ ⎪ ∂f ⎪ ∇f = grad( f ) = =⎪ ⎨ ⎬ ∂x ⎪ ∂x 2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ f ∂ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ⎪ ⎪ ∂ 3 ⎪ ⎪ ⎩ ⎭ Puede pensarse en el gradiente ∇ como un operador vectorial (no se añade la barra pues no ha lugar a confusión): 111 .3 Operaciones entre vectores y tensores • Producto ordinario.Introducción al análisis de estructuras con no linealidad geométrica 20. j 20.1. de un tensor de orden 2 por un vector.
Pueden emplearse varias notaciones: d = div(v)=tr(∇v) = ∇v : I = ∇ ⋅ v d = div(v)=∇T v 112 .k ∂D ∂x Tijk = ∂Dij ∂x k ∂Dij ∂x k ei ⊗ e j ⊗ ek 20. j ∂ vi ei ⊗ e j ∂x j • El gradiente de un tensor D de orden 2 es otro tensor T. es decir a la traza del gradiente del vector. T=∇D = grad(D) = T= = ∑ i . que se obtiene aplicando el operador gradiente a cada una de las componentes escalares del vector: grad(v) = ∇v = ∂v ∂x ∂v i ∂x j En notación de subíndices y matrices se representa como: (grad(v))ij = (∇v)ij = grad(v)= ∇ vT ( ) T Si consideramos tres vectores ei que definen una base ortogonal para las coordenadas x. su expresión es: grad(v)= ∇ ⊗ v ( ) T =∑ i. de orden 3. cada uno de cuyos términos es la derivada de las componentes del tensor respecto de las tres coordenadas. j .4 Divergencia • La divergencia de un vector v es un escalar d definido por: d=∑ i ∂vi ∂x i Esta suma de las tres derivadas direccionales de las componentes de vector corresponde a la suma de los términos de la diagonal del gradiente del vector.Introducción al análisis de estructuras con no linealidad geométrica ⎧ ⎪ ⎪ ∂ ⎫ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ x ∂ ⎪ 1⎪ ⎪ ⎪ ⎪ ⎪ ∂ ⎪ ⎪ ⎪ ⎪ ∇=⎨ ⎬ ⎪ ∂x 2 ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ∂x 3 ⎪ ⎪ ⎩ ⎭ • El gradiente de un campo vectorial v es un tensor de orden 2.
T div(D a)=∑ ai (div(DT ))i + ∑ ∑ Dij (grad(a))ij j j i ( ) En el primer sumando se identifica el producto escalar del vector a por la divergencia de DT. Se trata de hallar la divergencia del producto de ambos b=D a (que será un escalar): div(D a)=div(b)=∑ j ∂b j ∂x j ⎞ ⎛ ∂a i ⎞ ⎟ ⎟ ⎟ + ∑ ∑⎜ ⎟ ⎜ ai ⎟ D ⎟ ji ⎜ ⎟ ⎟ ⎜ x ∂ ⎠ j i ⎝ j ⎠ Sustituyendo el valor del término de bj: ⎛ ∂ ⎜ ⎜ div(D a)=∑ ⎜ ⎝ ∂x j j ⎜ ⎟ ⎜ ⎟= ⎜ ∑D a ⎠ ⎟ ⎜ ⎟ ∑ ∑⎝ ⎜ ∂x ji i i j i ⎞ ⎛ ∂Dji j Reordenando y empleando la traspuesta de D: T ⎞ ⎛ ∂Dij ⎛ ⎞ ⎟ T ∂a i ⎟ ⎜ ⎟ ⎟ ⎜ ⎜ div(D a)=∑ ai ∑ ⎜ D + ∑ ∑ ⎟ ⎟ ij ⎜ ⎜ ⎟ ⎟ ⎜ ∂x j ⎠ ⎜ ∂x j ⎠ j i ⎝ j i ⎝ En el primer término se identifica el término i de la divergencia de DT.5 Teoremas de integración • La integral a un volumen V del vector gradiente de una función escalar f es igual al flujo de dicho escalar en la superficie ∂V que rodea al volumen (n es el vector normal a la superficie): 113 . y en el segundo el producto contracto del ardiente de a por DT: div(D a)=a ⋅ div(DT ) + DT : grad(a) 20. Sea D un tensor de orden 2 y a un vector. Descomponiendo el gradiente en sus componentes simétrica y hemisimétrica ∇v = (∇v)sim + (∇v)hem div(v) = tr(∇v) = tr(∇v)sim + tr(∇v)hem = tr(∇v)sim + 0 • La divergencia de un tensor D de orden 2 es un vector d definido como la traza del gradiente de dicho tensor d=div(D)=tr(∇D) = (∇D) : I = ∇ ⋅ D Puesto en forma matricial su expresión es: d = div(D) = ∇T DT Sus componentes son: ( ) T {div(D)}i = ∑ j ∂Dij ∂x j • A continuación se desarrolla una propiedad de la divergencia que resulta útil para ciertos desarrollos teóricos.j del gradiente del vector a. y en el segundo se identifica el término i.Introducción al análisis de estructuras con no linealidad geométrica La divergencia afecta a la parte simétrica del gradiente del vector.
que establece que la divergencia del vector en un volumen V es igual al flujo de dicho vector en la superficie circundante de V. la expresión es: ∫ div(D) dV = ∫ D ⋅ n dS V ∂V ∫ ∇ ⋅ D dV = ∫ D ⋅ n dS V ∂V 114 .Introducción al análisis de estructuras con no linealidad geométrica ∫ ∇f dV = ∫ f n dS V ∂V • Aplicando la ecuación anterior a un vector v. • Para el caso de un tensor D de orden 2. ∫ div(v) dV = ∫ v ⋅ n dS V ∂V ∫ ∇ ⋅ v dV = ∫ v ⋅ n dS V ∂V Estas expresiones constituyen el teorema de la divergencia para un vector cualquiera. y el de la derecha es el flujo de dicho vector. Tomando la traza del tensor se obtiene una igualdad escalar: ∫ tr(∇v)dV = ∫ tr(v ⊗ n)dS V ∂V El término de la izquierda es la divergencia del vector v. componente a componente: ∫ ∇v dV = ∫ v ⊗ n dS V ∂V Cada integral proporciona un tensor.
PROCEDIMIENTOS MATLAB Procedimiento nolin Este procedimiento permite efectuar la simulación estática no lineal de estructuras planas compuestas por barras biarticuladas o vigas planas (no se pueden mezclar). 115 .Introducción al análisis de estructuras con no linealidad geométrica 21 ANEJO 3. Implementa el método de la longitud del arco y la iteración de Newton. Procedimiento dynim Este procedimiento permite efectuar la simulación dinámica no lineal de estructuras planas compuestas por barras biarticuladas o vigas planas (no se pueden mezclar). completo o modificado. Implementa el método explícito de diferencias centrales. Procedimiento dynex Este procedimiento permite efectuar la simulación dinámica no lineal de estructuras planas compuestas por barras biarticuladas o vigas planas (no se pueden mezclar). con paso fijo. Implementa el método implícito de Newmark. El proceso iterativo para alcanzar el equilibrio en cada paso de integración en el tiempo se efectúa mediante el método de Newton. completa o modificada. con paso fijo.
Documents Similar To No Lineal Continuo
Elasticidad - Luis Ortiz Berrocal
Calculo_Vectorial Usco 2015 A
fir 315 clase 4 u3 2011
Filtros-Traducido
Ensayo Triaxial.docx 11
07-Bordes
reportedemecanicadesuelostriaxial-140428004659-phpapp01.docx
Sistemas Coordenados(2014-U)Word93
More From t7acarra

References: RESOLUCIÓN 
 resolución 
 resolución 
 resolución 
 resolución 
 RESOLUCIÓN 
 resolución 
 resolución 
 resolución