MATLAB答案

runge kutta第四阶方法

17次观看（最近30天）
显示较旧的评论
瑞安·鲍曼（Ryan Bowman）
  • 0
回答： Meysam Mahooti 2021年5月5日
我已经开发了第四阶runge kutta方法，可以帮助我找到刚体的角速度。它应显示为3x1矩阵，以包括X，Y和Z速度。我的K2，K3和K4是不可分割的，是3x3矩阵，而我的K1是3x1矩阵。我相信所有矩阵都应该是3x1。这是我的代码：（我对我的函数进行了单独的文件。）
Clc
清除全部
清除变量
格式短e
%%问题1
i = [5.3 0 0; 0 7.4 0; 0 0 10.5];
W_0 = [0.1 0.2 0.3];
t = [0 0 0]';在这种情况下，％没有扭矩
DT = 0.025;％秒
t = 0：dt：10;％秒
[wdot] = getwdot（w_0，t，i）;
w（1，:) = w_0;
为了II = 1：长度（T）-1
k1 = dt*getwdot（w（ii，:)，t，i）;
k2 = dt*getwdot（w（ii，：）+0.5*k1，t+0.5*dt，i）;
k3 = dt*getwdot（w（ii，：）+0.5*k2，t+0.5*dt，i）;
k4 = dt*getwdot（w（ii，：）+k3，t+dt，i）;
w（ii+1，:) = w（ii，:)+（1/6）*（k1+2*k2+2*k3+k4）
结尾
功能[wdot] = getwdot（w_0，t，i）
imat = diag（i）;
imat_inv = diag（1./i）;
wdot = imat_inv。*cross（w_0'，imat。*w_0'）
任何解决此问题的建议都将不胜感激。

登录来回答这个问题。

答案（3）

吉姆·里格斯（Jim Riggs）
编辑：吉姆·里格斯（Jim Riggs） 2020年3月5日
我认为问题是T论点。
在K1中，您有一个资本t，这是一个3矢量
在K2，K3和K4中，这是一个较低的情况t，即0：dt：10。
那里有问题。
k1 = dt*getwdot（w（ii，:)，t，i）;
k2 = dt*getwdot（w（ii，：）+0.5*k1，t+0.5*dt，i）;
k3 = dt*getwdot（w（ii，：）+0.5*k2，t+0.5*dt，i）;
k4 = dt*getwdot（w（ii，：）+k3，t+dt，i）;
他们可能都是错误的。我认为您想要的是t（ii）（尽管我看不到t在功能中使用）
1条评论
吉姆·里格斯（Jim Riggs）
吉姆·里格斯（Jim Riggs） 2020年3月5日
另外，在您的功能中
imat_inv = diag（1./i）;
内部操作
（1./i）
在OFF对角项中产生3x3矩阵WTH除以零。
您最好使用
imat_inv = 1./imat;
由于IMAT已经是3个向量。该操作中没有零的分歧，您最终会得到相同的结果。
否则，您应该使用
imat_inv = diag（Inv（i））;
这会产生一个propper矩阵逆，您可以从中选择对角线。

登录发表评论。

詹姆斯·图萨（James Tursa）
编辑：詹姆斯·图萨（James Tursa） 2020年3月5日
在求解W_DOT的EOM时，您应该得到一个减号：
https://en.wikipedia.org/wiki/euler%27s_equations_%28rigid_body_dynamics%29
即，这个
wdot = imat_inv。*cross（w_0'，imat。*w_0'）
应该是这样：
wdot = -imat_inv。*cross（w_0'，imat。*w_0'）
另外，我想您正在尝试巧妙地计算惯性矩阵逆，并与所有这些诊断物相乘。何苦？只需使用BackSlash：
功能[wdot] = getwdot（w_0，t，i）
wdot = -i \ cross（w_0'，i*w_0'）
而且，您在上面的第三个论点中不一致。对于K1，您可以通过t，但是在K2，K3，K4中，您可以通过整个时间向量传递。您不会在派生型中使用此参数，因此不会引起任何问题，但是您的调用语句不正确……您可能是要将t（ii）而不是t用于该参数。
最后，尽管不是必需的，但我会选择让所有向量为列向量，因此您不需要将所有这些转置进行计算。但这只是一个nit。
Meysam Mahooti
Meysam Mahooti 2021年5月5日
//www.tatmou.com/matlabcentral/fileexchange/55430-runge-kutta-4th-core?s_tid=srchtitle

类别

标签

下载188bet金宝搏


释放

R2019A

社区寻宝

在Matlab Central中找到宝藏，发现社区如何为您提供帮助！

开始狩猎！