Re: Review ode csets

[c.]

On 18 Oct 2016, at 20:40, Rik <address@hidden> wrote:

> 10/18/16
> 
> Carlo,
> 
> Can you review this cset that I just committed? 
> http://hg.savannah.gnu.org/hgweb/octave/rev/7efa2d0e22c9.
> 
> I changed OutputFcn to behave in a more Matlab-compatible way.  However,
> part of the change was subtle.  When using the Refine option to generate
> intermediate timesteps, OutputFcn was called for each of the timesteps. 
> However, only the final call in the last timestep had the ability to halt
> the solver by returning a value of TRUE.  This seemed wrong, but maybe
> there was a reason for it.  The code is:
> 
> +            stop_solve = false;
>             for ii = 1:numel (approxtime)
> -              pltret = feval (options.OutputFcn, approxtime(ii),
> -                              approxvals(:, ii), [],
> -                              options.funarguments{:});
> +              stop_solve = feval (options.OutputFcn,
> +                                  approxtime(ii), approxvals(:, ii), [],
> +                                  options.funarguments{:});
> +              if (stop_solve)
> +                break;  # break from inner loop
> +              endif
>             endfor
> -            if (pltret)  # Leave main loop
> +            if (stop_solve)  # Leave main loop
> 
> There is also a patch I submitted to https://savannah.gnu.org/bugs/?49364
> that needs review.
> 
> Thanks,
> Rik


Rik,

I see the following problem with this changeset:

Let myofun be the following:

function val = myofun (t, y, s)
  disp ('t = ')
  disp (t)
  val = any ((t > 3.2e-2) & (t < .1));
  disp ('val =')
  disp (val)
end

and consider the following example:

>> myfun  = @(t, y) 1;
>> [t, y] = ode45 (myfun, [0 1], 1, odeset ('maxstep', 1, 'initialstep', .1, 
>> 'refine', 3, 'outputfcn', @myofun)); t

In Matlab this prints: 
------------------------------------------
t = 
     0     1

val =
     0

t = 
    0.0333    0.0667    0.1000

val =
     1

t = 
val =
     0


t =

         0
    0.0333
    0.0667
    0.1000
------------------------------------------

So integration stops after the call to myofun with t = [1/30 2/30 3/30] but all 
the results are returned

In Octave this produces an error:

------------------------------------------
t = 
   0
   1
val =
0
t = 
0
val =
0
error: approxvals(_,2): but approxvals has size 4x1
error: called from
    integrate_adaptive at line 249 column 24
    ode45 at line 222 column 12
------------------------------------------

c.