Local Stiffness Matrix Negative

Hey Everyone,

I am trying to construct a small 5 element, 2 noded finite element model of a tapered beam using matlab.

I currently have an acceptable answer, but I was wondering if there was any way that the clever people of the forum may have a more elegant solution.

Currently I am creating a global stiffness matrix using a connectivity matrix, describing the nodal connections. Since it is a beam, the nodes just line up with eachother simply.

In my nested loop I currently have the code at the bottom as a way of looping and checking if it is the corner of the sub-matrix 2x2. Do you all know if there is a better way of doing this so i dont have to use a if statement? It seems excessive.

B=[1 2;2 3;3 4;4 5;5 6]; %connectivity matrix
E=10e6 %youngs modulus, psi
L=24 %length of beam, inches
kg=zeros(6,6)
A=[11449/15625; 7921/15625; 5041/15625; 2809/15625; 49/625]; %area of cross section at each node
for e=1:5         %number of elements
   for i=1:length(B(1,:))      %number of nodes per element
       for j=1:length(B(1,:))  %number of nodes per element
           if i==2 & j==2
               k=A(e)*E/(L/5);
           elseif i==2 & j==1
               k=-A(e)*E/(L/5);
           elseif i==1 & j==2
               k=-A(e)*E/(L/5) ;    
           else
               k=A(e)*E/(L/5);
           end
           kg(B(e,i),B(e,j))=kg(B(e,i),B(e,j))+k;
       end
   end
end

Answer 1

[Alan Stevens]

How about replacing

           if i==2 & j==2
               k=A(e)*E/(L/5);
           elseif i==2 & j==1
               k=-A(e)*E/(L/5);
           elseif i==1 & j==2
               k=-A(e)*E/(L/5) ;    
           else
               k=A(e)*E/(L/5);
           end

by

            k=A(e)*E/(L/5)*((i==j) - (i~=j));