Elemento Finito Con Matlab

INSTITUTO TECNOLÓGICO DE CELAYA DEPARTAMENTO DE POSGRADO DE MECANICA

REPORTE Calculo de la deformación de estructuras mediante el método de elemento finito en Matlab

Alumno Ing. Miguel Angel Corzo Velázquez

Titular M.I. Raúl Lesso Arroyo

Celaya, Guanajuato, 24 de febrero de 2012

Método de Elemento Finito con Matlab

DESARROLLO

PROBLEMA 01- ELEMENTOS LINEALES La pieza solida de sección irregular, se encuentra bajo las condiciones mostradas en la Fig. 1. De un lado tenemos el lado empotrado y en el otro extremo la fuerza “F” que desconocemos pero debemos calcular con la única limitante de no sobrepasar la deformación de 3 mm. Para la solución se usara elementos lineales y cuadráticos, y el material de acero y aluminio

Extremo empotrado

F D=19

L=200

D=12

L=150

D=19

L=200

D=10

L=150

Fig. 1 Elemento bajo condiciones de carga (medidas en mm)

El problema cuenta con los siguientes datos para su resolución: Datos F=F E = 70 Gpa E = 210 Gpa

La primera etapa para la resolución del problema es la designación de los nodos que tomaremos para nuestro análisis.


4

3

1

2

1

3

5

4

2

Fig. 2 Disposición elegida para los nodos y elementos (circulo celeste)

la siguiente tabla: Elemento 1 2 3 4

Nodo 1 3 4 5

Nodo 3 4 5 2

Tabla 1. Conectividad de elementos

SOLUCION MEDIANTE MATLAB Ya definida la tabla de conectividad de elementos podemos proceder a la solución del sistema mediante Matlab. Ya que se trata de un elemento con un grado de libertad las subrutinas a utilizar son las siguientes:

%Entrada de datos Analisis de un elemento E=210e9; Acero A1=pi*(.019)^2/4; A2=pi*(.012)^2/4; A3=pi*(.019)^2/4; A4=pi*(.01)^2/4; L1=.2; L2=.15; L3=L1; L4=L2;


U=.003; %Evaluacion de las matrices locales de los elementos k1= LinearBarElementStiffness(E,A1,L1) k2= LinearBarElementStiffness(E,A2,L2) k3= LinearBarElementStiffness(E,A3,L3) k4= LinearBarElementStiffness(E,A4,L4) %Ensamblar matriz global del rigidez KG=zeros(5,5); KG=LinearBarAssemble(KG,k1,1,3) KG=LinearBarAssemble(KG,k2,3,4) KG=LinearBarAssemble(KG,k3,4,5) KG=LinearBarAssemble(KG,k4,5,2)

%Sustraccion de matrices de rigidez mediante la aplicaciones de la condiciones de %frontera KG1=[KG(3:5,1:2)] KG2=[KG(3:5,3:5)]

%Armar vector de desplazamientos u12=[0;-U]

u1 = u2 = 2 =


%Obtencion de valores de u3,u4,u5 u345=inv(KG2)*(KG1*u12)

u3 = u4 = u5 = u3 =

UT=[0;U;u345]; FT=KG*UT

Fx1 = Fx2 = Fx3 = Fx4 = Fx5 =

El procedimiento para el aluminio es similar, por lo tanto la fuerza para producir dicha deformación para los dos materiales es la siguiente:

FAcero

=

135.57 kN

FAluminio

=

45.19 kN


PROBLEMA 02- ELEMENTO CUADRÁTICO La estructura y el tipo de material es el mismo que el problema anterior con la diferencia que se utilizara elementos cuadráticos para su solución.

1

1 6

3

2

3 7

4

8

4

5

9

2


Elemento 1 2 3 4

Nodo 1 3 4 5

Nodo 3 4 5 2

Nodo 6 7 8 9


SOLUCION MEDIANTE MATLAB Ya definida la tabla de conectividad de elementos, podemos proceder a la solución del sistema mediante Matlab. Entrada de datos a Matlab %Entrada de datos Analisis de un elemento E=210e9; Acero A1=pi*(.014)^2/4; A2=pi*(.012)^2/4; A3=pi*(.019)^2/4; A4=pi*(.01)^2/4;


