decimal sqrt method doesn't use round-half-even

[mdickinson]

BPO	1725899
Nosy	@tim-one, @facundobatista, @mdickinson

^{Note: these values reflect the state of the issue at the time it was migrated and might not reflect the current state.}

Show more details

GitHub fields:

assignee = 'https://github.com/facundobatista'
closed_at = <Date 2007-08-02.03:15:41.000>
created_at = <Date 2007-05-25.22:52:08.000>
labels = ['library']
title = "decimal sqrt method doesn't use round-half-even"
updated_at = <Date 2007-08-02.03:15:41.000>
user = 'https://github.com/mdickinson'

bugs.python.org fields:

activity = <Date 2007-08-02.03:15:41.000>
actor = 'facundobatista'
assignee = 'facundobatista'
closed = True
closed_date = None
closer = None
components = ['Library (Lib)']
creation = <Date 2007-05-25.22:52:08.000>
creator = 'mark.dickinson'
dependencies = []
files = []
hgrepos = []
issue_num = 1725899
keywords = []
message_count = 8.0
messages = ['32126', '32127', '32128', '32129', '32130', '32131', '32132', '32133']
nosy_count = 3.0
nosy_names = ['tim.peters', 'facundobatista', 'mark.dickinson']
pr_nums = []
priority = 'normal'
resolution = 'fixed'
stage = None
status = 'closed'
superseder = None
type = None
url = 'https://bugs.python.org/issue1725899'
versions = ['Python 2.5']

[mdickinson]

According to version 1.66 of Cowlishaw's `General Decimal Arithmetic
Specification' the square-root operation in the decimal module should
round using the round-half-even algorithm (regardless of the rounding
setting in the current context). It doesn't appear to do so:

>>> from decimal import *
>>> getcontext().prec = 6
>>> Decimal(9123455**2).sqrt()
Decimal("9.12345E+6")

The exact value of this square root is exactly halfway between two
representable Decimals, so using round-half-even with 6 digits I'd
expect the answer to be rounded to the neighboring representable
Decimal with *even* last digit---that is,

Decimal("9.12346E+6")

This bug only seems to occur when the number of significant digits in
the argument exceeds the current precision (indeed, if the number of
sig. digits in the argument is less than or equal to the current
precision then it's impossible for the square root to be halfway
between two representable floats).

It seems to me that this is a minor bug that will occur rarely and is
unlikely to have any serious effect even when it does occur; however,
it does seem to be a deviation from the specification.

[tim-one]

Doesn't the spec also require that /inputs/ get rounded to working precision /before/ operations are performed? If so, the exact square 83237431137025 is rounded down to 83237400000000, and sqrt(83237400000000) is 9.12345E+6 rounded to 6 digits.

[mdickinson]

I don't think inputs should be rounded; note 1 near the start of the `Arithmetic Operations' section of the specification says:

"Operands may have more than precision digits and are not rounded before use."

[mdickinson]

Some more strange results for sqrt(): with Emax = 9, Emin = -9 and prec = 3:

1. In the following, should the Subnormal and Underflow flags be set? The result before rounding is subnormal, even though the post-rounding result is not, and the spec seems to say that those flags should be set in this situation. But Cowlishaw's reference implementation (version 3.50) doesn't set these flags. (If 9.998 is replaced by 9.99 then the flags *are* set, which seems inconsistent.)

>>> Decimal("9.998E-19").sqrt()
Decimal("1.00E-9")
>>> getcontext()
Context(prec=3, rounding=ROUND_HALF_EVEN, Emin=-9, Emax=9, capitals=1, flags=[Rounded, Inexact], traps=[DivisionByZero, Overflow, InvalidOperation])

2. As I understand the spec, the following result is incorrect:

>>> Decimal("1.12E-19").sqrt()
Decimal("3.4E-10")

(The true value of the square root is 3.34664...e-10, which should surely be rounded to 3.3E-10, not 3.4E-10?).
But again, Cowlishaw's reference implementation also gives 3.4E-10.

3. Similarly for the following

>>> Decimal("4.21E-20").sqrt()
Decimal("2.0E-10")

The answer should, I think, be 2.1E-10; here the reference implementation also gives 2.1E-10.

4. And I can't justify this one at all...

>>> Decimal("1.01001").sqrt()
Decimal("1.01")

Shouldn't this be 1.00? But again Python agrees with the reference implementation.

Either all this is pretty mixed up, or I'm fundamentally misunderstanding things. I have a patch that I think fixes all these bugs; if there's any agreement that they really are bugs then I'll post it. Can anyone shed any light on the above results?

[mdickinson]

One more result. This is definitely getting suspicious now--I'll submit my patch.

>>> from decimal import *
>>> getcontext().prec = 3
>>> Decimal(11772).sqrt()
Decimal("109")
>>> Decimal(11774).sqrt()
Decimal("108")

[tim-one]

Of course you're right, the spec does say inputs shouldn't be rounded. And that last example is horrendous: sqrt should definitely be monotonic (a floating function f is "monotonic" if it guarantees f(x) >= f(y) whenever x >= y; you found x and y such that x > y but sqrt(x) < sqrt(y) -- ouch!).

[mdickinson]

See patch 1741308. I'll contact Mike Cowlishaw to verify that the reference implementation really is buggy.

[facundobatista]

Fixed in the revisions 56654 to 56656, in the decimal-branch, using patches generated by Mark, thanks!