煜宸0922 发表于 2011-6-8 09:07

向各位大侠求助庞加莱投影图的画法

各位大侠,小弟刚入非线性不久,做两自由度的间隙碰撞系统,相图已经画出,但庞加莱投影图十分不会,好长时间没进展,请各位大侠指教,不胜感激。

kangarooli 发表于 2011-6-8 09:35

回复 1 # 煜宸0922 的帖子

方便的话你可以把系统发上来看一下,那样可能会更针对些,poincare图我现在参考的是论坛里频闪法思想,你可以搜搜看,好多的,看看有没有适合你的

煜宸0922 发表于 2011-6-8 09:41

回复 2 # kangarooli 的帖子

我要用的是截面法,系统的程序三部分:
1,系统方程的程序:
function dxdt=fun1(t,x)
global w
q1=0;
q2=w;   
q3=0.05;
q4=0;
q5=4;
q6=4;
q7=3;
dxdt=zeros(4,1);
dxdt(1)=x(3);
dxdt(2)=x(4);
dxdt(3)=(1-q1)*sin(q2*t+q4)-2*q3*x(3)+2*q3*x(4)-x(1)+x(2);
dxdt(4)=(q1*sin(q2*t+q4)+2*q3*x(3)-2*q3*(1+q5)*x(4)+x(1)-(1+q6)*x(2))/q7;
end
2.龙格库塔法:
functionx_new=rk144(x,h,T)
k1=h*feval('fun1',T,x(:,1));
k2=h*feval('fun1',T+h/2,x(:,1)+k1/2);
k3=h*feval('fun1',T+h/2,x(:,1)+k2/2);
k4=h*feval('fun1',T+h,x(:,1)+k3);
x_new=x(:,1)+(k1+2*k2+2*k3+k4)/6;
end
3.主程序
function =changerk344(N)
x(:,1)=;      % 初始状态
h=0.01;                % 初始步长,最大步长
AbsTol=1e-5;         % 容差
p=4;                        
T=zeros(1,N);          % 时间向量
T(1,1)=0;            % 时间起点
p21=2^p-1;
V =[];
tic
for i=1:N-1
   if abs(x(1,i))<0.09             % 定步长
      h=0.01;      
      x(:,i+1)=rk144(x(:,i),h,T(1,i));
      T(1,i+1) = T(1,i) + h;          % 更新时间
    else                                 % 变步长
      x1=rk144(x(:,i),h,T(1,i));
      x2=rk144(x(:,i),h/2,T(1,i));
      a=abs(x1-x2)/p21;
      while a(1,1)>AbsTol                % 选择合适步长
            h=h/10;
            x1=rk144(x(:,i),h,T(1,i));
            x2=rk144(x(:,i),h/2,T(1,i));
            a=abs(x1-x2)/p21;
      end
      x(:,i+1)=x1;      % 根据选择的步长更新状态
       if   abs(x(1,i+1))>=0.1
         V=x(3,i+1);
         x(3,i+1)=-0.8*x(3,i+1);
         x(1,i+1)=sign(x(1,i+1))*min(0.1,abs(x(1,i+1)));%限制位置在0.1处
       end
    T(1,i+1) = T(1,i) + h;   
      end
   h=0.01;
end
toc
T = T(:);   % 将时间转换为列向量
x = x(:,1:N);
选的截面是x(1)=0.1,1的速度等于碰后的速度

kangarooli 发表于 2011-6-8 09:46

回复 3 # 煜宸0922 的帖子

截面法没研究过,不过你程序既然有了,那现在问题是什么呢,出来的映射图对应不上还是怎么的

煜宸0922 发表于 2011-6-8 09:53

回复 4 # kangarooli 的帖子

程序只能画出来相图,庞加莱截面图我不会画呀

kangarooli 发表于 2011-6-8 10:18

回复 5 # 煜宸0922 的帖子

我运行了下你的程序老是有错误,看到截面我还以为你上面的程序是画截面图的呢,这么说你是根据荣格库塔思想自己编程序画的相图啊,我看你的模型也不是太复杂,怎么不用matlab自带的ode45什么的解呢,可能我还没看明白你的问题

煜宸0922 发表于 2011-6-8 10:29

回复 6 # kangarooli 的帖子

程序是对的,我这边运行没有问题,那三个程序要分开存储,然后调用。我们要求要用自己编的龙格库塔法。模型就是两自由度含间隙的碰撞振动系统,一个二阶常微分方程组。化成四状态方程。其中x1是质块1的位移,x3是质块1的速度,x2是质块2的位移,x4是它的速度。相图已经可以画出,但我要继续画庞加莱投影图和分岔图,但我试了试都不对,所以才求助

煜宸0922 发表于 2011-6-8 10:31

回复 7 # 煜宸0922 的帖子

主程序中用的变步长,另加碰撞条件

kangarooli 发表于 2011-6-8 10:52

回复 7 # 煜宸0922 的帖子

程序我是分开存储的,但老是提示N未定义,可能是我这里的问题,模型基本都是一个类型的,惭愧啊,以前都用现成的,还没做过变步长的呢,期待高手给你解答

煜宸0922 发表于 2011-6-8 16:40

回复 9 # kangarooli 的帖子

n是点数,可以直接赋值的。

liliangbiao 发表于 2011-6-9 15:44

频闪法不适合做含有间隙的碰撞动力系统,需要在碰撞面上建立poincare section. 建议看看罗冠炜的书籍和论文。

kezairenjian 发表于 2011-6-9 20:58

liliangbiao 发表于 2011-6-9 15:44 static/image/common/back.gif
频闪法不适合做含有间隙的碰撞动力系统,需要在碰撞面上建立poincare section. 建议看看罗冠炜的书籍和论文 ...

关于间隙振动系统的庞加莱截面法,是否可以不建立poincare setion,而直接在相图上找到截面呢?

liliangbiao 发表于 2011-6-9 22:57

建议看看罗冠炜的书籍和论文,我只能告诉你这么多。呵呵。至于poincare section 那是必须的!因为碰撞面就是速度相反的时刻。那么这个时刻就是所有轨线都被Poincare垂直相交的平面,难道你的相图上看不出来吗?

meiyongyuandeze 发表于 2011-6-12 08:32

可以截取碰撞面来最为庞加莱面吧,当初谢建华老师就是这么给我们讲的!

liliangbiao 发表于 2011-6-12 17:16

呵呵,你还上过谢教授的课,那么你应该会了,不妨给他说说!
页: [1] 2
查看完整版本: 向各位大侠求助庞加莱投影图的画法