Re: Root Finding Methods Gaurenteed to Find All Root Between (xmin, xmax)

• To: mathgroup at smc.vnet.net
• Subject: [mg112711] Re: Root Finding Methods Gaurenteed to Find All Root Between (xmin, xmax)
• From: "Ingolf Dahl" <ingolf.dahl at telia.com>
• Date: Tue, 28 Sep 2010 06:03:04 -0400 (EDT)

```It would be very good if the Semenov method can be transferred into a
program that it able to find all roots to multi-dimensional non-linear
equation systems. I have now also looked a little into this, but now have
three questions/statements:

1. The test T1 checks if the rectangle is free of roots. This is implemented
by comparing the center function values in the rectangle with upper limits
for the spatial derivatives times half the length of the rectangle. The
upper limits of the spatial derivatives are calculated by interval analysis.
Why is it done in that way? Why not apply interval analysis directly on the
functions as such? If this derivative method gives better estimates than an
direct interval analysis, then I think it would be a good idea to
incorporate such a formula in the interval analysis of functions.

2. The corresponding question is relevant for the T2 test. Here the task
seems to be to check if the Jacobian determinant is non-zero, and this test
is done by comparing the center value in the rectangle with upper limits for
the spatial derivatives of the Jacobian determinant times half the length of
the rectangle. The upper limits of the spatial derivatives are calculated by
interval analysis. Why is it done in that way? Why not apply interval
analysis directly on the Jacobian determinant as such? Then the program does
not need to evaluate the second-order derivatives of the functions in the
equations.

3. The reason to check if the Jacobian determinant is non-zero is that this
is said to check if there is at most a single root in the rectangle.
However, I looked in my textbook from my first year at the university
("Matematisk analys II" by Hylt=E9n-Cavallius and Sandgren, Lund 1965, in
emphasize that a non-zero Jacobian implies that the function values around a
point are unique in an environment, which is small enough, but we might have
several identical function values in a larger area. Then we have several
roots to the corresponding equation. I will repeat the counterexample here,
slightly modified, since I think that it is of interest to several of the
Suppose that we have -1/4 <== x<== 1/4 and 0 <== y <== 5*Pi/2 and the equations

Exp[x]*Sin[y] - 1/Sqrt[2] ==== 0
Exp[x]*Cos[y] - 1/Sqrt[2] ==== 0

Then the Jacobian is Exp[2x], which is larger than zero. The rectangle also
passes Test2, in spite of the two roots {x -> 0, y -> Pi/4} and {x -> 0, y
-> 9 Pi/4}. Maybe this is the error Daniel Lichtblau has noticed, but
anyhow, the paper by Semenov seems now a bit less impressive. And a plus
score for the mathematicians of the previous century. We thus need a better
test, and a valid mathematical proof of it, to be able to search for all
roots in a guaranteed way. My feeling is that one might test the elements of
the Jacobian matrix - if no one of them changes sign within the rectangle
(or hyper rectangle), and if the Jacobian determinant is non-zero of one
sign, then I guess that multiple roots are impossible. But that I have not
proved, and there might be better tests. Any suggestions? And what is meant
by "generalized jacobian" in this context?

Best regards

Ingolf Dahl
ingolf.dahl at telia.com

> -----Original Message-----
> From: Andrzej Kozlowski [mailto:akoz at mimuw.edu.pl]
> Sent: den 20 september 2010 11:42
> To: mathgroup at smc.vnet.net
> Subject: [mg112561] Re: Root Finding Methods Gaurenteed to Find All Root
Between
> (xmin, xmax)
>
> The method that I referred to does have a practical limitation in that
> it relies on an implementation of interval arithmetic or some equivalent
> method of finding upper bounds of the function(s) involved in our
> equations and their derivatives (or generalised jacobians in the
> multi-dimensional case). In other words, method will work as long as we
> can compute  f[Interval{{a,b}]], f'[Interval[{a,b}]] etc.
>
> The algorithm works as follows. It subdivides the interval (or
> rectangle) in which we are looking for roots into a finite number of
> sub-rectangles. It them performs two tests on them, let's call them T1
> and T2. If a rectangle passes the testT1  we no for sure that there are
> no roots in it and we eliminate it from further search. If a rectangle
> fails test T1 we apply to it test T2. If the rectangle passes T2 we know
> with certainty that there is at most a single root in the rectangle. We
> then apply Newton's method to the rectangle - it will either find the
> required root or will not converge. This will always determine whether
> there is a root in the rectangle or not. If a rectangle fails both tests
> it is subdivided further into smaller and smaller rectangles. As there
> are only finitely many roots in the original rectangle they will all be
> separated and we will now for sure that we have found all the roots.
> Since the process goes on until all roots have been found the example
> described below does not contradict the method.
>
> A problem occurs when there is a multiple root. A rectangle containing
> such a root will always fail both tests and the algorithm will never
> complete. We can of course always decide to stop search when the
> rectangle containing possibly multiple roots is sufficiently small. We
> will then know that there is possibly a multiple root in the remaining
> rectangle. If we wish we can of course let the method run longer. If the
> roots are really distinct the process will stop for sure, but if they
> are multiple it can go on for ever.
>
> Here may be the right place to add that Daniel Lichtblau has noticed
> that the method given by Semenov in the case systems with 2 or more
> equations and variables contains a mistake (it is fine for the one
> dimensional case). There is a mistake in the application of the mean
> value theorem for vector valued functions in Semenev's published
> article, which causes the algorithm to miss certain roots. However, we
> have been able to find the correct version of Semenev's method based on
> the concept of generalized jacobian. I am now writing up our corrected
> version of Semenev's method and will soon implement it in Mathematica.
> Eventually I will replace my demonstrations on the Demosntrations site
> with new ones using the new version, which we believe to be correct.
>
> Andrzej Kozlowski
>
>
>
> On 19 Sep 2010, at 16:50, Ted Ersek wrote:
>
> > In [1] I gave a long response to the thread [FindRoots]
> >
> > It seems in [2] Andrzej Kozlowski was saying there are ways to
> > obtain the complete set of roots (provably complete) over a
> > closed and bounded set provided the function is C2 continuous
> > and does not have multiple roots in the domain of definition.
> >
> > I suppose he is talking about numeric methods to find all the roots.
> > If that is the case I would expect these methods to have practical
> > limitations in what they can do?  I mean a function could be very
> > ill-behaved and still meet the conditions above. Consider f[x] below.
> >
> > pnts==Table[ {x,Log[x]}, {x,3.0,250.0,0.01}];
> > Part[ pnts, 16000, 2] == -0.05;
> > f==Interpolation[ pnts, Method:>Spline ];
> >
> > From the output of  Plot[f[x], {x,3,250}] it looks like f[x] is
monotonic,
> > well behaved and having no roots between (3, 250). In fact it is
> > C2 continuous, but it does have two root close to 18.99.  It actually is
> > well behaved except for the small interval between (18.992, 19.006).
> >
> > Suppose we wanted to find all the roots between (3,250) of a function
> > defined as  SomeNumericalAlgorithm[x] and the function always evaluated
to
> > the same value as the InterpolatingFunction above. How can a numeric
> > RootFinding
> > algorithm be gaurenteed to find these roots near 18.99 if it doesn't
have
> > the
> > good furtune of taking one or more sample between (18.992, 19.006), or
> > sampling higher order derivatives in a interval a bit larger than
(18.992,
> > 19.006)?
> >
> > (******* References ********)
> >
> > [1]   http://forums.wolfram.com/mathgroup/archive/2010/Sep/msg00255.html
> >
> > [2]  ==
> http://forums.wolfram.com/mathgroup/archive/2010/Sep/msg00262.html
> >
> >

```

• Prev by Date: Solve Question - 2 Non zero values
• Next by Date: Re: Interpolate in polar coordinates or cartesian
• Previous by thread: Re: Root Finding Methods Gaurenteed to Find All Root Between (xmin, xmax)
• Next by thread: Re: Root Finding Methods Gaurenteed to Find All Root Between (xmin, xmax)