Métodos Numéricos Aplicados EAP IH UNC 2012 -I §2. EL ELEMENTO DE BARRA LINEAL. 1. Ecuaciones b´ asicas. El elemen
Views 166 Downloads 2 File size 90KB
Métodos Numéricos Aplicados
EAP IH UNC
2012 -I
§2. EL ELEMENTO DE BARRA LINEAL.
1.
Ecuaciones b´ asicas.
El elemento de barras lineal es un elemento finito en una dimensi´on donde las coordenadas locales y lo globales coinciden. Es caracterizado por funciones de forma lineales y es id´entico al elemento de resorte elemental, salvo que la rigidez de la barra no se da directamente. El elemento lineal tiene m´ odulo de elasticidad E, ´area de la secci´ on transversal A y longitud L. Cada elemento de barra lineal tiene dos nodos, como se muestra en la Figura 1. En este caso la matriz de rigidez elemental est´a dada por
k=
i
"
EA L − EA L
− EA L
#
EA L
E, A
(1)
j x
L Figura 1: Elemento de barra lineal.
Es evidente que la matriz de rigidez elemental para el elemento de barra lineal es similar a la del resorte elemental con la rigidez reemplazada por EA/L. Es claro que el elemento de barra elemental tiene s´ olo dos grados de libertad - uno en cada nodo. Consecuentemente para una estructura con n nodos, la matriz de rigidez global K tendr´ a dimensi´on n × n (ya que tenemos un grado de libertad en cada nodo). La matriz de rigidez global K es ensamblada por la funci´on de MatLab LinearBarAssemble que est´a escrita espec´ıficamente para este prop´ osito. Este proceso se ilustra en detalle en los ejemplos. Una vez que ha sido obtenida la matriz de rigidez global K, tendremos la siguiente ecuaci´ on matricial:
[K] {U } = {F }
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
1
(2)
2012-1
donde U es el vector global de desplazamientos nodales y F es el vector de fuerzas global nodal. En este paso las condiciones de contorno se aplican manualmente a los vectores U y F . A continuaci´ on, el sistema (2) es resuelto por particionamiento seguido de eliminaci´on Gaussiana. Por u ´ltimo una vez que los desplazamientos y reacciones desconocidas han sido encontradas, las fuerzas de los elementos se obtienen para cada elemento de la siguiente manera:
{f } = [k] {u}
(3)
en donde f es un 2 × 1 vector elemental de fuerzas y u es un 2 × 1 vector elemental de desplazamientos. La tensi´on elemental es obtenida dividiendo las fuerzas elementales por el ´area de la secci´ on transversal A.
2.
Uso de las funciones de MatLab. Las cuatro funciones de MatLab usadas para el caso del elemento de barra lineal son: LinearBarElementStiffness(E,A,L) Esta funci´on calcula la matriz de rigidez elemental por cada barra elemental con m´ odulo de elasticidad E, ´area de secci´ on transversal A y longitud L. Retorna una matriz de rigidez elemental de 2 × 2. LinearBarAssemble(K, k, i, j) Estas funciones ensamblan la matriz de rigidez elemental k de la barra lineal, uniendo los nodos i (del extremo izquierdo) y j (del extremo derecho) en la matriz global de rigidez K. Se retorna una matriz global de rigidez K n × n cada vez que es ensamblado un elemento. LinearBarElementForces(k, u) Esta funci´on calcula el vector elemental de fuerzas usando la matriz de rigidez elemental k y el vector elemental de desplazamiento u. Se retorna un 2 × 1 vector elemental de fuerzas f . LinearBarElementStresses(k, u, A) Esta funci´on calcula el vector elemental de tensi´on usando la matriz elemental de rigidez k, el vector elemental de desplazamiento u y el ´area de secci´ on transversal A. Retorna un 2 × 1 vector elemental de tensi´on sigma o s.
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
2
2012-1
La siguiente es una lista de c´odigos fuente en MatLab para cada funci´on:
function y = LinearBarElementStiffness(E,A,L) % LinearBarElementStiffness Esta funci´ on retorna la matriz de rigidez % elemental para una barra lineal con m´ odulo de elasticidad E, ´ area de % secci´ on transversas A y longitud L. % La dimensi´ on de la matriz de rigidez elemental es 2 x 2. y = [E*A/L -E*A/L ; -E*A/L E*A/L];
function y = LinearBarAssemble(K,k,i,j) % % % %
LinearBarAssemble Esta funci´ on ensambla la matriz de rigidez elemental k de una barra lineal con nodos i y j en la matriz de rigidez global K. Esta funci´ on retorna una matriz de rigidez global K luego que hayan sido ensambladas las matrices de rigidez elemental.
K(i,i) K(i,j) K(j,i) K(j,j) y = K;
= = = =
K(i,i) K(i,j) K(j,i) K(j,j)
+ + + +
k(1,1) k(1,2) k(2,1) k(2,2)
; ; ; ;
function y = LinearBarElementForces(k,u) % LinearBarElementForces Esta funci´ on retorna un vector de fuerza % nodal elemental dada la matriz de rigidez elemental k y % el vector elemental de desplazamiento nodal u. y = k * u;
function y = LinearBarElementStresses(k, u, A) % LinearBarElementStresses Esta funci´ on retorna el vector de tensi´ on nodal % elemental dada la matriz de rigidez elemental k, el vector elemental de % desplazamiento nodal u y el a ´rea de secci´ on transversal A. y = k * u/A;
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
3
2012-1
Ejemplo 1 Considerar la estructura conformada por dos barras lineales como se muestra en la Fig. 2. Sea E = 210 GP a, A = 0,003 m2 , P = 10 kN y el nodo 3 es desplazado hacia la derecha en 0,002 m. Determinar: 1. la matriz de rigidez global para la estructura. 2. el desplazamiento en el nodo 2. 3. las reacciones en los nodos 1 y 3. 4. la tensi´ on en cada barra. P 3
2
1
1m
1.5 m
Figura 2: Estructura con dos barras para el ejemplo 1.
´ SOLUCION Paso 1. Discretizaci´ on del dominio: Este problema ya est´ a discretizado. El dominio ha sido dividido en dos elementos y tres nodos. Las unidades usadas en los c´ alculos del MatLab son kN y metros. La Tabla 1 muestra la conectividad elemental para este ejemplo: TABLA 1. Conectividad elemental para el ejemplo 1. N´ umero de elemento 1 2
Nodo i 1 2
Nodo j 2 3
Paso 2. Escritura de las matrices elementales de rigidez: Las dos matrices de rigidez elemental k1 y k2 son obtenidas por medio de la funci´ on de MatLab LinearBarElementStiffness. Cada matriz tiene dimensi´ on 2 × 2.
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
4
2012-1
>> E=210e6 E = 210000000 >> A=0.003 A = 0.0030 >> L1=1.5 L1 = 1.5000 >> L2=1 L2 = 1 >> k1=LinearBarElementStiffness(E,A,L1) k1 = 420000 -420000
-420000 420000
>> k2=LinearBarElementStiffness(E,A,L2) k2 = 630000 -630000
-630000 630000
Paso 3. Ensamble de la matriz de rigidez global. Dado que la estructura tiene tres nodos, el tama˜ no de la matriz de rigidez global es de 3×3. Por lo tanto para obtener K primeramente debemos construir una matriz nula de 3 × 3 luego hacer dos llamados a la funci´ on de MatLab LinearBarAssemble dado que tenemos dos barras lineales elementales en la estructura. Cada llamada a la funci´ on ensamblar´ a un elemento, MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
5
2012-1
>> K=zeros(3,3) K = 0 0 0
0 0 0
0 0 0
>> K=LinearBarAssemble(K,k1,1,2) K = 420000 -420000 0
-420000 420000 0
0 0 0
>> K=LinearBarAssemble(K,k2,2,3) K = 420000 -420000 0
-420000 105000 -630000
0 -630000 630000
Paso 4. Aplicaci´ on de las condiciones de contorno. El sistema (2) para este caso es obtenido utilizando la matriz de rigidez global del paso anterior: U1 F1 −420000 1050000 −630000 U F = 2 2 0 −630000 630000 U3 F3
420000
−420000
0
(4)
Las condiciones de contorno para este problema son dadas por: U1 = 0, Insertando ambas condiciones 420000 −420000 0
F2 = −10,
U3 = 0,002
(5)
en (4), obtenemos: 0 1050000 −630000 U2 −630000 630000 0,002 −420000
0
F 1 −10 = F3
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
6
(6)
2012-1
Paso 5. Resolviendo las ecuaciones. La soluci´ on del sistema de ecuaciones (6) se llevar´ a a cabo mediante un particionamiento (manual) y posterior eliminaci´ on Gaussiana (con MatLab). Primeramente particionamos (6) extrayendo la submatriz de la fila 2 y columna 2 que resulta ser una matriz 1 × 1. Debido al desplazamiento de 0.002 m efectuado por el nodo 3, tenemos que extraer la submatriz de la fila 2 y columna 3, que tambi´en resulta ser una matriz 1 × 1 Luego obtenemos:
[1050000] U2 + [−630000] (0,002) = {−10}
(7)
La soluci´ on de ambos sistemas es obtenido usando MatLab del siguiente modo. Notar que el operador \ (backslash) es usado para la eliminaci´ on Gaussiana.
>> k=K(2,2) k = 1050000 >> k0=K(2,3) k0 = -630000 >> u0=0.002 u0 = 0.0020 >> f=[-10] f = -10 >> f0=f-k0*u0 f0 = 1250 >> u=k\f0 MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
7
2012-1
u = 0.0012
Ahora est´ a claro que el desplazamiento en el nodo 2 es 0.0012 m. Paso 6. Post-procesamiento. En este paso obtenemos las reacciones en los nodos 1 y 3, y la tensi´ on en cada barra usando MatLab del siguiente modo. Primeramente creamos el vector global de desplazamientos nodales U y luego calculamos el vector global de fuerzas nodales F .
>> U=[0;u;u0] U = 0 0.0012 0.0020 >> F=K*U F = -500.0000 -10.0000 510.0000
As´ı, las reacciones en los nodos 1 y 3 son fuerzas de 500 kN (dirigida hacia la izquierda) y 510 kN (dirigida hacia la derecha), respectivamente. Queda claro que la fuerza de equilibrio se cumple. Luego establecemos los vectores elementales de desplazamiento nodal u1 y u2 ; luego calculamos el vector elemental de fuerzas f1 y f2 por medio de la funci´ on de MATLAB LinearBarElementForces. Finalmente dividimos cada fuerza elemental entre el a ´rea de la secci´ on transversal del elemento para obtener las tensiones elementales,
>> u1=[0;U(2)] u1 = 0 0.0012
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
8
2012-1
>> f1=LinearBarElementForces(k1,u1) f1 = -500.0000 500.0000 >> sigma1=f1/A sigma1 = 1.0e+005 * -1.6667 1.6667 >> u2=[U(2) ; U(3)] u2 = 0.0012 0.0020 >> f2=LinearBarElementForces(k2,u2) f2 = -510.0000 510.0000 >> sigma2=f2/A sigma2 = 1.0e+005 * -1.7000 1.7000
La tensi´ on en el elemento 1 es 1,667 × 105 kN/m2 (o 166,7 M P a de tensi´ on); y la tensi´ on en el 5 2 elemento 2 es 1,7 × 10 kN/m (o 170 M P a de tensi´ on). Alternativamente, podemos obtener la tensi´ on elemental directamente por medio de la funci´ on de MatLab LinearBarElementStresses. Esto se realiza de la siguiente manera, obtieni´endose los mismos resultados.
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
9
2012-1
>> s1=LinearBarElementStresses(k1,u1,A) s1 = 1.0e+005 * -1.6667 1.6667 >> s2=LinearBarElementStresses(k2,u2,A) s2 = 1.0e+005 * -1.7000 1.7000
´Indice 1. Ecuaciones b´ asicas.
1
2. Uso de las funciones de MatLab.
2
MÉTODOS NUMÉRICOS APLICADOS
EAP IH - UNC
10
2012-1