Problem with damping constant affecting displacement in a mass system displacement system

1 次查看(过去 30 天)
When no damping constant is present in my code it works as intended however when a damping constant is added the displacement increases which I believe shouldn't happen for my scenario. I cannot spot where the error is in my code.
% Define parameters
m = 10; % Mass (kg)
k = 160; % Spring constant (N/m)
b = 2; % Damping coefficient (N s/m)
F0 = 200; % Amplitude of forcing term (N)
v0 = [2, 3, 4, 5]; % Frequency of forcing term (rad/s)
x0 = -3; % Initial displacement (m)
v0_dot = 0; % Initial velocity (m/s)
% Define time span
tspan = linspace(0, 60, 1000); % from 0 to 20 seconds
% Define the function for ODE solver
for i = 1:length(v0)
omega = v0(i);
f = @(t, Y) [Y(2); (F0 * sin(omega * t) - k * Y(1) - b * Y(2)) / m];
% Solve using ODE45
[t_ode45{i}, Y_ode45{i}] = ode45(f, tspan, [x0, v0_dot]);
% Solve symbolically
syms x(t)
ode = diff(x, t, 2) == (F0 * sin(omega * t) - k * x - b * diff(x, t)) / m;
cond = [x(0) == x0, subs(diff(x), 0) == v0_dot]; % Corrected indexing for symbolic differentiation
x_sol = dsolve(ode, cond);
x_symbolic{i} = matlabFunction(x_sol);
% Plot results
subplot(length(v0), 1, i);
plot(t_ode45{i}, Y_ode45{i}(:, 1), 'b-', 'LineWidth', 1.5);
hold on;
t_symbolic = linspace(0, 60, 1000);
plot(t_symbolic, x_symbolic{i}(t_symbolic), 'r--', 'LineWidth', 1.5);
title(['Frequency = ', num2str(omega), ' rad/s']);
xlabel('Time (s)');
ylabel('Displacement (m)');
legend('ODE45', 'Symbolic');
grid on;
hold off;
end

采纳的回答

David Goodmanson
David Goodmanson 2024-2-9
Hi Daneil,
I don't see any signs that your calculation is incorrect, although in some cases you may have to run the calculation out further in time so that the solution can settle down. If you define
gamma = b/m
w0 = sqrt(k/m)
then for the underdamped case that you have, the requirement is gamma < 2*w0 and the expected amplitude for large times is
A = (F0/m)/sqrt((w^2-w0^2)^2 + gamma^2*w^2)
which agrees pretty well. With underdamping, the particular solutions fall off like exp((-gamma/2)*t) and a suitable wait time for the solution to down to a constant amplitude (so you get a correct result) is
waitime = 8/gamma
Here k/m = w0 = 4 and one of your frequencies is w=4, so in that case you would expect the solution to be quite large for small b, which is what happens. For fixed w, A decreases as b increases and I didn't see any indications otherwise.

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Symbolic Math Toolbox 的更多信息

标签

产品


版本

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by