Question

Problem 2: Calculate the force in all members of the truss shown in Fig. 2 by using MATLAB code. 50 kN 30 KN 4 m 5 m E 3 m 4 m 4 m Figure 2. Schematic of problem 2.

Answer 1

%Specify structure properties

nnodes=5; %specify the number of nodes

nmem=7; %specify the number of members

%Define nodal coordinates in order starting with node one

x(1)=0; y(1)=0;

x(2)=4; y(2)=3;

x(3)=12; y(3)=3;

x(4)=16; y(4)=0;

x(5)=8; y(5)=0;

%Define member connectivity start with member one

% Start node end node

mconn(1,1)=1; mconn(1,2)=2;

mconn(2,1)=2; mconn(2,2)=3;

mconn(3,1)=3; mconn(3,2)=4;

mconn(4,1)=5; mconn(4,2)=4;

mconn(5,1)=1; mconn(5,2)=5;

mconn(6,1)=5; mconn(6,2)=2;

mconn(7,1)=5; mconn(7,2)=3;

%input properties for each member

%For now E is constant.

E=200*10^6;

%A for each member

A(1)=0.0015;

A(2)=0.0015;

A(3)=0.0015;

A(4)=0.0015;

A(5)=0.0015;

A(6)=0.0015;

A(7)=0.0015;

%Specify supports

nsupp=3; %number of degrees of freedom with supports

sup(1)=1;%for the 1 to nsupp indicate which dof are supports

sup(2)=2;

sup(3)=8;

%sup(4)=8;

%sup(5)=1;

%sup(6)=2;

%Specify forces

%zero all the forces

for i=1:nnodes*2;

F(i)=0.0; %#ok<*SAGROW>

end;

nforce=7; %number of degrees of freedom with forces given

F(1)=0; %specify force at each degree of freedom where the external force is known, F(dof)

F(2)=0;

F(3)=0;

F(4)=1;

F(5)=0;

F(6)=0;

F(7)=0;

%***INPUT ENDS HERE ***

%Calculate the truss member lengths, and angles

for i=1:nmem;

dx=x(mconn(i,2))-x(mconn(i,1));

dy=y(mconn(i,2))-y(mconn(i,1));

L(i)=sqrt(dx^2+dy^2);

c(i)=dx/L(i);

s(i)=dy/L(i);

end;

%For each element create the global stiffness matrix and put it into the

%global stiffness matrix.

%zero the global stiffness matrix

for i=1:nnodes*2;

for j=1:nnodes*2;

kg(i,j)=0.0;

end;

end;

%Create each element global stiffness matrix as transpose[T][k][T]

for m=1:nmem; %loop over each element

%zero k and T

for i=1:7;

for j=1:7;

k(i,j)=0.0;

T(i,j)=0.0;

end;

end;

%create T

T(1,1)=c(m) ; T(1,2)=-1; T(1,3)=0; T(1,4)=0; T(1,5)=0; T(1,6)=-c(m); T(1,7)=0;

T(2,1)=s(m); T(2,2)=0; T(2,3)=0;T(2,4)=0; T(2,5)=0; T(2,6)=s(m); T(2,7)=0;

T(3,1)=0; T(3,2)=1; T(3,3)=-c(m); T(3,4)=0; T(3,5)=0; T(3,6)=0; T(3,7)=c(m);

T(4,1)=0; T(4,2)=0; T(4,3)=s(m); T(4,4)=0; T(4,5)=0; T(4,6)=0; T(4,7)=s(m);

T(5,1)=0; T(5,2)=0; T(5,3)=c(m); T(5,4)=1; T(5,5)=0; T(5,6)=0; T(5,7)=0;

T(6,1)=0; T(6,2)=0; T(6,3)=0;T(6,4)=-1; T(6,5)=1; T(6,6)=c(m); T(6,7)=-c(m);

T(7,1)=0; T(7,2)=0; T(7,3)=0; T(7,4)=0; T(7,5)=0; T(7,6)=-s(m); T(7,7)=-s(m);

%create k

k(1,1)=37500; k(1,2)=0; k(1,3)=0; k(1,4)=0; k(1,5)=0; k(1,6)=0; k(1,7)=0;

k(2,1)=0; k(2,2)=60000; k(2,3)=0; k(2,4)=0; k(2,5)=0; k(2,6)=0; k(2,7)=0;

k(3,1)=0; k(3,2)=0; k(3,3)=37500; k(3,4)=0; k(3,5)=0; k(3,6)=0; k(3,7)=0;

k(4,1)=0; k(4,2)=0; k(4,3)=0; k(4,4)=60000; k(4,5)=0; k(4,6)=0; k(4,7)=0;

k(5,1)=0; k(5,2)=0; k(5,3)=0; k(5,4)=0; k(5,5)=60000; k(5,6)=0; k(5,7)=0;

k(6,1)=0; k(6,2)=0; k(6,3)=0; k(6,4)=0; k(6,5)=0; k(6,6)=37500; k(6,7)=0;

k(7,1)=0; k(7,2)=0; k(7,3)=0; k(7,4)=0; k(7,5)=0; k(7,6)=0; k(7,7)=37500;

%transform k into the element global stiffness matrix

k=T'*k*T;

%put each element k into the global stiffness matrix kg

for i=1:2;

for j=1:2;

kg(mconn(m,i)*2-1,mconn(m,j)*2-1)=kg(mconn(m,i)*2-1,mconn(m,j)*2-1)+k(i*2-1,j*2-1);

kg(mconn(m,i)*2,mconn(m,j)*2)=kg(mconn(m,i)*2,mconn(m,j)*2)+k(i*2,j*2);

