在看MATLAB中ode45的应用
各位高手:在应用matlab中ODE45求解二阶微分方程时,如何输出状态方程本身随时间的时间历程。比如方程是二阶常微分方程,应用ODE45可以求出X和X的一次导数,如何获得X的二次导数呢?请高手们指导? 解析解:time=linspace(0,4,2000);
x=dsolve('5*D2x+2000*x+0.5*5*9.81*Dx=100*sin(30*t)','Dx(0)=0.1','x(0)=0.1');
xdot=diff(x);
x2dot=diff(xdot);
数值方法:
function solve_eq_examp
global k acc;
k=0; % 该计数器输出状态方程被调用的次数
acc=[]; %输出状态方程的时间历程
=ode45(@Slide_equ,,);
子函数
function xdot=Slide_equ(t,x)
global k acc;
k=k+1
xdot=[x(2)
-2000/5*x(1)-0.5*9.8*sign(x(2))+100/5*sin(30*t)];
acc=]; %这里输出状态方程的时间历程
数值方法中该计数器最后为3385,表明状态方程被调用了3385次。在应用全局变量返回状态方程每次被调用的是的结果后,得到的x二阶导数时间历程和解析解差距较大?
全局变量acc返回的矩阵大小为3385向,要不ODE放回的 大很多??
请问高手,有什么好的解决方法?不胜感激? 怎么一直没有人关注啊!:@L 雨人 发表于 2011-3-26 12:33 static/image/common/back.gif
数值方法中该计数器最后为3385,表明状态方程被调用了3385次。在应用全局变量返回状态方程每次被调用的是的 ...
说明一下你是如何求二阶导数的?怀疑可能是数值求导方面的问题
最好把代码给全,并附上结果图 {:{39}:}! 请问X何X一阶导数和解析值相差大吗?ODE45也不是所有的微分方程都是用! 回复 8 # meiyongyuandeze 的帖子
ode45是用来求解二阶微分方程,具体的实现是通过变量代换,将二价的转换为一阶的进行求解的,方程的输出只有X和X的倒数,你如果想得到X的二阶倒数,你可以根据方程本身,将X的二阶倒放到一边,其他的放到另一边,直接代入数值就可以得到X的二阶倒,也不知道这样可不可行。 E:\matlab\Steerability\单自由度方程的不同解法、解析数值解对比.bmp
感谢大家的回复,以上代码是完整的,可以使用,这里我把结果贴上,以供大家讨论。 E:\matlab\Steerability\单自由度方程的不同解法\解析数值解对比.bmp
我上传不了图片!:@L
对于之前发的帖子,现在有点结果了。做这个简单的方程是为了验证一下方法的可行性。由于ODE45在积分时,多次调用了被积函数,直接使用全局变量返回被积函数的每次迭代值是不准确的。图片上的加速度项可以看出误差的波动项。
之前之所以想使用二阶微分项的原因是,由于动力系统各个自由度之间存在着惯性耦合,(同一个微分方程中包含着其它自由度上微分方程变量的二阶微分项),所以进行数值求解时遇到了困难,如果使用SIMULINK去左,出现代数环,同样很难解决。现在的决定从方程本身着手,手动解除微分方程组的二阶导数之间的耦合项,一遍求解方便。
再次谢谢大家的支持 回复 6 # gghhjj 的帖子
图片上传出错了,程序是完整的,感谢参与讨论 关注,希望能分享的进展! 回复 13 # meiyongyuandeze 的帖子
个人感觉使用ODE系列函数进行数值求解,也就是对状态方程的积分,获得各个状态变量的时间序列。由于ODE系列方程只能返回各个状态变量积分后的时间序列,我的问题有什么方法可以把各个时间点对应的状态也输出出来。这个问题我一直没有解决。考虑过全局变量,但是结果不理想。也使用过用积分后的状态变量重新计算带入状态方程,计算各个时间点的状态。感觉有点麻烦。请问各位有什么好的方法。
http://forum.vibunion.com/data/attachment/album/201104/01/145100k7jmk9k6k826im2m.jpg
页:
[1]
2