# \[Ann\] IntervalRootFinding.jl for finding all roots of a multivariate function

**URL:** <https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515>\
**Category:** Community\
**Tags:** intervals, roots, scientific-computing\
**Created:** [March 5, 2018, 6:38pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515 "2018-03-05T18:38:42Z")\
**Posts on this page:** 13\
**Page:** 3

<div class="post-metadata">

**Author:** ![foobar\_lv2](https://avatars.discourse-cdn.com/v4/letter/f/ee59a6/32.png) [@foobar\_lv2](https://discourse.julialang.org/u/foobar_lv2)\
**Post date:** [March 13, 2018, 1:41pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/41 "2018-03-13T13:41:22Z")

</div>

> [@dpsanders](#):
>
> What we have shown is that there are no simultaneous roots of g(I\_k) and h(I\_k).

Oops. Yes, of course.

> Is there a reference I can look at for this stuff?

Hmm. I would myself use google and wikipedia mainly. “Computational Homology” by Kaczynski, Mischaikow, Mrozek (available on ligben) contains some context, and I know that Mischaikow’s group is also in the business of producing code, not just theorems. That being said, I failed to find any references at all that explicitly talks about the connection between Lefschetz index and interval-math / rectangle based root finding. So this might turn out to become a research project yet 😉

* * *

The way I thought about what you are doing:

You decompose the domain into rectangles which intersect at their boundaries (so you leave no “gap” between `x` and `nextfloat(x)`). For each rectangle, and each coordinate function, write down a `+1` if the coordinate is positive in the rectangle, a `-1` if it is negative, and a `0` if it might be zero.

The rectangles that might contain roots are labeled with `zeros(dim)`.

The intersection of R,R’ gets label zero if both R,R’ are labeled zero, plus if one is labeled plus, minus if one is labeled minus (for each coordinate separately). The plus/minus case cannot happen, due to the fact that the rectangles intersect at their boundaries.

Now, every rectangle R’s boundary decomposes into a union of intersections between R and some R’. The same holds for each subset M of rectangles; the boundary of the union decomposes into a union of intersections between one M-rectangle and one M-complement-rectangle.

Hence, you obtain (without extra calculation) a map from the boundary of the union of maybe-root-containing rectangles to the space of non-zero labels. Even better, this map is continuous, in the sense that the labels are non-contradictory at the transitions, where you move from one bounding rectangle to the next; non-contradictory means no jump from plus to minus (must go via zero).

This map contains the Lefshetz-index! And, by finding connected components of your set of rectangles, you can compute the Lefshetz indices separately (so: each connected maybe-root-containing region gets its own personal Lefshetz index, and each connected component of the boundary of such a region also gets its own index, e.g. if you are left with a maybe-zero-containing double-annulus, like a disc with two disjoint smaller discs punched out).

* * *

Is this picture of your procedure roughly correct? If this really is the information we can go on, then I’ll think a little about how to get a proper homology solver for this (with the goal of almost no additional f-evaluations, up to the difference between `x` and `nextfloat(x)`).

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [March 15, 2018, 9:19pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/42 "2018-03-15T21:19:55Z")

</div>

Also arbitrary functions of complex numbers, e.g. `sin`, seem to be problematic currently.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [April 30, 2018, 7:53pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/43 "2018-04-30T19:53:18Z")

</div>

This article is probably relevant:

> **[An efficient degree-computation method for a generalized method of bisection...](https://link.springer.com/article/10.1007/BF01404868)**
>
> LetP be ann-dimensional polyhedron and let $$b(P) = \\sum\\limits\_{q = 1}^m {\\langle X\_1^q , \\ldots ,X\_n^q \\rangle } $$ be the oriented boundary ofP in terms of the oriented (n−1)-simplexesS q =〈X 1 q ,...,X n q 〉,q=1,...,m. LetF=(f 1,...,f n):P→R n,...

---

<div class="post-metadata">

**Author:** ![abulak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/abulak/32/28314_2.png) [@abulak](https://discourse.julialang.org/u/abulak)\
**Post date:** [July 21, 2018, 9:30pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/44 "2018-07-21T21:30:20Z")

</div>

A different (already implemented) reference:

[https://arxiv.org/abs/1207.6331](https://arxiv.org/abs/1207.6331)

---

<div class="post-metadata">

**Author:** ![iwelch](https://avatars.discourse-cdn.com/v4/letter/i/8c91f0/32.png) [@iwelch](https://discourse.julialang.org/u/iwelch)\
**Post date:** [September 3, 2018, 6:55pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/45 "2018-09-03T18:55:05Z")

</div>

hi david—is Interval\*.jl working with julia 1.0? planned to become working again?

regards,

/iaw

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 3, 2018, 7:14pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/46 "2018-09-03T19:14:36Z")

</div>

There are open pull requests updating IntervalArithmetic.jl and IntervalRootFinding.jl to full compatibility with Julia 1.0. I hope to tag the former today.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 9, 2018, 2:53pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/47 "2018-09-09T14:53:22Z")

</div>

Versions 0.15 of IntervalArithmetic.jl and version 0.4 of IntervalRootFinding.jl have been released. These are fully compatible with both Julia 0.7 and Julia 1.0 (and, in fact, the latest developement version of Julia, and presumably all 1.x versions).

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 25, 2022, 2:26pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/49 "2022-09-25T14:26:35Z")

</div>

What code are you running?

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 25, 2022, 2:37pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/50 "2022-09-25T14:37:56Z")

</div>

Actually I don’t think that function has isolated roots, which is why you’re getting singular Jacobians. The zero set should be a three-dimensional surface in 4d space.

---

<div class="post-metadata">

**Author:** ![vinod\_kumar](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/vinod_kumar/32/42952_2.png) [@vinod\_kumar](https://discourse.julialang.org/u/vinod_kumar)\
**Post date:** [September 25, 2022, 9:07pm UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/53 "2022-09-25T21:07:32Z")

</div>

if I make b=0, do I able to find roots then?

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 26, 2022, 3:17am UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/54 "2022-09-26T03:17:48Z")

</div>

You are not calculating roots of the function `f`, but rather critical points, i.e. roots of the gradient.

Nonetheless, I think that what I said holds, and there are whole curves (?) of critical points satisfying e.g.

```julia
a - b == pi
b - c == pi
d == anything

```

You can try with `Bisection` instead of `Newton`, since that does not use a gradient, but even with a tolerance of `1e-1` it takes a very long time to compute, for this reason: it will return any box that contains any part of the zero set:

```julia
A = roots(∇f, (0..6.28) × (0..6.28) × (0..6.28) × (0..6.28), Bisection, 1e-1)

```

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 26, 2022, 3:18am UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/55 "2022-09-26T03:18:27Z")

</div>

Fixing `b = 0` may help, yes.

By the way, in the future, please do not reply to an old post like this; rather make a new one.

---

<div class="post-metadata">

**Author:** ![dpsanders](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dpsanders/32/3573_2.png) [@dpsanders](https://discourse.julialang.org/u/dpsanders)\
**Post date:** [September 26, 2022, 5:12am UTC](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515/56 "2022-09-26T05:12:37Z")

</div>

Try computing the derivative symbolically, e.g. using `Symbolics.jl`, and then finding the constraints that follow from setting gradient = 0. I think there will be too few independent constraints for the number of variables.

[Previous page](https://discourse.julialang.org/t/ann-intervalrootfinding-jl-for-finding-all-roots-of-a-multivariate-function/9515.md?page=2)