ASKSAGE: Sage Q&A Forum - RSS feedhttps://ask.sagemath.org/questions/Q&A Forum for SageenCopyright Sage, 2010. Some rights reserved under creative commons license.Wed, 31 Oct 2012 13:04:34 +0100find roots with arbitrary precisionhttps://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/It seems that `find_root` will always convert its parameter range to `float`, even if `a` and `b` were originally given in some arbitrary precision real number field. Is there a variant of this algorithm that can find roots to arbitrary precision?
Even greater would be some algorithm that can make use of interval arithmetic based on `RealIntervalField`, which in particular will return an interval that is guaranteed to contain the zero. I have written some code along these lines:
def bisect_interval(f, i, d):
# find zero of function f in interval i with desired diameter d
# f must be monotonically increasing and must contain a zero in i
d2 = d/2
zero = i.parent()(0)
hi = i
while hi.absolute_diameter() > d2:
l, r = hi.bisection()
c = l.intersection(r)
fc = f(c)
if fc > zero:
hi = l
else:
hi = r
lo = i
while lo.absolute_diameter() > d2:
l, r = lo.bisection()
c = l.intersection(r)
fc = f(c)
if fc < zero:
lo = r
else:
lo = l
return lo.union(hi)
Am I reproducing functionality that's already available somewhere in Sage? If not, do you consider this functionality worth adding? Should it use some better algorithm than simple bisection? Is the algorithm even correct? You don't *have* to answer all of these questions, as my core question remains the one about an arbitrary precision version of `find_root`. But additional answers would be welcome.Wed, 31 Oct 2012 05:44:29 +0100https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/Comment by kcrisman for <p>It seems that <code>find_root</code> will always convert its parameter range to <code>float</code>, even if <code>a</code> and <code>b</code> were originally given in some arbitrary precision real number field. Is there a variant of this algorithm that can find roots to arbitrary precision?</p>
<p>Even greater would be some algorithm that can make use of interval arithmetic based on <code>RealIntervalField</code>, which in particular will return an interval that is guaranteed to contain the zero. I have written some code along these lines:</p>
<pre><code>def bisect_interval(f, i, d):
# find zero of function f in interval i with desired diameter d
# f must be monotonically increasing and must contain a zero in i
d2 = d/2
zero = i.parent()(0)
hi = i
while hi.absolute_diameter() > d2:
l, r = hi.bisection()
c = l.intersection(r)
fc = f(c)
if fc > zero:
hi = l
else:
hi = r
lo = i
while lo.absolute_diameter() > d2:
l, r = lo.bisection()
c = l.intersection(r)
fc = f(c)
if fc < zero:
lo = r
else:
lo = l
return lo.union(hi)
</code></pre>
<p>Am I reproducing functionality that's already available somewhere in Sage? If not, do you consider this functionality worth adding? Should it use some better algorithm than simple bisection? Is the algorithm even correct? You don't <em>have</em> to answer all of these questions, as my core question remains the one about an arbitrary precision version of <code>find_root</code>. But additional answers would be welcome.</p>
https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/?comment=18778#post-id-18778I think this is because `find_root` uses SciPy, which probably doesn't have arbitrary precision.Wed, 31 Oct 2012 08:12:45 +0100https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/?comment=18778#post-id-18778Answer by achrzesz for <p>It seems that <code>find_root</code> will always convert its parameter range to <code>float</code>, even if <code>a</code> and <code>b</code> were originally given in some arbitrary precision real number field. Is there a variant of this algorithm that can find roots to arbitrary precision?</p>
<p>Even greater would be some algorithm that can make use of interval arithmetic based on <code>RealIntervalField</code>, which in particular will return an interval that is guaranteed to contain the zero. I have written some code along these lines:</p>
<pre><code>def bisect_interval(f, i, d):
# find zero of function f in interval i with desired diameter d
# f must be monotonically increasing and must contain a zero in i
d2 = d/2
zero = i.parent()(0)
hi = i
while hi.absolute_diameter() > d2:
l, r = hi.bisection()
c = l.intersection(r)
fc = f(c)
if fc > zero:
hi = l
else:
hi = r
lo = i
while lo.absolute_diameter() > d2:
l, r = lo.bisection()
c = l.intersection(r)
fc = f(c)
if fc < zero:
lo = r
else:
lo = l
return lo.union(hi)
</code></pre>
<p>Am I reproducing functionality that's already available somewhere in Sage? If not, do you consider this functionality worth adding? Should it use some better algorithm than simple bisection? Is the algorithm even correct? You don't <em>have</em> to answer all of these questions, as my core question remains the one about an arbitrary precision version of <code>find_root</code>. But additional answers would be welcome.</p>
https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/?answer=14216#post-id-14216You can compare your function with that in mpmath:
sage: import mpmath
sage: mpmath.findroot?Wed, 31 Oct 2012 05:58:58 +0100https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/?answer=14216#post-id-14216Comment by MvG for <p>You can compare your function with that in mpmath:</p>
<pre><code>sage: import mpmath
sage: mpmath.findroot?
</code></pre>
https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/?comment=18777#post-id-18777The mpmath package and its findroot function look interesting. I see no mention of any support for interval arithmetic, so this code does something different from my own code even for `solver='bisect'`. And it does not integrate well with Sage functions, as even a simple expression like `2*arccos(x)` causes a `TypeError`. I just wrote [a follow-up question about this](http://ask.sagemath.org/question/1937/).Wed, 31 Oct 2012 13:04:34 +0100https://ask.sagemath.org/question/9486/find-roots-with-arbitrary-precision/?comment=18777#post-id-18777