声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1288|回复: 2

[非线性振动] 相图和幅值谱的问题

[复制链接]
发表于 2008-10-10 17:25 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
这是我做的时间历程,相图和幅值谱,可是前面两个和频谱对不上,从相图可以看到theta的运动是4倍周期运动,可频谱出现了5种频率成分,这应该怎么解释呀,高手帮忙给解决一下。

                               
登录/注册后可看大图



下面是我的程序
function dy=moxing3(t,y,flag)
global c;
theta0=pi*0.1;
g=10;
p=1.29;
T=26000;
ml=1.7;h=10/1000;
a=0.4;
dd=0.028;
mi=dd*h*0.9;
D=dd+h*sin(theta0);
L=244;
U=17;
ky=4*T*3/L;
kt=864;
kx=4*T*3/L;
S=3*mi*dd/2*L*(2/pi)^2;
m=(2/pi)^2*L*3*(ml+mi);
J=(2/pi)^2*L*(ml+mi)*a^2;
xiy=0.04;
xix=0.04;
xit=0.06;
omegay=sqrt(ky/m);
omegax=sqrt(kx/m);
omegat=sqrt(kt/J);
a1=-4.5;
a2=1;
a3=2;
b1=-1;
b2=-1;
b3=1;
c1=-5.2;
c2=3;
c3=1;
theta1=y(5)+theta0-y(2)/U-(dd/2/U)*y(6);
cy=a1*theta1+a2*theta1^2+a3*theta1^3;
cx=b1*theta1+b2*theta1^2+b3*theta1^3;
cm=c1*theta1+c2*theta1^2+c3*theta1^3;
Fy=1/2*p*U^2*D*cy*L;
Fx=1/2*p*U^2*D*cx*L;
M=1/2*p*U^2*D^2*cm*L;
func1=y(2);
func2=-(-(2*m*xiy*omegay*y(2)+ky*y(1)-Fy)*S^2*sin(theta0)^2+S^2*cos(theta0)*sin(theta0)*(2*m*xix*omegax*y(4)+kx*y(3)-Fx)-S*cos(theta0)*(2*J*xit*omegat*y(6)+kt*y(5)-M)*m+(2*m*xiy*omegay*y(2)+ky*y(1)-Fy)*J*m)/(-S^2+J*m)/m;
func3=y(4);
func4=(-S^2*sin(theta0)*cos(theta0)*(2*m*xiy*omegay*y(2)+ky*y(1)-Fy)-S^2*sin(theta0)^2*(2*m*xix*omegax*y(4)+kx*y(3)-Fx)+S*sin(theta0)*(2*J*xit*omegat*y(6)+kt*y(5)-M)*m+(2*m*xix*omegax*y(4)+kx*y(3)-Fx)*S^2-(2*m*xix*omegax*y(4)+kx*y(3)-Fx)*J*m)/(-S^2+J*m)/m;
func5=y(6);
func6=-(-S*cos(theta0)*(2*m*xiy*omegay*y(2)+ky*y(1)-Fy)-S*sin(theta0)*(2*m*xix*omegax*y(4)+kx*y(3)-Fx)+(2*J*xit*omegat*y(6)+kt*y(5)-M)*m)/(-S^2+J*m);
dy=[func1;func2;func3;func4;func5;func6];



clear;
step=0.01;
tstop=1000;
y0=rand(1,6);
pt=500*1/step;
[t,y]=ode45('moxing3',[0:step:tstop] ,y0,[]);
set(gcf,'color','w')
subplot(2,1,1);
plot(t,y(:,5));
xlabel('t')
  ylabel('theta')
  axis([900 1000 -Inf Inf])
  grid;
subplot(2,2,3);
plot(y(end-pt:end,5),y(end-pt:end,6));
xlabel('theta')
  ylabel('dtheta')
  grid;
  subplot(2,2,4);
yy=fft(y(end-pt:end,5));
N=length(yy);
power=abs(yy);
freq=(1:N-1)*1/step/N;
plot(freq(1:N/2),power(1:N/2));
xlabel('f(theta)')
  ylabel('')
  axis([0 1.5 0 Inf])
回复
分享到:

使用道具 举报

发表于 2008-10-11 21:06 | 显示全部楼层

回复 楼主 zhanghb 的帖子

这个帖子你在数学版发了,那里面已经有回复了
发表于 2008-10-15 08:03 | 显示全部楼层
从你的频率上看,都是倍数关系,所以个人认为相图看不出来其周期数
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-18 21:45 , Processed in 0.059990 second(s), 17 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表