3/24/25 9:36 AM C:\Users\hp\OneDr...\FRCTR_MECH_CODE_2.
m 1 of 2
clc; clear; close all;
% GIVEN PARAMETERS
L = 100; % Beam Length (m)
E = 200e9; % Young's Modulus (Pa)
I = 1e-6; % Moment of Inertia (m^4)
w = 1; % UDL (N/m)
num_elems = 100; % Number of finite elements
num_nodes = num_elems + 1;
% MESH GENERATION
le = L / num_elems;
x_nodes = linspace(0, L, num_nodes);
% INITIALIZATION
K = zeros(2*num_nodes);
F = zeros(2*num_nodes, 1);
% ELEMENT STIFFNESS MATRIX
ke = (E * I / le^3) * [12 6*le -12 6*le;
6*le 4*le^2 -6*le 2*le^2;
-12 -6*le 12 -6*le;
6*le 2*le^2 -6*le 4*le^2]; % Element stiffness matrix
% GLOBAL STIFFNESS MATRIX ASSEMBLY
for i = 1:num_elems
indices = [2*i-1, 2*i, 2*i+1, 2*i+2]; % DOFs for element
K(indices, indices) = K(indices, indices) + ke; % Assembly into global matrix
end
% APPLYING DISTRIBUTED LOAD AS NODAL FORCES
for i = 1:num_elems
fe = (w * le / 2) * [1; le/6; 1; -le/6]; % Equivalent nodal force vector
indices = [2*i-1, 2*i, 2*i+1, 2*i+2];
F(indices) = F(indices) + fe; % Add to global force vector
end
% APPLYING BOUNDARY CONDITIONS
fixed_dofs = [1, 2]; % Fixed end at x = 0 (zero displacement & zero rotation)
roller_dofs = 2*num_nodes-1; % Roller at x = L (zero displacement, free rotation)
free_dofs = setdiff(1:2*num_nodes, [fixed_dofs, roller_dofs]); % Remaining DOFs
% ================== SOLVING FOR DISPLACEMENTS ==================
U = zeros(2*num_nodes, 1); % Initialize displacement vector
U(free_dofs) = K(free_dofs, free_dofs) \ F(free_dofs); % Solve reduced system
% ================== EXTRACTING DEFLECTION RESULTS ==================
y_fem = U(1:2:end); % Extract vertical displacement
3/24/25 9:36 AM C:\Users\hp\OneDr...\FRCTR_MECH_CODE_2.m 2 of 2
% ================== ANALYTICAL SOLUTION ==================
x_exact = linspace(0, L, 1000);
y_exact = (w*x_exact)/(24*E*I) .* (L^3 - 2*L*x_exact.^2 + x_exact.^3);
% ================== PLOTTING RESULTS ==================
figure;
plot(x_nodes, y_fem, 'ro-', 'LineWidth', 2); hold on;
plot(x_exact, y_exact, 'b-', 'LineWidth', 2);
xlabel('Beam Length (m)');
ylabel('Deflection (m)');
title('FEM vs Analytical Beam Deflection' );
legend('FEM Solution', 'Analytical Solution');
grid on;