CODEWITHSUNDEEP
matlabdifferential-equationsode
Oliver Amundsen
Oliver Amundsen

Reputation: 1511

how to solve a system of Ordinary Differential Equations (ODE's) in Matlab

I have to solve a system of ordinary differential equations of the form:

dx/ds = 1/x * [y* (g + s/y) - a*x*f(x^2,y^2)]
dy/ds = 1/x * [-y * (b + y) * f()] - y/s - c

where x, and y are the variables I need to find out, and s is the independent variable; the rest are constants. I've tried to solve this with ode45 with no success so far:

y = ode45(@yprime, s, [1 1]);

function dyds = yprime(s,y)
    global g a v0 d
    dyds_1 = 1./y(1) .*(y(2) .* (g + s ./ y(2)) - a .* y(1) .* sqrt(y(1).^2 + (v0 + y(2)).^2));
    dyds_2 = - (y(2) .* (v0 + y(2)) .* sqrt(y(1).^2 + (v0 + y(2)).^2))./y(1) - y(2)./s - d;
   dyds = [dyds_1; dyds_2];
return

where @yprime has the system of equations. I get the following error message:

YPRIME returns a vector of length 0, but the length of initial conditions vector is 2. The vector returned by YPRIME and the initial conditions vector must have the same number of elements.

Any ideas? thanks

Upvotes: 0

Views: 3979

Answers (2)

WYSIWYG
WYSIWYG

Reputation: 494

the syntax for ODE solvers is [s y]=ode45(@yprime, [1 10], [2 2])

and you dont need to do elementwise operation in your case i.e. instead of .* just use *y(:,1) plotted against s with initial conditions [2 2] and parameter values: g=1;a=2;v0=1;d=1.5;

y(:,2) vs s

Upvotes: 0

s.bandara
s.bandara

Reputation: 5664

Certainly, you should have a look at your function yprime. Using some simple model that shares the number of differential state variables with your problem, have a look at this example.

function dyds = yprime(s, y)
    dyds = zeros(2, 1);
    dyds(1) = y(1) + y(2);
    dyds(2) = 0.5 * y(1);
end

yprime must return a column vector that holds the values of the two right hand sides. The input argument s can be ignored because your model is time-independent. The example you show is somewhat difficult in that it is not of the form dy/dt = f(t, y). You will have to rearrange your equations as a first step. It will help to rename x into y(1) and y into y(2).

Also, are you sure that your global g a v0 d are not empty? If any one of those variables remains uninitialized, you will be multiplying state variables with an empty matrix, eventually resulting in an empty vector dyds being returned. This can be tested with

assert(~isempty(v0), 'v0 not initialized');

in yprime, or you could employ a debugging breakpoint.

Upvotes: 2

Related Questions