[Numpy-discussion] Odd-looking long double on windows 32 bit

[Matthew Brett]

Hi,

On Sun, Nov 13, 2011 at 8:21 AM, Charles R Harris
<charlesr.harris at gmail.com> wrote:
>
>
> On Sun, Nov 13, 2011 at 12:57 AM, Matthew Brett <matthew.brett at gmail.com>
> wrote:
>>
>> Hi,
>>
>> On Sat, Nov 12, 2011 at 11:35 PM, Matthew Brett <matthew.brett at gmail.com>
>> wrote:
>> > Hi,
>> >
>> > Sorry for my continued confusion here.  This is numpy 1.6.1 on windows
>> > XP 32 bit.
>> >
>> > In [2]: np.finfo(np.float96).nmant
>> > Out[2]: 52
>> >
>> > In [3]: np.finfo(np.float96).nexp
>> > Out[3]: 15
>> >
>> > In [4]: np.finfo(np.float64).nmant
>> > Out[4]: 52
>> >
>> > In [5]: np.finfo(np.float64).nexp
>> > Out[5]: 11
>> >
>> > If there are 52 bits of precision, 2**53+1 should not be
>> > representable, and sure enough:
>> >
>> > In [6]: np.float96(2**53)+1
>> > Out[6]: 9007199254740992.0
>> >
>> > In [7]: np.float64(2**53)+1
>> > Out[7]: 9007199254740992.0
>> >
>> > If the nexp is right, the max should be higher for the float96 type:
>> >
>> > In [9]: np.finfo(np.float64).max
>> > Out[9]: 1.7976931348623157e+308
>> >
>> > In [10]: np.finfo(np.float96).max
>> > Out[10]: 1.#INF
>> >
>> > I see that long double in C is 12 bytes wide, and double is the usual 8
>> > bytes.
>>
>> Sorry - sizeof(long double) is 12 using mingw.  I see that long double
>> is the same as double in MS Visual C++.
>>
>> http://en.wikipedia.org/wiki/Long_double
>>
>> but, as expected from the name:
>>
>> In [11]: np.dtype(np.float96).itemsize
>> Out[11]: 12
>>
>
> Hmm, good point. There should not be a float96 on Windows using the MSVC
> compiler, and the longdouble types 'gG' should return float64 and complex128
> respectively. OTOH, I believe the mingw compiler has real float96 types but
> I wonder about library support. This is really a build issue and it would be
> good to have some feedback on what different platforms are doing so that we
> know if we are doing things right.

Is it possible that numpy is getting confused by being compiled with
mingw on top of a visual studio python?

Some further forensics seem to suggest that, despite the fact the math
suggests float96 is float64, the storage format it in fact 80-bit
extended precision:

On OSX 32-bit where float128 is definitely 80 bit precision we see the
sign bit being flipped to show us the beginning of the number:

In [33]: bigbin(np.float128(2**53)-1)
Out[33]: '00000000000000000000000000000000101111111111111101000000001100111111111111111111111111111111111111111111111111111111100000000000'

In [34]: bigbin(-np.float128(2**53)+1)
Out[34]: '00000000000000000000000000000000000000000000000011000000001100111111111111111111111111111111111111111111111111111111100000000000'

I think that's 48 bits of padding followed by the number (bit 49 is
being flipped with the sign).

On windows (well, wine, but I think it's the same):

bigbin(np.float96(2**53)-1)
Out[14]: '000000000000000001000000001100111111111111111111111111111111111111111111111111111111100000000000'
bigbin(np.float96(-2**53)+1)
Out[15]: '000000000000000011000000001100111111111111111111111111111111111111111111111111111111100000000000'

Thanks,

Matthew

<bigbin-definition>
import sys
LE = sys.byteorder == 'little'

import numpy as np

def bigbin(val):
    val = np.asarray(val)
    nbytes = val.dtype.itemsize
    dt = [('f', np.uint8, nbytes)]
    out = [np.binary_repr(el, 8) for el in val.view(dt)['f']]
    if LE:
        out = out[::-1]
    return ''.join(out)