软弹簧duffing方程不能用ode45么？

lihaitao123

这是the duffing equations：nonlinear oscillator and their behaviors 里的原图。
但是这种soften duffing 方程我用matlab 中的ode45（自适应步长）和ode4（定步长）算出
的结果全是发散的。

描述以及方程

F=0.35，388，393

分岔以及局部放大

lihaitao123

clear all
%数据的初始化
global  gamma omega f
gamma=0.2;
omega=0.8;
f=.35;
t0=0;
tfinal=40;
tspan=t0:0.01:tfinal;
y0=[0,0,0];
[t,y]=ode45(@softduffing,tspan,y0);
figure;
plot(t,y(:,1))
xlabel('t'),ylabel('y(1)'),title('时间历程图')%时间历程图
figure
plot(y(end-5000:end,1),y(end-5000:end,2))%像轨迹
xlabel('y(1)'),ylabel('y(2)')

function dy=softduffing(t,y)
global gamma omega f
dy=zeros(3,1);
dy(1)=y(2);
dy(2)=-2*gamma*y(2)-y(1)+y(1)^3+f*cos(omega*y(3));
dy(3)=1;

gghhjj

lihaitao123 发表于 2014-4-14 14:41
clear all
%数据的初始化
global  gamma omega f

尝试了一下你的程序，从程序上来看可能是方程奇异性造成的

如果方程中三阶项减小点，比如：

1. dy(1)=y(2);
   
2. dy(2)=-2*gamma*y(2)-y(1)+0.8*y(1)^3+f*cos(omega*y(3));
   
3. dy(3)=1;

复制代码

则计算没有什么问题，因此建议换换算法看看，比如用差分法，newmark法等试试看

zzqmessi

lihaitao123 发表于 2014-4-14 14:41
clear all
%数据的初始化
global  gamma omega f

你也可以用ode15s 试试

zzqmessi

lihaitao123 发表于 2014-4-14 14:41
clear all
%数据的初始化
global  gamma omega f

确实很怪异 你把cos换成sin

gghhjj

zzqmessi 发表于 2014-4-17 16:50
你也可以用ode15s 试试

这个方程并不是刚性的

gghhjj

zzqmessi 发表于 2014-4-17 17:16
确实很怪异 你把cos换成sin

确实邪门，呵呵

zzqmessi

gghhjj 发表于 2014-4-18 08:37
确实邪门，呵呵

我还是想问问你 有没有经典解决 duffing方程 这样骨干曲线的 用数值法或者 解析法

lihaitao123

ode15s 和 ode23 等等都是不行的！

zzqmessi

lihaitao123 发表于 2014-4-18 11:25
ode15s 和 ode23 等等都是不行的！

cos 改为 sin 就可以  我想问问你 你有没有做duffing方程的幅频特性曲线呀

lihaitao123

cos 改为 sin 频率就会发生变化了吧，系统应该变化了吧。没有，论坛里有很多duffing分叉图的应该可以借鉴一下思想，（频率作为分叉参数）吧

lihaitao123

这个系统无论sin还是cos如果做分岔图的话，都会碰到解不收敛的情况。。

lihaitao123

图可以调出来了，不是诡异，非线性系统本身具有初值敏感的特点。初值奉上（0，-0.5），（0,0.05），（0.11,0）。还是自已以前学业不精呀。

gghhjj

lihaitao123 发表于 2014-4-18 15:28
图可以调出来了，不是诡异，非线性系统本身具有初值敏感的特点。初值奉上（0，-0.5），（0,0.05），（0.11, ...

谢谢分享经验！

车萌

你好，楼主，我做了一下你这个软刚度非线性，用的是和你这个一样的数据，带入了你说的初值，但是还是出现警告，我把我的程序复制上，楼主你给帮忙看一下和你的区别。