L1=.2 L2=.15 L3=L1 L4=L2 U=.003 %Evaluacion de las matrices locales de los elementos

k1= QuadraticBarElementStiffness(E,A1,L1) k2= QuadraticBarElementStiffness(E,A2,L2) k3= QuadraticBarElementStiffness(E,A3,L3) k4= QuadraticBarElementStiffness(E,A4,L4)

%Ensamblar matriz global del rigidez KG=zeros(9,9); KG=QuadraticBarAssemble(KG,k1,1,3,6) KG=QuadraticBarAssemble(KG,k2,3,4,7) KG=QuadraticBarAssemble(KG,k3,4,5,8) KG=QuadraticBarAssemble(KG,k4,5,2,9)

%Sustraccion de matrices de rigidez mediante la aplicaciones de la condiciones de %frontera KG1=[KG(3:9,1:2)] KG2=[KG(3:9,3:9)]

%Armar vector de desplazamientos u12=[0;-U]

%Obtencion de valores de u u3456789=inv(KG2)*(KG1*u12)


u3 = u4 = u5 = u6= u7 = u8 = u9=

UT=[0;U;u3456789]; %Obtencion del vector de fuerzas totales

FT=KG*UT

Fx1 = Fx2 = Fx3 = Fx4 = Fx5 = Fx6 = Fx7 = Fx8 = Fx9 =

El procedimiento para el aluminio es similar, por lo tanto la fuerza para producir dicha deformación para los dos materiales es la siguiente:

FAcero

=

135.57 kN

FAluminio

=

45.19 kN


PROBLEMA 03Para la siguiente estructura. Los miembros tiene una sección transversal de 10 in2 y módulo de elasticidad de E= 29x106 lb/in2. Determinar la defección de cada junta, el esfuerzo en cada miembro y las reacciones en la base.

Fig. 4 Elemento bajo condiciones de carga (medidas en in)

La primera etapa para la resolución del problema es la designación de los nodos que tomaremos para nuestro análisis.


21

22 19

25

15

28 27

18

14

20

24

23

29

16

26

13

17

11

10

12 9

7

8

6 5

3 4

2

1


Elemento 1 2 3 4 5 6 7 8 9

Longitud 15 15 21.21 15 15 10.31 16 10.31 10

Nodo 1 1 1 2 3 3 4 4 5

Nodo 2 3 4 4 4 5 5 6 6

 0 90 45 90 0 75.96 141.34 204.03 0


10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29

10 14.14 10 10 14.14 10 10 11.18 11.18 5 11.18 10 10 10 11.18 10 11.18 5 10 11.18

5 5 6 7 8 13 7 7 9 9 9 11 12 8 10 14 8 10 15 10

7 8 8 8 13 14 13 9 10 12 11 12 13 14 14 15 10 15 16 16

90 45 90 0 135 0 90 153.43 26.56 90 153.43 0 0 90 153.43 0 26.56 90 0 26.56


SOLUCION MEDIANTE MATLAB Ya definida la tabla de conectividad de elementos podemos proceder a la solución del sistema mediante Matlab. Ya que se trata de un elemento con un grado de libertad las subrutinas a utilizar son las siguientes: %Dimensiones y material E=29e6;P=1e3; L1=15/12;L2=L1;L3=((15^2+15^2)^.5)/12;L4=L1;L5=L1;L6=((2.5^2+10^2)^.5)/12; L7=16;L8=L6;L9=10/12;L10=L9;L11=(200^.5)/12;L12=L9;L13=L9;L14=L11;L15=L9 L16=L9;L17=(125^.5)/12;L18=L17;L19=5/12;L20=L17;L21=L9;L22=L9;L23=L9; L24=L17;L25=L9;L26=L17;L27=L19;L28=L9;L29=L17; T1=0;T2=90;T3=45;T4=90;T5=0;T6=75.96;T7=141.34;T8=104.04;T9=0;T10=90; T11=45;T12=90;T13=0;T14=135;T15=0;T16=90;T17=153.44;T18=26.56;T19=90; T20=153.44;T21=0;T22=0;T23=90;T24=153.44;T25=0;T26=26.56;T27=90;T28=0; T29=26.56; A=10;


