圆盘非线性微分方程组使用RK4计算的结果
这样的方程难道MATLAB也无能?方程组见附件1。
附件2,3是我写的MATLAB程序。使用MATLAB的ode45函数求解此方程组,在我的机器上运行长达4小时,始终不能计算出结果。
检查程序,应该没有什么问题。欢迎各位高手和群里的同行指点。
我的QQ:635359794
E-mail:kpg19850521@126.com
[ 本帖最后由 kpgkpg 于 2008-11-28 16:24 编辑 ] main
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%系统方程(圆盘不平衡)
% x'' + 2ξx' + x +γ( x.^2 + y.^2) x = e*ω_bar.^2cosτ
% y'' + 2ξy' + y +γ( x.^2 + y.^2) y = e*ω_bar.^2sinτ - G
% G =g/(ω.^2*δ) ;τ=ω*t
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%初始参数
% 2ξ = 0. 24
% ω_bar = 2. 1
% γ = 0. 5
%偏心量e :
%e = 0. 1175、0. 1384 和0. 1750
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
keci=0.12;
omega_bar=2.1;
gama=0.5;
delta=0.1;
omega=2.0;
G=9.8/(omega.^2*delta);
% e=0.08:0.01:0.18;
t=0:500;
tao=omega*t;
initial=;
cf=[];
ccf=[];
e=[ 0.1175, 0.1384 ,0.1750];
for i=1:length(e)
% if e(i)==0.1175 || e(i)==0.1384||e(i)==0. 1750
coef=;
% for t=0:0.1:500;%500个周期(omega)
=ode45('right',tao,initial,[],coef);
cf=;
% ccf=;
% end
end
figure
title('PointCare 图');
plot(cf(2,:),f(1,:))
xlabel('位移')
ylabel('速度')
function f=right(t,x,flag,coef)
f=[x(2);
-x(2)*coef(1)-coef(2)*x(1)-coef(3)*(x(1).^2+x(3).^2)*x(1)+coef(4)*cos(t);
x(4);
-x(4)*coef(1)-coef(2)*x(3)-coef(3)*(x(1).^2+x(3).^2)*x(1)+coef(4)*sin(t)-coef(5)];
[ 本帖最后由 无水1324 于 2008-11-28 16:55 编辑 ] 程序值得稍加修改,问题取得初步进展。修改的程序,画出的图形参见附件。
图片总是太大,裁剪了好几次,看来只能上传这一张了,其他的大家运行下MATLAB自己可以画出来的
[ 本帖最后由 kpgkpg 于 2008-11-28 19:06 编辑 ] 无水兄,请问怎么样把图片贴出来呢?
回复 板凳 kpgkpg 的帖子
下面有一个“发表新回复”你点那个之后就知道了有一个增加附件的就可以了轴心轨迹图
轴心轨迹图学习学习
多谢二位的奉献回复 6楼 kpgkpg 的帖子
是不是你的步长太大了,这个轴心轨迹很粗糙的其他图形如何作出
文章中的其他图形又是如何作出来呢,比如波形图,最大lyapunow指数图,相平面图,轴心轨迹图,时域波形图和幅值谱图等。回复 楼主 kpgkpg 的帖子
文章中的其他图形又是如何作出来呢,比如波形图,最大lyapunow指数图,相平面图,轴心轨迹图,时域波形图和幅值谱图等。 学习了~~~~~ t的步长取小一些,试试0.001. 无水1324 发表于 2008-11-28 16:52 static/image/common/back.gifmain
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%系统方程(圆盘不平衡)
请问下,就在这个模型,如果我想画无量纲振动随转速变化的分岔图,要求转速是50:1900,怎么修改下程序呢,谢谢 回复 13 # fanhuasijin 的帖子
那要看你的模型了,参数,初值,方程,取值等等都可能要修改 回复 14 # hsfy919 的帖子
Jeffcott刚性转子轴承系统的动力学模型。其方程如下:
mx’’=-fx+mew2sin(wt)
my’’=-fy+mew2cos(wt)+mg
上采用无量纲处理,无量纲轴颈坐标:X=x/c,Y=y/c;X’=x’/(cw), Y’=y’/(cw^2);X’’=x’’/(cw), Y’’=y’’/(cw^2);无量纲时间:tau=wt;无量纲油膜力分量:Fx=fx/…, Fy=fy/…其余无量纲量(略).
设z1=X,Z2=Y,Z3=X’,Z4=Y’. 其中“’”表示d/dtau, 则(式-1)用状态变量表示的无量纲形式为z1’= Z3,z2’=Z4,z3’=-Fx/M+p*sin(tau),z4=-Fy/M+p*cos(tau)+G1
我要画出无量纲振动随转速变化的分岔图,那是不是就意味着要有时间t和转速w这两个自变量了,如果有两个自变量还能算吗,如果可以的话可不可以指定些参考资料,
问题比较大,希望先给我指个方向,我再琢磨琢磨,谢谢
页:
[1]
2