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

4 visualizaciones (últimos 30 días)
Donghun Lee
Donghun Lee el 16 de Mzo. de 2020
Comentada: Ameer Hamza el 16 de Mzo. de 2020
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.

Respuesta aceptada

Ameer Hamza
Ameer Hamza el 16 de Mzo. de 2020
Editada: Ameer Hamza el 16 de Mzo. de 2020
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')
  2 comentarios
Donghun Lee
Donghun Lee el 16 de Mzo. de 2020
Omg.. Thank you sooo much. You've just saved my life! Thanks!!
Ameer Hamza
Ameer Hamza el 16 de Mzo. de 2020
Glad to be of help.

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Numerical Integration and Differential Equations en Help Center y File Exchange.

Community Treasure Hunt

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

Start Hunting!

Translated by