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 that 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.
