Re: fzero initial bracketing [WAS: bug #31070]

[Jaroslav Hajek]

On Mon, Sep 20, 2010 at 5:24 PM, Carlo de Falco <address@hidden> wrote:
> Dear all,
>
> This discussion started as bug #31070 in the tracker,
> but Jaroslav convinced me that this is not really a bug
> but rather a suggestion for improvement so I am moving
> the discussion here.
>
> The algorithm used by fzero [1] for finding x s.t. f(x) = 0,
> requires that the user provide an interval [a, b] such that f(a) * f(b) < 0
>
> for compatibility with Matlab, though, fzero allows the user to just
> provide an initial guess x0 for x and attempts at internally
> constructing the initial
> bracketing by looking for a and b in some 'neighborhood' of x0
>
> To my knowledge there is no algorithm which is guaranteed to find such
> initial bracketing,

Such an algorithm may in fact be constructed given some assumptions, e.g.

x0 = initial guess;
n = 1;
fx0 = fun (x0);

while (1)
  for k = 1:2:n
    d = k*2^(k-n);
    x1 = x0 + d;
    if (sign (fun (x1)) * sign (fx0) <= 0)
      break;
    endif
    x1 = x0 - d;
    if (sign (fun (x1)) * sign (fx0) <= 0)
      break;
    endif
  endfor
  n++;
endwhile

Lemma:
In infinite precision arithmetics, if fun is a continuous function
such that fun(x0)*fun(x1) < 0 for at least one x1, the above algorithm
will always find a valid bracketing.

Proof:
Left to the reader :)

The question here is whether something like this is practical to use.

> so it is possible that the bracketing is not found in which case fzero
> will throw an error and die.
>
> Nonetheless, the use of fzero with an initial guess rather than an
> initial bracketing as input is quite
> common and there is a lot of Matlab code out there that relies on it,
> even though other functions
> would be probably more appropriate for this purpose.
>
> Indeed, Matlab reportedly seems to be able to construct the initial
> bracketing more often than Octave,

Reportedly? I'd like to see the reports first.


-- 
RNDr. Jaroslav Hajek, PhD
computing expert & GNU Octave developer
Aeronautical Research and Test Institute (VZLU)
Prague, Czech Republic
url: www.highegg.matfyz.cz