麻烦大家帮忙看看这个振动微分方程，不知道如何去解......
大家好，我的这个振动微分微分方程里面由于质量是耦合的，在用ode求解时，由于赋初值时只有y0，再代入方程迭代时，就算不了，不知道ode是不不能解这类方程，还请大师指导！谢谢大家，以下是方程及程序代码，期待得到大家的指导，真心求教......函数程序：
function dy=xuanjiafun(t,y)
m=1380;I=2440;J=380;e=0;i=0;l1=1.25;l2=1.51;
b1=0.74;b2=0.74;k11=17;k21=17;k31=17;k41=17;c1=1.5;c2=1.5;c3=1.5;c4=1.5;A=0.05;w=8;t1=0.13;
dy(1)=y(2);
dy(3)=y(4);
dy(5)=y(6);
dy(2)=(-(k11*(y(1)-l1*y(3)-b1*y(5))+k21*(y(1)-l1*y(3)+b2*y(5))+k31*(y(1)+l2*y(3)-b1*y(5))+k41*(y(1)+l2*y(3)+b2*y(5))+c1*(dy(1)-l1*dy(3)-b1*dy(5))+c2*(dy(1)-l1*dy(3)+b2*dy(5))+c3*(dy(1)+l2*dy(3)-b1*dy(5))+c4*(dy(1)+l2*dy(3)+b2*dy(5)))+2*A*sin(w*t)+2*A*sin(w*(t+t1)))/m-e*dy(4)-i*dy(6);
dy(4)=((-((b1^2+l1^2)^1/2+(b2^2+l1^2)^1/2)*A*sin(w*t)+((b2^2+l2^2)^1/2+(b1^2+l2^2)^1/2)*A*sin(w*(t+t1))-(l2*(k31*(y(1)+l2*y(3)-b1*y(5))+k41*(y(1)+l2*y(3)+b2*y(5)))-l1*(k11*(y(1)-l1*y(3)-b1*y(5))+k21*(y(1)-l1*y(3)+b2*y(5)))-l1*(c1*(dy(1)-l1*dy(3)-b1*dy(5))+c2*(dy(1)-l1*dy(3)+b2*dy(5)))+l2*(c3*(dy(1)+l2*dy(3)-b1*dy(5))+c4*(dy(1)+l2*dy(3)+b2*dy(5))))-I*dy(4))/(e*m)-dy(2)-i*dy(6))/e;
dy(6)=(((b2^2+l1^2)^1/2-(b1^2+l1^2)^1/2)*A*sin(w*t)+((b2^2+l2^2)^1/2-(b1^2+l2^2)^1/2)*A*sin(w*(t+t1))-(b2*(k21*(y(1)-l1*y(3)+b2*y(5))+k41*(y(1)+l2*y(3)+b2*y(5)))-b1*((k11*(y(1)-l1*y(3)-b1*y(5))+k31*(y(1)+l2*y(3)-b1*y(5)))-b1*(c1*(dy(1)-l1*dy(3)-b1*dy(5))+c3*(dy(1)+l2*dy(3)-b1*dy(5)))+b2*(c2*(dy(1)-l1*dy(3)+b2*dy(5))+c4*(dy(1)+l2*dy(3)+b2*dy(5))))-J*dy(6))/(i*m)-dy(2)+e*dy(4))/i;
t(end,1)
dy=;
主程序：
clc
clear
hold on;
box on;
Y1=[];
w=8;k=0;w_vc =w/(2*pi);period=1/w_vc;
step=period/100;
tspan=0:step:500*period;
y0=;
options=odeset('RelTol',1.0e-9); % 10-9 for ode15s; 10-5 for ode45
=ode45('xuanjiafun',tspan,y0,options); %计算积分
%=ode45('bevelchaos_fun',tspan,y0,options,b); %计算积分
Y100=;
Y200=;
Y300=;
Y400=;
Y500=;
Y600=;
save Y100.dat Y100 -ascii
save Y200.dat Y200 -ascii
save Y300.dat Y300 -ascii
save Y400.dat Y400 -ascii
save Y500.dat Y500 -ascii
save Y600.dat Y600 -ascii
y=
save datay.dat y -ascii
同求同求！质量耦合还能不能用龙格库塔法？？求大神 希望你用公式说明，或者文字。。看程序很累 dw04116 发表于 2013-3-13 11:39 static/image/common/back.gif
希望你用公式说明，或者文字。。看程序很累
我将没用的删掉，方程如下：
function dy=fun(t,y)
dy(1)=y(2);
dy(3)=y(4);
dy(5)=y(6);
dy(2)=。。。。。。-e*dy(4)-i*dy(6);
dy(4)=。。。。。。。。
dy(6)=。。。。。。。。
t(end,1)
dy=
在用龙格库塔计算时，赋初值是y0=，在算dy(2)的时候出现了dy(4)和dy(6)由于没有初值算不下去，还是我在龙格库塔解法的解法上理解有问题？ 问题解决啦！小小的问题......还是要认真学习，哈哈...... 冬夜 发表于 2013-3-14 11:03 static/image/common/back.gif
问题解决啦！小小的问题......还是要认真学习，哈哈......
怎么解决的？求指导！谢谢
页:
[1]