ODE method issue in matlab

Dear all, I have a legacy simulation model written in Matlab using fixed step Euler integration. This model has 4 control inputs and 10 state variables. I am rewriting the same program in Matlab. I have switched to ODE23 method (ODE23 is a variable step size method). In order to validate the new implementation, i passed the 4 control inputs as a function of time to both models (One is fixed step, other is variable). I am not getting same results (I suspect that both the models do not receive same control inputs at a time t). I wish to know the following

1. Is there a way i can make ODE23 or any other ODE method fixed step in Matlab?

2. If not, then in order to validate my new implementation, how do i access time variable in ODE23 method and write a routine to pass on similar inputs at same integration time as that of earlier implementation. I had currently used time variable available in ODEFUN.

regards,

Tarun

Answer 1

[Mike Hosea]

If the derivative odefun(t,y) depends on the four "control" functions of t (or even of both t and y), then these functions would be incorporated into the definition of yprime = odefun(t,y). I don't see how the fact that the solver is fixed-step or variable-step should matter, except that fixed-step simulations are numerically suspect to begin with. Even small amounts of instability will be magnified when integrating for any length using the (Forward) Euler method with fixed step sizes.

ODE23 provides a dizzying array of options, including event functions, nonnegativity contraints, and output functions. You can even have it return a struct (use one output) and then pass that struct to ODEXTEND to extend the integration. Using that facility, you could conceivably implement a pseudo-fixed-size integration with a loop, but the only benefit I can think of to that is if the control functions change in a non-smooth fashion at discrete points. In that case it would be best not to try to integrate over these points as if they weren't there, rather to step on, and subsequently from, them precisely.

[IISc]

Dear Mike,

Thanks for your answer. Meanwhile, i got a fixed step solver ODE1.m from earlier message trail on the this question/answers forum. In order to overcome my doubt over variable step integration i decided to use it. I implemented two seperate simulations: One in which euler method is applied my me, the second in which i use ODE1.m function to integrate. Now even the step size is also same!! The results of the two match only for 150 secs and then the second one diverges. The dynamical system which is being simulated is of fast dynamics and open loop unstable. I printed state variables in both cases and found that both match for a very short time(1 secs only) and then later start differing. I am pretty sure that vector u is same in both cases and is being read properly.

Can there be issue with the way feval function works...(The ODE1.m uses feval function for returning the value to state variable for the current input. regards,

[Mike Hosea]

Every time you take a step with a numerical method there are TWO sources of error that enter the solution. One is round-off error, and the other is discretization error. So the "noise" introduced in each step is a mixture of these two sources. If the problem being integrated is stupid-simple, then there may be no discretization error in some cases, but round-off error is hard to avoid. Now, accepting that we always have some noise, whether it is the mixture or round-off error only, the dynamical system may have limited stability, in which case a small perturbation may be amplified from one time step to the next. To complicate the matter further, the numerical integration method itself also may be stable or unstable in a given step.

The attraction of variable-step solvers is that they 1) control the discretization error (at least locally) and 2) keep the numerical method stable. A fixed-step method cannot do this except in very limited circumstances where some of the properties of the physical system are incorporated into, and assumed by, the numerical method. I am thinking of symplectic methods here. The bottom line is, that unless one does all the exact same calculations, even to the point of using the exact same expressions to calculate everything, even the next time value t (e.g. t = t + h in a loop and t = t0 + n*h would give slightly different results), the smallest perturbations are both inevitable and likely to be magnified. It is something of a numerical miracle when two different integrations match for very long.

The issue is that fixed-step integration is even more vulnerable than variable-step. What you must do to get comfortable with a fixed-step integration is repeat it with a smaller step sizes and verify that convergence is occurring. That is to say, you need to verify that the step size h is small enough that the method has any accuracy going at all. Sometimes, when the solution is sufficiently smooth, you can even use Richardson Extrapolation to get an error estimate and improve the accuracy of your results.