kg(mconn(m,i)*2-1,mconn(m,j)*2)=kg(mconn(m,i)*2-1,mconn(m,j)*2)+k(i*2-1,j*2);

kg(mconn(m,i)*2,mconn(m,j)*2-1)=kg(mconn(m,i)*2,mconn(m,j)*2-1)+k(i*2,j*2-1);

end;

end; %now kg should be complete

end; %end loop over all the elements

%Put kg into kgs and modify kgs to solve for the unknown displacements

%we save kg so that we can calculate the reactions also.

kgs=kg;

%modify kgs based on the support conditions. Zero the entire row for each

%dof sup(i) except for entry kgs(sup(i),sup(i)) set it equal to 1. When we solve

%for the displacements we have essentially set them equal to zero since

%they are supports.

for i=1:nsupp; %loop through all the specified supports

for j=1:nnodes*2;

kgs(sup(i),j)=0.0; %zero the whole row

end;

kgs(sup(i),sup(i))=1.0; %put a 1 in position kgs(sup(i),sup(i))

end;

kgs %#ok<*NOPTS>

%solve for the global nodal displacements

display('The global nodal displacements are:')

G=inv(kgs)*F' %#ok<*MINV>

%Solve for the global external forces, here we use kg that we set aside

%above.

display('The global external forces are:')

F=kg*G

%Calculate the forces for each member

for m=1:nmem; %loop over each element

%zero k and T

for i=1:7;

for j=1:7;

k(i,j)=0.0;

T(i,j)=0.0;

end;

end;

%create T

T(1,1)=c(m) ; T(1,2)=-1; T(1,3)=0; T(1,4)=0; T(1,5)=0; T(1,6)=-c(m); T(1,7)=0;

T(2,1)=s(m); T(2,2)=0; T(2,3)=0;T(2,4)=0; T(2,5)=0; T(2,6)=s(m); T(2,7)=0;

T(3,1)=0; T(3,2)=1; T(3,3)=-c(m); T(3,4)=0; T(3,5)=0; T(3,6)=0; T(3,7)=c(m);

T(4,1)=0; T(4,2)=0; T(4,3)=s(m); T(4,4)=0; T(4,5)=0; T(4,6)=0; T(4,7)=s(m);

T(5,1)=0; T(5,2)=0; T(5,3)=c(m); T(5,4)=1; T(5,5)=0; T(5,6)=0; T(5,7)=0;

T(6,1)=0; T(6,2)=0; T(6,3)=0;T(6,4)=-1; T(6,5)=1; T(6,6)=c(m); T(6,7)=-c(m);

T(7,1)=0; T(7,2)=0; T(7,3)=0; T(7,4)=0; T(7,5)=0; T(7,6)=-s(m); T(7,7)=-s(m);

%create k

k(1,1)=37500; k(1,2)=0; k(1,3)=0; k(1,4)=0; k(1,5)=0; k(1,6)=0; k(1,7)=0;

k(2,1)=0; k(2,2)=60000; k(2,3)=0; k(2,4)=0; k(2,5)=0; k(2,6)=0; k(2,7)=0;

k(3,1)=0; k(3,2)=0; k(3,3)=37500; k(3,4)=0; k(3,5)=0; k(3,6)=0; k(3,7)=0;

k(4,1)=0; k(4,2)=0; k(4,3)=0; k(4,4)=60000; k(4,5)=0; k(4,6)=0; k(4,7)=0;

k(5,1)=0; k(5,2)=0; k(5,3)=0; k(5,4)=0; k(5,5)=60000; k(5,6)=0; k(5,7)=0;

k(6,1)=0; k(6,2)=0; k(6,3)=0; k(6,4)=0; k(6,5)=0; k(6,6)=37500; k(6,7)=0;

k(7,1)=0; k(7,2)=0; k(7,3)=0; k(7,4)=0; k(7,5)=0; k(7,6)=0; k(7,7)=37500;

%Extract the global displacements associated with element m

u(1)=G(mconn(m,1)*2-1);

u(2)=G(mconn(m,1)*2);

u(3)=G(mconn(m,2)*2-1);

u(4)=G(mconn(m,2)*2);

u(5)=G(mconn(m,2)*2-1);

u(6)=G(mconn(m,2)*2);

u(7)=G(mconn(m,2)*2-1);

%calculate the local element end forces

f=k*T*u';

%save the important end forces

fm(m)=f(3); %positive values mean truss member in tension, negative means compression

end;%end loop over each element

display('The member forces are:')

fm'

%plot the truss on a graph

plot(0,0,'w')

hold

%loop over all the elements and plot each element, one at a time

magnify=2000; %magnification factor for deflections

for m=1:nmem;

xm(1)=x(mconn(m,1));

ym(1)=y(mconn(m,1));

xm(2)=x(mconn(m,2));

ym(2)=y(mconn(m,2));

xu(1)=G(mconn(m,1)*2-1)*magnify;

yu(1)=G(mconn(m,1)*2)*magnify;

xu(2)=G(mconn(m,2)*2-1)*magnify;

yu(2)=G(mconn(m,2)*2)*magnify;

if fm(m) < 0.0;

plot(xm,ym,'k'); %undisplaced truss members in compression will be blue

else

plot(xm,ym,'b'); %undisplaced truss members in tension will be red

end;

plot(xm+xu,ym+yu,'r--'); %displaced truss members

end;

xlabel('xcoordinates')

ylabel('ycoordinates')

title('Undeflected(solid) and Magnified Deflected(dashed) Truss')

hold