%Evaluacion de matrices por elemento k1 = PlaneTrussElementStiffness(E,A,L1, T1); k2 = PlaneTrussElementStiffness(E,A,L2, T2); k3 = PlaneTrussElementStiffness(E,A,L3, T3); k4 = PlaneTrussElementStiffness(E,A,L4, T4); k5 = PlaneTrussElementStiffness(E,A,L5, T5); k6 = PlaneTrussElementStiffness(E,A,L6, T6); k7 = PlaneTrussElementStiffness(E,A,L7, T7); k8 = PlaneTrussElementStiffness(E,A,L8, T8); k9 = PlaneTrussElementStiffness(E,A,L9, T9); k10 = PlaneTrussElementStiffness(E,A,L10, T10); k11 = PlaneTrussElementStiffness(E,A,L11, T11); k12 = PlaneTrussElementStiffness(E,A,L12, T12); k13 = PlaneTrussElementStiffness(E,A,L13, T13); k14 = PlaneTrussElementStiffness(E,A,L14, T14); k15 = PlaneTrussElementStiffness(E,A,L15, T15); k16 = PlaneTrussElementStiffness(E,A,L16, T16); k17 = PlaneTrussElementStiffness(E,A,L17, T17); k18 = PlaneTrussElementStiffness(E,A,L18, T18); k19= PlaneTrussElementStiffness(E,A,L19, T19); k20= PlaneTrussElementStiffness(E,A,L20, T20); k21= PlaneTrussElementStiffness(E,A,L21, T21); k22= PlaneTrussElementStiffness(E,A,L22, T22); k23= PlaneTrussElementStiffness(E,A,L23, T23); k24= PlaneTrussElementStiffness(E,A,L24, T24); k25= PlaneTrussElementStiffness(E,A,L25, T25); k26= PlaneTrussElementStiffness(E,A,L26, T26); k27= PlaneTrussElementStiffness(E,A,L27, T27); k28= PlaneTrussElementStiffness(E,A,L28, T28); k29= PlaneTrussElementStiffness(E,A,L29, T29); %Ensamble de la matriz goblal de rigidez K=zeros(32,32); KG = PlaneTrussAssemble(K,k1,1,2); KG = PlaneTrussAssemble(KG,k2,1,3); KG = PlaneTrussAssemble(KG,k3,1,4); KG = PlaneTrussAssemble(KG,k4,2,4); KG = PlaneTrussAssemble(KG,k5,3,4); KG = PlaneTrussAssemble(KG,k6,3,5); KG = PlaneTrussAssemble(KG,k7,4,6); KG = PlaneTrussAssemble(KG,k8,4,6); KG = PlaneTrussAssemble(KG,k9,5,6); KG = PlaneTrussAssemble(KG,k10,5,7); KG = PlaneTrussAssemble(KG,k11,5,8); KG = PlaneTrussAssemble(KG,k12,6,8);


KG = PlaneTrussAssemble(KG,k13,7,8); KG = PlaneTrussAssemble(KG,k14,8,13); KG = PlaneTrussAssemble(KG,k15,13,14); KG = PlaneTrussAssemble(KG,k16,7,13); KG = PlaneTrussAssemble(KG,k17,7,9); KG = PlaneTrussAssemble(KG,k18,9,10); KG = PlaneTrussAssemble(KG,k19,9,12); KG = PlaneTrussAssemble(KG,k20,9,11); KG = PlaneTrussAssemble(KG,k21,11,12); KG = PlaneTrussAssemble(KG,k22,12,13); KG = PlaneTrussAssemble(KG,k23,8,14); KG = PlaneTrussAssemble(KG,k24,10,14); KG = PlaneTrussAssemble(KG,k25,14,15); KG = PlaneTrussAssemble(KG,k26,8,10); KG = PlaneTrussAssemble(KG,k27,10,15); KG = PlaneTrussAssemble(KG,k28,15,16); KG = PlaneTrussAssemble(KG,k29,10,16) %Sustraer la matriz de desplazamientos

KG1=[KG(1:1,1:1) KG(1:1,5:32);KG(5:32,1:1) KG(5:32,5:32)] %Estableciendo vector de fuerzas Fu=[0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;0;-P;0;0;0;0;0;0;0;0;0;-P] %Obtención de desplazamientos U=inv(KG1)*Fu %Obtencion de reacciones UT=[U(1:1);0;0;0;U(2:29)] RT=KG*UT


