MATLAB Answers

Assembly of Global Stiffness Matrix

120 views (last 30 days)
hamzah khan
hamzah khan on 9 Nov 2019
Edited: hamzah khan on 9 Nov 2019
Hey Guys,
I am developing the code for a special 1D FEM problem for beam elements where each element has 2 nodes but each node has 6 degree of freedom.
So I have per element stiffness matrix of 12x12.
If I want to assemble its global stiffness matrix, I am having trouble how to do it.
I have everything else ready, element connectivity matrix, nodal list and also indivual element stiffness matrix, but I dont know how to do it using 2 nodes.
I know it can be done if my element has 12 nodes each having 1 degree of freedom, but with 2 nodes and 12 degree of freedom. I dont know how to do it.
I am attaching my code, If it is possible, kindly tell me how to do it?
% number of elements in the beam
nE =2
% generate number of nodes
nN=nE+1;
% DOFs per node
dpn = 6;
% Nodes per element
npe = 2;
% step size
step=(xb-xa_1)/nE;
if nE<1
error('Number of elements <1');
end
%initiliasation ( speeds up calculations)
Nod= zeros(nN,1);
Elm = zeros(nE,2);
%generate list of nodes
for i=1:1:nN
Nod(i,1) = xa + (i-1)*(xb-xa)/nE;
end
%generate connectivity array for elelments:
for kk=1:1:nE
Elm(kk,1)=kk;
Elm(kk,2)=kk+1;
end
%nodelabel and mark if node is a degree of dfreedon ( nodeLabel=0)
%or a boubdary node ( nodeLabel>0) and assign index of boundaryValue
nodeLabel=zeros(size(Nod));
nodeLabel(1) = 1 ;
nodeLabel(nN)=2;
%initialisation of A matrix
A=sparse(nN,nN);
disp('matrix assemly')
tic % start stopwatch
for n=1:1:nE
Se=alpha*local_element_matrix_K % Element Stiffness Matrix
Te=beta*local_element_matrix_M % Element Mass Matrix
Ge= local_element_matrix_G % Element Gyroscopic Matrix
Ce= (Se + Te)+Ge; % Element Damping Matrix
% Total Element Stiffness Matrix for 1 element
Ke=Se+Te+Ce; % size(Ke) = 12x12(Linear interpolation in 1D, multiple DOF per node)
for i=1:1:size(Ke,1)
glo_i= Elm(n,i);
if(nodeLabel(glo_i)>0)%Dirichlet boundary?%
A(glo_i,glo_i)=1.0;
else
for j = 1:1:size(Ke,2)
glo_j=Elm(n,j);
A(glo_i,glo_j)= A(glo_i,glo_j)+ Ke(i,j);
end
end
end %i
end %k
disp('')
spy(A) %

  0 Comments

Sign in to comment.

Answers (0)


Translated by