Matlab code: deflection of simply supported beam using ode45

9 次查看（过去 30 天）
Camryn
Camryn 2023-11-7
编辑： hamawandy 2025-3-7
Ran in:
Hi, my code keeps returning a function for a cantilever beam rather than a simply supported. How do I fix this??
My code:
clear;
clc;
x0 = 0;
Y0 = [0; 0; 0; 0];
L = 5000; % Length of the beam in mm
xRange = [x0, L];
[x, y] = ode45(@(x, y) beamDeflection(x, y, L), xRange, Y0);
plot(x, y);
xlabel('Length (mm)');
ylabel('Deflection (mm)');
title('Simply Supported Beam Deflection');
function dydx = beamDeflection(x, y, L)
E = 69000; % Young's modulus of aluminum in MPa
b = 500; % Length of cross-sectional base in mm
h = 750; % Height of cross-sectional height in mm
I = (b * h^3) / 12; % Moment of inertia about x
q = 10000; % Load in newtons per mm^2
if x == 0 || x == L
dydx = [0; 0; 0; 0];
else
dydx = [y(2); y(3); y(4); (q / (E * I))];
end
end

请先登录，再回答此问题。

回答（2 个）

Torsten
Torsten 2023-11-7
编辑：Torsten 2023-11-7
Ran in:
Hi, my code keeps returning a function for a cantilever beam rather than a simply supported. How do I fix this??
By using bvp4c instead of ode45. You have a boundary value problem, not an initial value problem because you try to impose boundary condition at x = 0 and x = L.
Here is a symbolic solution:
syms x y(x)
E = 69000; % Young's modulus of aluminum in MPa
b = 500; % Length of cross-sectional base in mm
h = 750; % Height of cross-sectional height in mm
I = (b * h^3) / 12; % Moment of inertia about x
q = 10000; % Load in newtons per mm^2
L = 5000;
Dy = diff(y,x);
D2y = diff(y,x,2);
D3y = diff(y,x,3);
D4y = diff(y,x,4);
eqn = D4y == q/(E*L);
conds = [y(0)==0,y(L)==0,D2y(0)==0,D2y(L)==0];
sol = dsolve(eqn,conds)
sol = 
fplot(sol,[0 L])
  7 个评论
显示 5更早的评论
Torsten
Torsten 2025-3-5
Ran in:
syms x y(x)
E = 69000; % Young's modulus of aluminum in MPa
b = 500; % Length of cross-sectional base in mm
h = 750; % Height of cross-sectional height in mm
I = (b * h^3) / 12; % Moment of inertia about x
q = 10000; % Load in newtons per mm^2
L = 5000;
Dy = diff(y,x);
D2y = diff(y,x,2);
D3y = diff(y,x,3);
D4y = diff(y,x,4);
eqn = D4y == q/(E*I);
conds = [y(0)==0,y(L)==0,D2y(0)==0,D2y(L)==0];
sol = dsolve(eqn,conds);
xn = 0:0.1:L;
yn = subs(sol,x,xn);
area(xn,yn)
●

请先登录，再进行评论。

hamawandy
hamawandy 2025-3-6
编辑：hamawandy 2025-3-6
Dear Sir
I am facing a problem in sketching the shear force diagram from the span 4m till 8m.......the line should be horizontal from 4m to 8m.
Can you help me to fix this error ?
This is the code:
R1=6;
R2=2;
X= linspace (0, 8, 100);
V= (R1-2*X).*(X>=0) - (R2)*(X>=4);
subplot (211)
plot (X, V)
area (X,V, 'FaceColor', 'r', 'LineWidth', 2);
  4 个评论
显示 2更早的评论
Torsten
Torsten 2025-3-7
Which function for M do you want to plot in a mathematical notation (not code) ?
hamawandy
hamawandy 2025-3-7
编辑：hamawandy 2025-3-7
Ran in:
Thank you so much I have fixed the error just now.
The final code for both shear force and bending moment diagrams is the following:
R1=6;
R2=2;
X= linspace (0, 8, 100);
V= (R1-2*X).*(X>=0).*(X<=4) - (R2)*(X>=4);
subplot (211)
plot (X, V)
area (X,V, 'FaceColor', 'r', 'LineWidth', 2);
M=(R1*X-X.^2 ).*(X>=0).*(X<4) + (R1*X-8*(X-2)).* (X>=4).* (X<=8);
subplot(212)
plot (X,M)
area(X,M, 'FaceColor', 'b', 'LineWidth', 2);
xline(0, 'Color', 'k', 'LineWidth', 2); % Draw line for Y axis
yline(0, 'Color', 'k', 'LineWidth', 2); % Draw line for X axis
●

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

产品


版本

R2023b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by