Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

ode15i for second order implicit ODE

888 views
Skip to first unread message

EL

unread,
Mar 16, 2011, 11:18:04 AM3/16/11
to
Hello all,

I have a implicit ODE System (two large coupled equations) of second-order:

f(x,x',x'') = 0

and need to solve them for some initial values x(0) = x0, x'(0) =x'0. Matlab provides the numerical solver "ode15i" for implicit ODEs, but it is only capabale to solve first-order ones:

f(x,x') = 0.

I am familiar of reducing higher-order explicit ODEs to a system of first-order ODEs and solve them with i.e. ode45. Now, my ODE system is very large and I do not want to convert it to explicit form x'' = f(x,x',t). I am not sure how to write a function for the implicit case, so that ode15i is able to work with it. I tried to do this for the following simple example first, but I got stuck.

2nd-order "implicit" ODE:
m x'' - kx = 0

Transformation to 1st-order (?!) :
(1) x'- v = 0
(2) m v' - kx =0

Now, how do function code and the ode15i-call have to look?
numsolution = ode15i(@myode, tspan, ..., ..., options);

function f = myode(t, q, dq)
f = [..., ...];
end

Thanks you.

Torsten

unread,
Mar 16, 2011, 11:54:36 AM3/16/11
to

f = [dq(1)-q(2),m*dq(2)-k*q(1)];

Best wishes
Torsten.

EL

unread,
Mar 16, 2011, 1:15:22 PM3/16/11
to
Torsten <Torsten...@umsicht.fraunhofer.de> wrote in message <55f7d63d-e413-47fc...@dr5g2000vbb.googlegroups.com>...


Thank you for your quick reply. Your solution works.
The call should look like:

ode15i(@myode, tspan, [x0, dx0]' , [dx0, 0], options);

where the last value dx0 of the first initial conditions vector matches the first entry of the second initial conditions vector. The last entry should be zero. At least in most cases, I guess.

Btw, I intended to use the one-mass-oscillator ODE for the test example above. There is a sign error, which has kept me to get the right state plots. The correct equation is:

m x'' + kx = 0

which leads to: f = [dq(1)-q(2),m*dq(2) + k*q(1)];

Thanks again!

Regards.

sijo...@gmail.com

unread,
Jan 2, 2019, 12:10:52 PM1/2/19
to
I tried this code, but I am getting an error called MYODE must return a column vector.

Could you please help me with this ?
0 new messages