FEM ANALYSIS OF 1-D BAR ELEMENT WITH MATLAB
QUESTION:
Determine the nodal displacement of 1-D bar with uniform cross-section as shown in the above figure using elemental matrix equation.The length of the bar is 100mm and the area of cross-section is 150mm^2 .The value of Young's Modulus is 200KN/mm^2.The bar is fixed at one end and the force of 1000N is applied at the other free end.
AIM:
To find the displacement at the nodal points of a 1-D bar element fixed at one end and force applied at the free end by elemental matrix equation(FEM technique) in matlab.
THEORY:
What is Finite Element Method(FEM)?
- The finite element method (FEM) is a numerical technique used to perform finite element analysis (FEA) of any given physical phenomenon.
It is necessary to use mathematics to comprehensively understand and quantify any physical phenomena, such as structural or fluid behavior, thermal transport, wave propagation, and the growth of biological cells.
GOVERNING EQUATION:
K = `(A*E)/(Le)` `[[1,-1],[-1,1]]`
where,
K = elemental stiffness matrix
A = Area of cross section(`mm^2`)
E = Young's Modulus`(N/(mm)^2)`
Le = Length of each element(mm)
OBJECTIVE:
The main objective of the project is to find the displacement at the nodal points of the 1-D bar element fixed at one point by using elemental matrix equation.
WORKFLOW OF THE PROJECT:
- Initially we take the input of length,area,young's modulus and the number of nodes.
- Then we generate the zeros ,force and elemental stiffness matrix one by one.
- After that we obtain the global stiffness matrix.
- Then we obtain the displacement matrix by solving the global stiffness matrix and force matrix
- At last we make the plot of displacement at the nodal points.
PROGRAMMING LANGUAGE USED - MATLAB
MAIN CODE:
%A program for FEM analysis of 1-D bar with uniform cross section
%By using elemental matrix equation
%Bar is fixed at one end and free at other
clear all
close all
clc
disp(' FEM ANALYSIS OF 1-D BAR WITH UNIFORM CROSS-SECTION ')
L = input('Enter the length of the element(in mm) :');
A = input('Enter the area of the element(in mm^2) :');
E = input('Enter the value of Young Modulus(in N/mm^2) :');
nodes = input('Enter the number of nodes ');
elements = nodes - 1;
%length of each element
l_e = L/elements;
%Generating zeros matrix
Z = zeros(nodes);
%generating force matrix
F = input('Enter the value of force in column matrix(in N) : ') %Here we have to enter the same number of forces as the number of nodes
C = [1 -1 ;-1 1];
%Generating elemental matrix
K = ((A*E)/l_e)*C;
%Global stiffness matrix generation
for i =1:elements
range = [i:i+1];
Z(range,range) = Z(range,range) + K;
end
disp('Global Stiffness Matrix is')
disp(Z)
%Getting the number of fixed and free nodes
fixed_nodes = 1;
free_nodes = setxor(1:nodes,fixed_nodes); %setxor is used to get all the nodes except the fixed nodes
%Generating matrix of free nodes
Z_p = Z(free_nodes ,free_nodes);
F_p = F(free_nodes ,1);
%Getting the value of displacement
U_p = Z_p\F_p;
U = zeros(nodes,1);
disp('Resultat Displacement matrix at nodal points')
U(free_nodes) = U_p
%For getting deflection in X and Y direction
X = zeros(nodes,1);
Y = zeros(nodes,1);
for k = 2:nodes
X(k) = X(k-1) + l_e;
end
X_def = X+U; %Deflection in X direction
Y_def = Y; %Deflection in Y direction
%there is no force in Y direction so Y deflection is zero
%Plotting the displacement at the nodal points
patch(X_def,Y_def,U,'EdgeColor','interp','linewidth',10);
title('Resultant Displacement Plot');
colormap jet
cb = colorbar; %setting the colorbar
CODE EXPLAINATION:
- Initially we take the input of length ,area,Young's modulus and number of nodes by using input command.
- Then we generate the zeros ,force and elemental matrix.
- After that we obtain Global stiffness matrix by using the for loop.
- After that we obtain the displacement matrix by solving the force and global stiffness matrix.
- After that we plot the graph of displacement by using patch command and used colorbar to show the displacement range.
OUTPUT:
The below plot shows the displacement of the 1-D bar by applying the force in X-direction:
The below screenshot shows the values of output in the command window:
CONCLUSION:
Thus we have calculated the global stiffness matrix and the displacement matrix at the different nodal points.Also the plot of displacement was obtained as shown above.
ERRORS OBTAINED AND IT'S SOLUTION:
Initially i didn't used backslash(\) as F_p is 2`×`1 matrix and Z_p is 2`×`2 matrix so we cannot perform the operation as it is against matrix rule.This was rectified after using the backslash sign.
REFERENCES:
https://www.mathworks.com/help/matlab/visualize/introduction-to-patch-objects.html