[<prev] [next>] [<thread-prev] [day] [month] [year] [list]
```Date: Sat, 11 Nov 2023 11:26:58 +0100
From: Szabolcs Nagy <nsz@...t70.net>
To: Damian McGuckin <damianm@....com.au>
Cc: musl@...ts.openwall.com
Subject: Re: Double Rounding in ....

* Damian McGuckin <damianm@....com.au> [2023-11-11 19:11:41 +1100]:
>
> 	hypot, hypotf, hypotl
> and
> 	scalbn, scalbnf, scalbnl
>
> Does anybody have a good/clear explanation for how double rounding occurs in
> the case of squaring-of/multiplication-by small numbers and how the action
> taken in the code addresses that? Not sure if it needs to mention

n * 2^n is normally exact, but in the subnormal range there
are less precision bits so rounding can occur.

if you do the scaling in two steps (e.g. because the scale
factor has such a big exponent that it is not representible
in a single float) then you can get double rounding.

e.g.

scalbn(0x1.000000000000bp-1, -1024) == 0x1.0000000000008p-1025

normally the scaling is done in two steps: x*0x1p-1022*0x1p-2
because up to 0x1p-1022 the scaling is easy to construct from an int,
while 0x1p-1024 is subnormal.

with naive implementation the bit pattern at the end is
..001011 x
..00110  x*0x1p-1022 (rounded)
..010    x*0x1p-2 (rounded)

with the fix
..001011 x*0x1p-969 (969 = 1022-53)
..001    x*0x1p-55 (rounded)

so simply split the scale factor such that the first scaling can
only do rounding if the second scaling (<0x1p-53) completely rounds
the entire result away.

>
> 	FLT_EVAL_METHOD
>
> because certainly the float and double versions of those routines are
> written in terms of float_t and double_t respectively.
>
> Thanks - Damian
>
> Pacific Engineering Systems International ..... 20D Grose St, Glebe NSW 2037
> Ph:+61-2-8571-0847 .. Fx:+61-2-9692-9623 | unsolicited email not wanted here
> Views & opinions here are mine and not those of any past or present employer
```

Confused about mailing lists and their use? Read about mailing lists on Wikipedia and check out these guidelines on proper formatting of your messages.