110 likes | 227 Views
Project I Tutorial. MECH 300H Spring 2008. Project Description. Write a FEM code using MATLAB for stress analysis of a bridge as described below.
E N D
Project I Tutorial MECH 300H Spring 2008
Project Description Write a FEM code using MATLAB for stress analysis of a bridge as described below. A small railroad bridge is constructed of steel members, all of have a cross-sectional area of 3250 mm2, Young’s modulus of 1GPa. A train stops on the bridge, and the loads applied to the truss on one side of the bridge are as show in the following figure. Determine the nodal displacements and element stresses.
Project Description What to submit? (1) Both the hard copy and the electronic copy of your code with sufficient comments. (Send the electronic copy to honggang@ust.hk) In your code, the following major components must be included (a) formulation of element stiffness matrix, (b) stiffness matrix assembling, (c) construction of force vector, (d) condensed matrix, (e) linear system solver, and (f) post computation. (2) A short report – description of your implementation, such as number of nodes, elements, results and discussion. Notes: Independent work is required! Grading will be based on (a) the correctness of your results, (b) report writing (complete yet concise) and (c) the generality of your code.
Project task • Scheme is the most important! • You should know how to solve this problem by hand. • Project task: writing the code • Tell the computer how to solve this problem automatically.
Program steps The program consists of these steps. • InputData Input the information of nodes and elements • FormEquations Form the local and global Equation • BoundaryCondition Add the B.C. • CondenseEquations Condense the system equations • SolveEquations Solve the equations • Postcalculation Calculation the displacements and stress
Program steps [Main program] % 1D Bar ELEMENT ON 2D PLANE TRUSS ANALYSIS USING FEM function [] = bar(); clear all; close all; %The program consists of these steps. InputData; FormEquations; BoundaryCondition; CondenseEquations; SolveEquations; Postcalution;
Input data • The number of nodes and elements • NN=7,NE=11%dimension of global matrix: 2NN* 2NN • The information of node • x(NN), y(NN) % -x , y coordinate • The information of element • E(NE), A(NE), L(NE), Angle(NE) %properties of each bar • Noc=[1 2; 1 3; 3 2; 2 4; 3 4; 3 5; 5 4; 4 6;5 6; 5 7; 7 6] • %connect joints of each bar • The information of B. C • BC_U=[0 0 1 1 1 1 1 1 1 1 1 1 1 0 ]'; • %displacement B.C. 0-known; 1-unknown • BC_F=1*10e3[0 -280 0 0 0 -210 0 0 0 -280 0 0 0 -360]'; • %force B.C.
FormEquations Task: The local stiffness matrix is assembled into the global matrix. For the Nth Element N1=Noc(N,1)%Index of the first node N2=Noc(N,2) %Index of the second node N1*2-1 N1*2 N2*2-1 N2*2 Global index 1 2 3 4 Local index 1 N1*2-1 2 N1*2 3 N2*2-1 4 N2*2
FormEquations Pseudo code: Do N = 1, NE N1=Noc(N,1) N2=Noc(N,2) KG(N1*2-1,N1*2-1)=KG(N1*2-1,N1*2-1)+KL(1,1) KG(N1*2-1,N1*2 )=KG(N1*2-1,N1*2 )+KL(1,2) KG(N1*2-1,N2*2-1)=KG(N1*2-1,N2*2-1)+KL(1,3) KG(N1*2-1,N2*2 )=KG(N1*2-1,N2*2 )+KL(1,4) KG(N1*2 ,N1*2-1)=KG(N1*2 ,N1*2-1)+KL(2,1) KG(N1*2 ,N1*2 )=KG(N1*2 ,N1*2 )+KL(2,2) KG(N1*2 ,N2*2-1)=KG(N1*2 ,N2*2-1)+KL(2,3) KG(N1*2 ,N2*2 )=KG(N1*2 ,N2*2 )+KL(2,4) KG(N2*2-1,N1*2-1)=KG(N2*2-1,N1*2-1)+KL(3,1) KG(N2*2-1,N1*2 )=KG(N2*2-1,N1*2 )+KL(3,2) KG(N2*2-1,N2*2-1)=KG(N2*2-1,N2*2-1)+KL(3,3) KG(N2*2-1,N2*2 )=KG(N2*2-1,N2*2 )+KL(3,4) KG(N2*2 ,N1*2-1)=KG(N2*2 ,N1*2-1)+KL(4,1) KG(N2*2 ,N1*2 )=KG(N2*2 ,N1*2 )+KL(4,2) KG(N2*2 ,N2*2-1)=KG(N2*2 ,N2*2-1)+KL(4,3) KG(N2*2 ,N2*2 )=KG(N2*2 ,N2*2 )+KL(4,4) Hint: Due to the symmetry of the Matrix, you can write the code more short and efficient.
CondenseEquations Task: The global systems are condensed in order to solve it. BC_U 0 01 … 1 0 0 0 1 … 1 0 Scheme: If BC_U(N) ==1, then the element in the Nth line and Nth column should be reserved.
CondenseEquations Pseudo code: Subscript: I ,J - the index of KG(I ,J) CI,CJ - the index of KC(CI,CJ) CI=0 Do I=1 , NN*2 If (BC_U(I) == 1) Then CI=CI+1 CJ=0 Do J=1 , NN*2 If (BC_U(J) == 1) Then CJ=CJ+1 KC(CI,CJ)=KG(I,J) FC(CI)=BC_F(I)