Error trying to solve 2 Second ODE

3 次查看（过去 30 天）
Scott Angus
Scott Angus 2019-11-15
评论： Scott Angus 2019-11-15
Hi,
I have a function that solves a differential equation using the RK method:
function [tvec, xvec] = rksolve(f, t0, tf, x0, dt)
%Solve f with initial conditions x0 between t0 -> tf
%with increments dt
%Set initial conditions
t = t0;
x = x0;
%Answer output
tvec = t;
xvec = x;
while t < tf
%RK method of solving
k1 = f(x, t);
k2 = f(x+0.5*k1*dt, t+0.5*dt);
k3 = f(x+0.5*k2*dt, t+0.5*dt);
k4 = f(x+k3*dt, t+dt);
x = x + ((k1 + 2*k2 + 2*k3 +k4)*dt)/6;
t = t+dt;
%Store outputs
tvec = [tvec t];
xvec = [xvec x];
end
But when I try to run it with the equations defined in func.m, it comes up with a horzcat error in the line:
xvec = [xvec x];
func.m:
function a = coords(x, t);
%Initial Conditions
xpos = x(1);
vx = x(2);
ypos = x(3);
vy = x(4);
a = [(-xpos - 2*xpos*ypos); (-ypos - (xpos^2) +ypos^2)];
%dxdt = vx
%dvxdt = -x -2*x*y;
%dydt = vy;
%dvydt = -y -(x^2) + y^2;
Any help would be greatly appreciated!
  1 个评论
darova
darova 2019-11-15
What about preallocation?
n = round((tf-t0)/dt+1);
tvec1 = zeros(1,n);
i = 1;
while t < tvec
%% ...
tvec1(i) = t;
i = i + 1;
end
tvec = tvec1;

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

James Tursa
James Tursa 2019-11-15
编辑：James Tursa 2019-11-15
Just looking at func.m, it appears you pass in a 4-element x vector, but you only return a 2-element a vector. You need to return the derivative for all the states. So something like this instead based on your comments:
a = [ vx; (-xpos - 2*xpos*ypos); vy; (-ypos - (xpos^2) +ypos^2) ];
  2 个评论
Scott Angus
Scott Angus 2019-11-15
Thanks for your response. I tried changing the output of func.m to
a = [ vx; (-xpos - 2*xpos*ypos); vy; (-ypos - (xpos^2) +ypos^2) ];
but still come up with an error:
>> [ts, xmat] = rksolve(@func, 0, 200, [0, 0.3, 0, 0], 0.2)
Error using horzcat
Dimensions of arrays being concatenated are not consistent.
Error in rksolve (line 25)
xvec = [xvec x];
I'm not sure why because the x and xvec dimensions should be the same?
Scott Angus
Scott Angus 2019-11-15
I found the issue. I was passing the function rksolve with a horizontal vector for x0 instead of a vertical vector. Adding semicolons between the components of the vector when I passed it solved the problem.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Ordinary Differential Equations 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by