sssssxxxxx921 发表于 2007-9-18 10:38

那你怎么看出两周期的?我实在看不出怎么 会有个两周期

sssssxxxxx921 发表于 2007-9-18 10:41

function dq=Rotors_System_Sub_Func(t,q)
global BN Nb1 Nb2 w1 w2 Ro1 Ri1 Ro2 Ri2
%-------------------------------------------------------------------------
r01=0.00002;
m1=7.86;
m2=11.93;
W1=m1*9.8;
W2=m2*9.8;
Kb1=4.6346e10;
Kb2=3.4953e9;
%-------------------------------------------------------------------------
F=1485;
F1 = 1485/2;
delta_i=1.2443*w2^2*1e-11;
delta_e=2.2183*w1^2*1e-11;
Dm=(117+133.3)/2;
wc=(w2*(1-8/Dm)+(1+8/Dm)*w1)/2;
delta_o=4.912*wc^2*1e-14;
delta_F=8.10*F1^0.925/7.8^0.85*1e-8;
r02=0.00005-(delta_i-delta_o-delta_e-delta_F);
%-------------------------------------------------------------------------
Fx11=0;
Fy11=0;   
Fx22=0;
Fy22=0;
for i=1:Nb1
   
    sita(i)=2*pi/Nb1*(i-1)+BN/Nb1*w1*2*pi*t;      
    Deformation(i,1)=1e-6*q(1,1)*cos(sita(i))+1e-6*q(2,1)*sin(sita(i))-r01;
    if Deformation(i)<=0
      Deformation(i)=0;
    end
   
    fx11=Kb1*(Deformation(i))^1.5*cos(sita(i));
    fy11=Kb1*(Deformation(i))^1.5*sin(sita(i));      
   
    Fx11=Fx11+fx11;   
    Fy11=Fy11+fy11;
end
Fx12=Fx11;
Fy12=Fy11;
%-------------------------------------------------------------------------
Fx21=0;
Fy21=0;
for j=1:Nb2
    sita(j)=2*pi/Nb2*(j-1)+(Ro2/(Ro2+Ri2)*w1+Ri2/(Ro2+Ri2)*w2)*2*pi*t;
    Deformation1(j,1)=1e-6*q(3)*cos(sita(j))+1e-6*q(4)*sin(sita(j))-r02;
    if Deformation1(j)<=0
       Deformation1(j)=0;
    end   
    fx21=Kb2*(Deformation1(j))^1.1*cos(sita(j));
    fy21=Kb2*(Deformation1(j))^1.1*sin(sita(j));
   
Fx21=Fx21+fx21;
Fy21=Fy21+fy21;
end
%-------------------------------------------------------------------------
for k=1:Nb1
   sita(k)=2*pi/Nb1*(k-1)+BN/Nb1*w2*2*pi*t;
    Deformation2(k,1)=1e-6*q(3,1)*cos(sita(k))+1e-6*q(4,1)*sin(sita(k))-r01;
    if Deformation2(k)<=0
       Deformation2(k)=0;
   end
   
    fx22=Kb1*(Deformation2(k))^1.5*cos(sita(k));
    fy22=Kb1*(Deformation2(k))^1.5*sin(sita(k));
   
    Fx22=Fx22+fx22;
    Fy22=Fy22+fy22;
end
%-------------------------------------------------------------------------
P=1470;
Cx11=1500;Cy11=1500;         
Cx12=2500;Cy12=2500;
Cx21=2000;Cy21=2000;
Cx22=7000;Cy22=7000;
%-------------------------------------------------------------------------
dq(5:8,1)=[-1e6/m1*((Cx11+Cx12)*q(5,1)*1e-6+Fx11+Fx12-(W1+Fx21+Cx21*q(7,1)*1e-6));...
         -1e6/m1*((Cy11+Cy12)*q(6,1)*1e-6+Fy11+Fy12-(Fy21+Cy21*q(8,1)*1e-6));...
         -1e6/m2*((Cx21+Cx22)*q(7,1)*1e-6+Fx21+Fx22-(P+W2));...
         -1e6/m2*((Cy21+Cy22)*q(8,1)*1e-6+Fy21+Fy22)];
dq(1:4,1)=q(5:8,1);         

主程序:
function Rotors_System_Func         
clear
clc
global BN Nb1 Nb2 w1 w2 Ro1 Ri1 Ro2 Ri2
%-------------------------------------------------------------------------
n_one_T=100;         
n_T=400;            
      
Ro1=4.885e-2;      
Ri1=6.865e-2;
Nb1=10;      

