Info

此问题已关闭。 请重新打开它进行编辑或回答。

How to simplify this and make it efficient?

1 次查看(过去 30 天)
Donghun Lee
Donghun Lee 2020-5-6
关闭: Ameer Hamza 2020-5-6
lc, clear all
A= 0.06;
k_l = 26400; %Linear stiffness
m = 483; %Mass
l =0.5;
d =-0.005;
f_n = sqrt(k_l/m)/(2*pi); %Natural frequency
%%
Om_array = linspace(0,20,20); %in rad/s-1
l_array = linspace(0,1,20);
[om_array, L_array] = meshgrid(Om_array, l_array);
Response_amp = zeros(size(Om_array));
T = 150;
x0 = [0,0];
for i=1:numel(Om_array)
for j=1:numel(l_array)
Om = om_array(i,j);
l = L_array(i,j);
k_s = -(k_l*(l-d))/(4*d); %Spring stiffness
f = @(t,x) [ x(2); ...
-(2*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,[100,T],x0);
Response_amp(i,j) = (max(x(:,1)) - min(x(:,1)))/2;
% xval(i) = Om/(2*pi) ;
end
end
%% plot
figure(1);
ax = axes();
view(3);
hold(ax);
view([30 33]);
grid on
mesh(om_array/(2*pi),L_array,Response_amp) ;
xlabel('Frequency (Hz)')
ylabel('Length of the spring (m)')
zlabel('Response Amplitude (m)')
set(gca,'FontSize',15)
% set(gca,'xtick',[])
% set(gca,'ytick',[])
% set(gca,'ztick',[])
%%
%l = linspace(0,1,40);
%b = max(max(Response_amp));
hold on
d =-0.01;
Om_array = linspace(0,20,20); %in rad/s-1
l_array = linspace(0,1,20);
[om_array, L_array] = meshgrid(Om_array, l_array);
Response_amp = zeros(size(Om_array));
T = 150;
x0 = [0,0];
for i=1:numel(Om_array)
for j=1:numel(l_array)
Om = om_array(i,j);
l = L_array(i,j);
k_s = -(k_l*(l-d))/(4*d); %Spring stiffness
f = @(t,x) [ x(2); ...
-(2*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,[100,T],x0);
Response_amp(i,j) = (max(x(:,1)) - min(x(:,1)))/2;
% xval(i) = Om/(2*pi) ;
end
end
%% plot
figure(1);
ax = axes();
view(3);
hold(ax);
view([30 33]);
grid on
mesh(om_array/(2*pi),L_array,Response_amp) ;
xlabel('Frequency (Hz)')
ylabel('Length of the spring (m)')
zlabel('Response Amplitude (m)')
set(gca,'FontSize',15)
mesh(om_array/(2*pi),L_array,Response_amp) ;
hold off
Hi, all. This code shows 2 different graphs in the same figure when d = -0.005 and -0.01. However, I wish to simplify this code as it seems to be quite long-winded. Also, I want to vary d from -0.005 to -0.03.
Thanks for reading and I would appreciate your time for solving this problem.

回答(0 个)

此问题已关闭。

Community Treasure Hunt

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

Start Hunting!

Translated by