u1 = v1 = u2 = v2 = u3 = v3 = u4 = v4 = u5 = v5 = u6 = v6 = u7 = v8 = u8 = v8 = u9 = v9 = u10 = v10= u11 = v11= u12 = v12= u13 = v13= u14 = v14= u15 = v15= u16 = v16=

%Encontrando los esfuerzos u1 =[UT(1:4)];

Fx1 = Fy1 = Fx2 = Fy2 = Fx3 = Fy3 = Fx4 = Fy4 = Fx5 = Fy5 = Fx6 = Fy6 = Fx7 = Fy7 = Fx8 = Fy8 = Fx9 = Fy9 = Fx10 = Fy10 = Fx11 = Fy11 = Fx12 = Fy12 = Fx13 = Fy13 = Fx14 = Fy14 = Fx15 = Fy15 = Fx16 = Fy16 =


S1 = PlaneTrussElementStress(E,L1,T1,u1) u2 =[UT(1:2);UT(5:6)]; S2 = PlaneTrussElementStress(E,L2,T2,u2) u3 =[UT(1:2);UT(7:8)]; S3 = PlaneTrussElementStress(E,L3,T3,u3) u4 =[UT(3:4);UT(7:8)]; S4 = PlaneTrussElementStress(E,L4,T4,u4) u5 =[UT(5:6);UT(7:8)]; S5 = PlaneTrussElementStress(E,L5,T5,u5) u6 =[UT(5:6);UT(9:10)]; S6 = PlaneTrussElementStress(E,L6,T6,u6) u7 =[UT(7:8);UT(9:10)]; S7 = PlaneTrussElementStress(E,L7,T7,u7) u8 =[UT(7:8);UT(11:12)]; S8 = PlaneTrussElementStress(E,L8,T8,u8) u9 =[UT(9:10);UT(11:12)]; S9 = PlaneTrussElementStress(E,L9,T9,u9) u10 =[UT(9:10);UT(13:14)]; S10 = PlaneTrussElementStress(E,L10,T10,u10) u11 =[UT(9:10);UT(15:16)]; S11 = PlaneTrussElementStress(E,L11,T11,u11) u12 =[UT(11:12);UT(15:16)]; S12 = PlaneTrussElementStress(E,L12,T12,u12) u13 =[UT(13:14);UT(15:16)];


S13 = PlaneTrussElementStress(E,L13,T13,u13) u14 =[UT(15:16);UT(25:26)]; S14 = PlaneTrussElementStress(E,L14,T14,u14) u15 =[UT(25:28)]; S15 = PlaneTrussElementStress(E,L15,T15,u15) u16 =[UT(13:14);UT(25:26)]; S16 = PlaneTrussElementStress(E,L16,T16,u16) u17 =[UT(13:14);UT(17:18)]; S17 = PlaneTrussElementStress(E,L17,T17,u17) u18 =[UT(17:20)]; S18 = PlaneTrussElementStress(E,L18,T18,u18) u19 =[UT(17:18);UT(23:24)]; S19 = PlaneTrussElementStress(E,L19,T19,u19) u20 =[UT(17:18);UT(21:22)]; S20 = PlaneTrussElementStress(E,L20,T20,u20) u21 =[UT(21:24)]; S21 = PlaneTrussElementStress(E,L21,T21,u21) u22 =[UT(23:26)]; S22 = PlaneTrussElementStress(E,L22,T22,u22) u23 =[UT(15:16);UT(27:28)]; S23 = PlaneTrussElementStress(E,L23,T23,u23) u24 =[UT(19:20);UT(27:28)]; S24 = PlaneTrussElementStress(E,L24,T24,u24)


u25 =[UT(27:30)]; S25 = PlaneTrussElementStress(E,L25,T25,u25) u26 =[UT(15:16);UT(19:20)]; S26 = PlaneTrussElementStress(E,L26,T26,u26) u27 =[UT(19:20);UT(29:30)]; S27 = PlaneTrussElementStress(E,L27,T27,u27) u28 =[UT(29:32)]; S28 = PlaneTrussElementStress(E,L28,T28,u28) u29 =[UT(19:20);UT(31:32)]; S29 = PlaneTrussElementStress(E,L29,T29,u29)

Elemento Finito Con Matlab

Recommend Documents