Ro2=66.515e-3;
Ri2=58.5e-3;         
Nb2=34;      
w2_min=100;   
w2_max=10400;
w2_step=100;         
w1_rpm=11000;   
q_initial(1:8,1)=1e-11;
BN=Ri1/(Ri1+Ro1)*Nb1;   
%-------------------------------------------------------------------------
tic
    w2_rpm=10000;   
    w2=w2_rpm/60;
    w1=w1_rpm/60;            
    r=Ro2/(Ro2+Ri2);
   
    w_vc_1=w1*BN   
      
    w_cage_3=w1*r+w2*(1-r);      
    w_vc_3=w_cage_3*Nb2
   
   
    w_cage=w2*BN;
    w_vc=w_cage
   
    T_vc=2*pi/w_vc;      
    dt=T_vc/n_one_T;   
    time=n_T*T_vc;         
    n=round(time/dt);      
    t_span(1:n)=linspace(0,time,n);   
   =ode45('Rotors_System_Sub_Func',t_span,q_initial);
%-------------------------以下为绘制系统频谱图程序--------------------------
    xc =q(end-8192:4:end,1);
    yc =q(end-8192:4:end,2);
    xcc=q(end-8192:4:end,3);
    ycc=q(end-8192:4:end,4);
X1 = fft(xc, 2048);
X2 = fft(xcc,2048);
X1(1) = [];
X2(1) = [];
Pxx = X1.* conj(X1) / 2048;
Pyy = X2.* conj(X2) / 2048;
f = 0.25*1/dt*(0:1024)/2048;
subplot(221);plot(t_span(end-500:end),q(end-500:end,1))
xlabel('t (s)')
ylabel('x1 (m)')
subplot(222);plot(f(1:1024),Pxx(1:1024))
xlabel('频率(Hz)')
ylabel('功率')
subplot(223);plot(t_span(end-500:end),q(end-500:end,3))
xlabel('t (s)')
ylabel('x2 (m)')
subplot(224);plot(f(1:1024),Pyy(1:1024))
xlabel('频率(Hz)')
ylabel('功率')
%-------------------------------------------------------------------------
toc

sssssxxxxx921 发表于 2007-9-18 10:42

不好意思    上边没加任何注释

咕噜噜 发表于 2007-9-18 14:47

具体来说,从图一中来看有两个周期,小周期就是你如图2中标的那样,而大周期是你如图1标的A到C

咕噜噜 发表于 2007-9-18 14:47

程序看的有点糊涂,很长,不怎么明白,你能不能说一下你的系统还有你系统类型

sssssxxxxx921 发表于 2007-9-18 16:14

我的是一个转子系统,至于模型不是一两句话能说清的
不过我觉得这对与结果没什么影响吧, 我很奇怪   出来的图这么直观,程序和系统模型难道能说明是几周期的吗?

秋月 发表于 2007-9-18 16:50

回复 #36 sssssxxxxx921 的帖子

你好!请问你有转子系统分岔图的程序吗?
能否给个例子?
谢谢!

sssssxxxxx921 发表于 2007-9-18 17:04

http://forum.vibunion.com/forum/viewthread.php?tid=50058&extra=&page=1
你看看她的吧她也是做转子系统的, 而且最后分岔图也出来了。我的分岔图暂时做不出来

sssssxxxxx921 发表于 2007-9-18 17:11

还有一个关于功率谱的问题,知道的话帮忙看下
http://forum.vibunion.com/forum/thread-51961-1-1.html

秋月 发表于 2007-9-18 18:22

回复 #38 sssssxxxxx921 的帖子

谢谢了!
我刚开始接触还不太懂,谢谢指点.
请多高手多指教!

sssssxxxxx921 发表于 2007-9-18 19:07

回复 #34 咕噜噜 的帖子

另外从相图上看   应该是个极限环才对啊,那是几个闭合的环线,
关于有毛刺现象说明极限环正在变得不稳定,开始有破裂的迹象,这样分析是否合理

咕噜噜 发表于 2007-9-18 19:12

回复 #41 sssssxxxxx921 的帖子

:@L :@L 那个不是极限环,对于极限环可不是这么画的
毛刺也不是极限环破裂造成的

sssssxxxxx921 发表于 2007-9-18 19:32

那极限环是怎么样画的
极限环:运动微分方程的解在相平面上所确定的相轨迹是一条孤立的封闭曲线,它所对应的周期运动由系统的物理参数唯一确定,与初始运动状态无关。这种孤立的封闭相轨迹称为极限环

sssssxxxxx921 发表于 2007-9-18 19:48

看一下这个图,他在判断周期几运动时好像不是像你那样说的A-B 和 A-C两周期来判断的

octopussheng 发表于 2007-9-18 20:54

他这个判断应该是去掉瞬态后确定的吧!

其实方法就是小咕说的那种

小咕,我理解的对不对?
页: 1 2 [3] 4 5
查看完整版本: 帮忙看看这个相图