FEA.zip_FEA_finite element_matlab有限元_有限元分析_有限元模拟

clc clear nn=6; % nn=number of nodes ne=nn-1; %ne=number of elements nc=[0 0; %np=node coordinates 20 0; 40 0; 60 0; 80 0; 100 0]; ele=[1 2; %element connections 2 3; 3 4; 4 5; 5 6]; Lb=100; %length of beam L=Lb/ne; %length of element I=1; %moment of interia E=69e9; % young's modolus F=[0 0 -400 0 -400 0 -400 0 -400 0 -400 0]'; %force matrix %-------------------------------------------------- %plot the initial figure for e=1:ne x=[nc(ele(e,1),1) nc(ele(e,2),1)]; % For X Coordinates% y=[nc(ele(e,1),2) nc(ele(e,2),2)]; % For Y Coordinates% plot(x,y,'b') h=text(x,y, num2str(e)); hold on end %--------------stiffness matrix-------------- k=(E*I/(L*L*L))*[12 6*L -12 6*L; 6*L 4*L*L -6*L 2*L*L; -12 -6*L 12 -6*L; 6*L 2*L*L -6*L 4*L*L]; %-------calculation of global matrix------------- Kgl=zeros(2*nn); %initialization of global matrix for e=1:ne %for each element for i=ele(e,1) %i=first node of element for j=ele(e,2) %j=second node of element dof(1)=2*i-1; % x coordinate of node 1 dof(2)=2*i; dof(3)=2*j-1; dof(4)=2*j; for n1=1:4 %for row for n2=1:4 %for coloum Kgl(dof(n1),dof(n2))=Kgl(dof(n1),dof(n2))+k(n1,n2); %formula for global matrix end end end end end %---------------deflection----------------------- K=Kgl(3:12,3:12); f=F(3:12,1); D=inv(K)*f; % D=deflection disp('displacement for each degree of freedom'); for i=1:5 u(i+1,1)=D(2*i-1,1); end x=0:20:100; plot(x,u); % %------------------------------------------------------------------ % % %strain % for e=1:nn dx1=D(1,(2*e-1)); dy1=D(1,(2*e)); De=[dx1;dy1]; B=[-1 1 ]; Strain=De*B; disp('strain for each element');disp(e); disp(Strain) %--------------------------------------------------------------- %stress % v=0.33; Es=(E/1-v^2)*[1 v;v 1]; Stress=Es*Strain; disp('stress for each node');disp(e); disp(Stress); end