Re: 'InitialStep' and 'MaxStep' in ode45

[Thomas Treichl]

John B. Thoo schrieb:
> Hi, everyone.  I have a question about 'InitialStep' and 'MaxStep' in  
> ode45.
>
> I tried to solve an equation using
>
>      [T, u] = ode45 (@F, [t0 tf], [u0, u0_t]);
>
> and this is what I got:
>
> octave-3.2.0:4> run_twoway
> warning: Option "RelTol" not set, new value 0.000001 is used
> warning: Option "AbsTol" not set, new value 0.000001 is used
> warning: Option "InitialStep" not set, new value 1.200000 is used
> warning: Option "MaxStep" not set, new value 1.200000 is used
> error: Solving has not been successful. The iterative integration  
> loop exited at time t = 5.818179 before endpoint at tend = 12.000000  
> was reached. This may happen if the stepsize grows smaller than  
> defined in vminstepsize. Try to reduce the value of "InitialStep" and/ 
> or "MaxStep" with the command "odeset".
> octave-3.2.0:4>
>
> So I halved tend and set 'InitialStep' and 'MaxStep' as follows:
>
>      options = odeset ('InitialStep', 1.0e-3, 'MaxStep', 1.0e-3, ...
>                        'RelTol', 1.0e-6, 'AbsTol', 1.0e-6);
>      [T, u] = ode45 (@F, [t0 tf], [u0, u0_t], options);
>
> This time I got:
>
> octave-3.2.0:4> run_twoway
> error: Solving has not been successful. The iterative integration  
> loop exited at time t = 5.553544 before endpoint at tend = 6.000000  
> was reached. This may happen if the stepsize grows smaller than  
> defined in vminstepsize. Try to reduce the value of "InitialStep" and/ 
> or "MaxStep" with the command "odeset".
> octave-3.2.0:4>
>
>
> My question is this:
>
>    If the first time 'InitialStep' and 'MaxStep' were set at 1.2  
> resulted in ode45 exiting prematurely at time t = 5.818179, then why  
> does decreasing 'InitialStep' and 'MaxStep' result in the solver  
> exiting prematurely at a *shorter* time t = 5.553544?

Hi John,

this strongly depends on your ODE system. I think you have a singularity 
somewhere at t=5.5..6.0, ie. that this is the third possibility for your problem 
that ode23..ode78 can not solve for you however you choose InitialStep and/or 
MaxStep. If your problem is that kind of stiff then you maybe should use a 
solver that can solve stiff problems, cf. description about the different 
solvers in odepkg.pdf and the description about all the different options that 
can be set.

InitialStep and MaxStep in ode23..ode78 are not computed (as like the solvers in 
Matlab do) but must always be given by the user. If no value is given (eg. your 
first trial) then simply (tf-t0)/10 is taken - this might be completely wrong 
and depends on which ODE system should be solved.

Now, if you choose another InitalStep then solving can completely different (in 
your first trial we have t1_1=(tf-t0)/10, in your second trial t1_2=1e-3, you 
see that t1_1 and t1_2 are already different in your first InititialStep, it 
would be pure luck if you would hit tx_1 = 5.818179 and tx_2 = 5.553544).

Just for info, you normally should hit tx_1=tx_2 exactly if you take fixed step 
sizes, eg.

  [T, u] = ode45 (@F, [t0:tf], [u0, u0_t]);
                   ## --> <-- pass a vector with time values

> Another question:
>
>    How should I choose to set 'InitialStep' and 'MaxStep'?  And  
> should they always be set to the same value?

Ok, I hope I could explain InitialStep, MaxStep is different: There are errors 
computed every time you solve one delta time (AbsTol, RelTol). The solver tries 
to increase the dt automatically every time one time step was solved 
successfully (ie. computed error < AbsTol and/or RelTol). Let's say we had 
dt_old = 1 and the solver thinks dt_new = 2. If you set MaxStep = 1.5 then the 
user decreases dt to 1.5 and takes control over the solver's estimation.

*If* ode45 is the right choice for your problem (because of stiffness) then not 
necessarily InitialStep and MaxStep have the same values. It depends once again 
what kind of system you have - might there be stability problems in your ODE 
system right at the beginning of solving (make InitialStep smaller) or somewhere 
 in the middle or the end (make MaxStep smaller).

Best regards

  Thomas