How can I plot a graph with a time dependent ode45?

4 次查看(过去 30 天)
clc
clear all
A = 0.05; %Excitation Amplitude
l_r = 2; %Wave length of the road
v = 45; %Speed(m/s)
P = l_r/v; %Period
Om = 1/P*2*pi; %Forcing Frequency
%Om = %0.07m,2m,45m/s
k_l = 26400; %Linear stiffness
m = 483; %Mass
d = -0.1; %Stretching condition
l = 0.5; %Length of the spring
k_s = -(k_l*(l-d))/(4*d); %Spring stiffness
f_n = sqrt(k_l/m)/(2*pi); %Natural frequency
%%
f = @(t,x) [ x(2); ...
-(4*k_s*(x(1)-A*sin(Om*t))*(sqrt((l-d)^2 + (x(1)-A*sin(Om*t))^2) - l))/ ...
(m*(sqrt((l-d)^2 + (x(1)-A*sin(Om*t))^2))) ];
T = 50;
x0 = [0,0];
[t,x] = ode45(f,[0,T],x0);
Response_amp = max(x(:,1)) - min(x(:,1));
plot(t,x(:,1))
xlabel('Time (s)')
ylabel('Amplitude (m)')
title('When d=-0.1', 'fontsize', 20)
set(gca,'FontSize',15)
Hi, all. This is my ode45 code. If I chagne A (Excitation Amplitude), Response_amp will have a different value. I wish to plot this relationship as A vs Response_amp.

采纳的回答

Ameer Hamza
Ameer Hamza 2020-3-16
编辑:Ameer Hamza 2020-3-16
You can visualize such relation with 3D plots
l_r = 2; %Wave length of the road
v = 45; %Speed(m/s)
P = l_r/v; %Period
Om = 1/P*2*pi; %Forcing Frequency
%Om = %0.07m,2m,45m/s
k_l = 26400; %Linear stiffness
m = 483; %Mass
d = -0.1; %Stretching condition
l = 0.5; %Length of the spring
k_s = -(k_l*(l-d))/(4*d); %Spring stiffness
f_n = sqrt(k_l/m)/(2*pi); %Natural frequency
%%
fig = figure();
ax = axes();
hold(ax);
view([-53 33]);
grid on
A_array = 0.05:0.05:0.3; %Excitation Amplitude
T = 15;
x0 = [0,0];
for i=1:numel(A_array)
A = A_array(i);
f = @(t,x) [ x(2); ...
-(4*k_s*(x(1)-A*sin(Om*t))*(sqrt((l-d)^2 + (x(1)-A*sin(Om*t))^2) - l))/ ...
(m*(sqrt((l-d)^2 + (x(1)-A*sin(Om*t))^2))) ];
[t, x] = ode45(f,[0,T],x0);
Response_amp = max(x(:,1)) - min(x(:,1));
plot3(t, A*ones(size(t)), x(:,1), 'LineWidth', 1);
end
xlabel('Time (s)')
ylabel('A (m)')
zlabel('Amplitude (m)')
title('When d=-0.1')

更多回答(0 个)

类别

Help CenterFile Exchange 中查找有关 Numerical Integration and Differential Equations 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by