User defined function within ode45

I have the following ODE system:

If f(t) is a constant, then I can easily write the code, but I have f(t) as a piecewise function which I coded up in two ways, using the piecewise command, and a function that has lots of loops.

When I use the piecewise command I get the following error:

Array indices must be positive integers or logical values.

Error in sym/subsref (line 898)

R_tilde = builtin('subsref',L_tilde,Idx);

Any idea what I'm doing wrong?

Answer 1

[Torsten]

Symbolic functions like the one you get if you use the "piecewise" command are not allowed with ODE45.

Either define the right-hand side of your ODE system purely numerical in a function (and f(t) e.g. using if-statements to cope with its piecewise definition) or use matlabFunction to convert a symbolic into a numerical function.

[Walter Roberson]

matlabFunction generating an anonymous function cannot handle piecewise. matlabFunction asked to write to file can handle piecewise but not (usually) vectorization. If you write to file be sure to turn off optimization (unless they have finally fixed it, it has been broken in recent releases)

But piecewise or if/elseif or logical indexing still have the same problem that the derivatives need to be continuous. Any one call to ode45 can only involve one branch of the conditions unless the derivatives have been matched between the conditions. For example no interp1 'linear' but interp1 'spline' is ok

Answer 2

[Walter Roberson]

You have some line in which you have VARIABLE(SymbolicVariable). For example you might have

A = rand(1,5);
syms N
symsum(exp(-A(N).^2), N, 1, 5)

symbolic variables can never be used to for indices.

Note that the mathematics of ode45 requires that the first two derivatives of the function are continuous. If the first derivative is not continuous you might possibly see an error message but not necessarily; violating this condition typically just gets you wrong results.

Most of the time piecewise functions have discontinuous derivatives, unless they have been carefully matched. So most of the time you need to use an event function to detect the change of conditions and have it stop integration at the boundary, and then a new ode45 call to continue. When the boundaries are known times the setup is slightly easier.

[Matthew Hunt]

In the end, I just used a for loop to loop over the vector. The function itself is continuous.

Matlab needs to reduce the use of syms.

[Walter Roberson]

For the purpose of ode45, it is not enough that the function is continuous: the first two derivatives must be continuous as well. This is true no matter how you calculate the values inside the function.

[Walter Roberson]

MATLAB makes available the Symbolic toolbox to make it more convenient to create your calculations in a way that is more readable. It is common, for example, for people to make mistakes when applying the Chain Rule, so the Symbolic Toolbox can take care of that detail, making it more likely that your code does what you want. The only times that you "required" to use the Symbolic Toolbox are:

• you need more precision or greater range than is supported by double precision (for example you need factorial of a number greater than 15); or
• you need one of the specialized functions such hypergeometric implemented by the Symbolic Toolbox for which no double precision version is available