Ask Your Question
1

Finding root with difficulties

asked 2013-08-08 20:46:22 -0500

MvG gravatar image

I'm working with this example:

sage: var('x a t', domain='real')
sage: assume(a > 0, a < 1)
sage: xt(t) = solve(integrate(-4*sqrt(17)/(9*x^2 + 6*x - 16),
...                 (x, 0, a)).subs(a = x) == t, x)[2].rhs()
sage: find_root(xt - 1, 2, 3)
0.0
sage: plot(xt - 1, (t, 2, 3))
verbose 0 (2424: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 137 points.
verbose 0 (2424: plot.py, generate_plot_points) Last error message: ''

The plot reports errors, but it looks reasonable enough. The existance of a root appears very clear, but sage is unable to find it. Furthermore, sage does not even properly report its inability to find a root, but instead returns a value as if there were no problems. The fact that this value lies outside the interval is an indication that there might be something fishy happening here.

Is there any way to get at my root numerically, or perhaps even symbolically? I know that the function xt is an ugly beast, but some form of robust bisection should be possible, I think.

edit retag flag offensive close merge delete

1 answer

Sort by ยป oldest newest most voted
2

answered 2013-08-08 22:16:08 -0500

Eviatar Bach gravatar image

updated 2013-08-09 12:15:26 -0500

Unfortunately, the current root-finding only works with real numbers (see #14991), and find_root is encountering complex numbers.

With mpmath's root-finding you get the following answer, although the syntax is extremely clumsy at the moment:

from sage.libs.mpmath.utils import call
import mpmath
def mpmath_xt(t, prec):
    t = ComplexField(prec)(t.real, t.imag)
    ft = xt(t).n(prec)
    return mpmath.mpc(ft.real(), ft.imag()) - 1

Then,

sage: call(mpmath.findroot, lambda t: mpmath_xt(t, 100), [2, 3], parent=ComplexField(100))
2.4629947809517837884000550617 + 2.0631855558523985234600213842e-88*I

Change 100 in the above to whatever precision you want.

Hope this helps!

edit flag offensive delete link more

Comments

Hmmm. *Why* is my way of passing the function not appropriate? As far as I understand things, the lambda is a python construct, which sage will treat as a black box. The expression, on the other hand, is something sage can work with, e.g. by turning it into some optimized fastcall thingy. So why would a black box be more appropriate than a symbolic expression?

MvG gravatar imageMvG ( 2013-08-09 08:00:32 -0500 )edit

Okay, I looked into it, and I think you're right; it's not the wrong way to do it. The problem is precisely that it is converting the expression to a `fast_float`, but then chokes when it encounters complex numbers. Look what happens: sage: (xt - 1)._fast_float_(t)(1) nan

Eviatar Bach gravatar imageEviatar Bach ( 2013-08-09 12:17:05 -0500 )edit

I should also mention that you need to be careful with the precision; your expression is so big that there's probably roundoff error, and the answer is not correct to the precision you specify. So you should work at a higher precision than the one you want and then discard some of the trailing digits.

Eviatar Bach gravatar imageEviatar Bach ( 2013-08-09 12:21:48 -0500 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

Stats

Asked: 2013-08-08 20:46:22 -0500

Seen: 297 times

Last updated: Aug 09 '13