clear all
n=0.0005;
ro=2600;
E0=2.3*10^5;
E1=5.3*10^4;
E2=-1.6*10^4;
A=20*10^(-6);
pe=0.005;
l=0.1;
P=0.5;
T=10;
mu=n/ro/A;
alpha=(1-pe*T)*E1*pi^2/ro/l^3/3;
beta=3*(1-pe*T)*E2*pi^4/ro/l^4/64;
omega=sqrt((1-pe*T)*E0*pi^2/ro/l^2/4);
F=2*P/ro/A;
%w1=1/3*omega-0.2;
x0_1=0.00001;x0_2=0.00001;
w=0.3:0.005:2;
hold on
for n=1:length(w)
T=2*pi/w(n);
t_final=500*T;
x0=[x0_1;x0_2];
% options=odeset('RelTol',1e-3,'AbsTol',1e-4);
[t,x]=ode45('poin',[10:T/10:t_final],x0,[],w(n));
plot(w(n),x(10:10:end,1),'k.','markersize',5),hold on;
end
grid
%xlabel('轴向力 Excitation P/N');% ,'fontsize',14);
ylabel('位移 Displacement a/m');
%ylabel('速度 Speed V/m/s');
xlabel('频率 Frequency \Omega/rad/s');
%ylabel(' 位移 Displacement a/m'); |