|
马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
這是我已經做好的函數
有11各微分方程,我是做機械的齒輪系統分析,分析其震動位移與加速度的情況
初始條件都是[0],時間從0~1000
function y=ivp4(t,y)
yp=zeros(0); % note: column vector
yp=yp';
yp(1)=y(2);
yp(2)=(6.15*10.^6*(y(3)-y(1)))/13.2;
yp(3)=y(4);
yp(4)=((8.1*10.^7*(y(5)-y(3)))-(6.15*10.^6*(y(3)-y(1))))/(1.3);
yp(5)=y(6);
yp(6)=((-6.03*10.^6*(y(5)-y(3)))+(1.545*10.^4*0.176*(0.176*y(7)-0.176*y(5))))/0.76;
yp(7)=y(8);
yp(8)=((1.545*10.^4*0.176*(0.176*y(7)-0.176*y(5)))-(1.545*10.^4*0.176*(0.176*y(9)-0.176*y(7))))/0.76;
yp(9)=y(10);
yp(10)=((-1.545*10.^4*0.176*(0.176*y(9)-0.176*y(7)))+(6.03*10.^6*(y(11)-y(9))))/0.76;
yp(11)=y(12);
yp(12)=((-6.03*10.^6*(y(11)-y(9)))+(1.618*10.^4*0.150*(0.150*y(13)-0.150*y(11))))/0.56;
yp(13)=y(14);
yp(14)=((-1.618*10.^4*0.150*(0.150*y(12)-0.150*y(11)))+(2.01*10.^6*(y(15)-y(13))))/0.56;
yp(15)=y(16);
yp(16)=(-2.01*10.^6*(y(15)-y(13)))/3.66-500;
yp(17)=y(18);
yp(18)=(1.159*10.^4*0.156*(0.084*y(19)-0.156*y(17)))/0.25;
yp(19)=y(20);
yp(20)=(2.01*10.^6*(y(21)-y(19))-1.159*10.^4*0.156*(0.084*y(19)-0.156*y(17)))/0.03;
yp(21)=y(22);
yp(22)=(-2.01*10.^6*(y(21)-y(19)))/1.5-3000;
執行程式碼
y0=[0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0; 0] ;
>> tspan = [0 1000];
>> [t, y] = ode45('ivp4', tspan, y0);
>> plot(t,y(:,2))
可是當我執行ODE45之後,其值就發散,不知道是那出現問題了
請大俠邦幫忙
[ 本帖最后由 eatnche594 于 2007-4-13 09:59 编辑 ] |
|