Incorrect rounding at low precision

asked 2022-05-12 14:28:25 +0100

I'm trying to use Sage to implement a slightly modified bisection algorithm on a certain function, and at each step, displaying the interval endpoints as well as the function at these endpoints. I don't need much precision for the latter numbers; I just want an indication of the size and sign of these numbers. However, the numerical_approx() function is giving me weird responses when I specify the number of digits. For example, the following calculation (which is the function value at the greater endpoint, in the final step of the algorithm)

numerical_approx(1/2427673471984599040*(1530627*sqrt(4630419677705) - 4630419677705)*(sqrt(4630419677705) - 482051) + 482051/524288)

evaluates to

1.64597827478907e-7

which appears to be "reasonable", though I don't know for certain. However, when I write,

numerical_approx(1/2427673471984599040*(1530627*sqrt(4630419677705) - 4630419677705)*(sqrt(4630419677705) - 482051) + 482051/524288,digits=3)

it evaluates to

-0.0000160

which is clearly unreasonable: it's the wrong sign and it's too large. Further more, these endpoints come back to this same number a lot, despite how the input changes. At one point, in the middle of this string of -0.0000160s, Sage estimates 0.000, just once, and goes back to -0.0000160.

Changing the precision to 4 and 5 does seem to improve the situation (correct sign, more reasonable magnitude), but the significant digits are all different!

So, my questions are:

What is going wrong here?
Can I fix it?
What methods could I use to be certain that these approximations have 3 correct significant figures?

answered 2022-05-12 16:23:46 +0100

tmonteil
27413 ●31 ●204 ●517 http://wiki.sagemath.o...

updated 2022-05-12 16:27:32 +0100

To get a reliable interval, you can use the real interval field RIF:

sage: a = 1/2427673471984599040*(1530627*sqrt(4630419677705) - 4630419677705)*(sqrt(4630419677705) - 482051) + 482051/524288
sage: RIF(a).endpoints()
(1.64597823149037e-7, 1.64597830476510e-7)

So, your first result seems correct.

Now, regarding your low precision result, let us observe the following:

sage: b = a.numerical_approx(digits=3)
sage: b
-0.0000610
sage: R = b.parent()
sage: R
Real Field with 14 bits of precision

The computation of the numerical value of the symbolic expression is done within the "field" R of floating-point numbers with only 14 bits of precision, meaning that every number apperaing in the symbolic expression is converted in R and then the operations are done within R, in particular, numbers with large precision get rounded:

sage: R(2427673471984599040)
2.43e18

Since your symbolic expression has differences, and several coefficients get rounded, there is no reason why the sign must be preserved.

If you do the computation by hand, you can check:

sage: 1/R(2427673471984599040)*(R(1530627)*sqrt(R(4630419677705)) - R(4630419677705))*(sqrt(R(4630419677705)) - R(482051)) + R(482051)/R(524288)
-0.0000610

Now, if you only want to keep the 3 digits of precision, you can first do the computation within a field of large precision, and round afterwards:

sage: a.numerical_approx().numerical_approx(digits=3)
1.65e-7

sage: R(a.numerical_approx())
1.65e-7

sage: R(RR(a))
1.65e-7

edit flag offensive delete link

Comments

Thank you for the answer! The RIF command strikes me as particularly handy. There is something unsettling about this though: how am I supposed to detect when the precision of the underlying field is lacking precision and giving me an erroneous answer? Suppose I'm not so astute next time to notice that the result is wrong, or I'm not specifically casting an eye over the figures. Is there a way to ask sage to grant me precise results? The RIF function suggests to me that sage was perfectly capable of giving the number to 3 correct significant figures, but the numerical_approximation function did not do this. Is there a built-in function (possibly with a larger footprint) that guarantees the figures it produces (or warns me when it cannot guarantee the figures)?

Theo Bendit ( 2022-05-13 11:39:19 +0100 )edit

Actually, I think I just answered my own question.

Theo Bendit ( 2022-05-13 11:42:45 +0100 )edit

Indeed, you can replace RIF with RealIntervalField(prec) with a higher precision.

tmonteil ( 2022-05-13 18:40:17 +0100 )edit

add a comment

Your Answer

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

Incorrect rounding at low precision

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

Incorrect rounding at low precision edit

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

Incorrect rounding at low precision