Note

Go to the end to download the full example code.

Float Conversion

I came up with the following question (float64)x < (float64)y \Longrightarrow (float32) x < (float32)y? What is the probability this holds?

Probability (float64)x == (float32)x

Let’s evaluate how many time we draw a random double number equal to its float conversion.

import random
import numpy
import pandas
import matplotlib.pyplot as plt


rnd = numpy.random.random(100000000)
rnd.shape, rnd.dtype

((100000000,), dtype('float64'))

rnd32 = rnd.astype(numpy.float32).astype(numpy.float64)
equal = (rnd == rnd32).sum()
equal

1

It is very low. Let’s check the reverse is true.

rnd32b = rnd32.astype(numpy.float64).astype(numpy.float32)
equal = (rnd32b == rnd32).sum()
equal

100000000

Let’s study the distribution of the difference.

delta = rnd - rnd32
numpy.min(delta), numpy.max(delta)

(-2.9802321277472288e-08, 2.9802320833383078e-08)

numpy.min(rnd), numpy.max(rnd)

(7.923312783653103e-09, 0.9999999805091377)

plt.hist(delta, bins=1000)
plot float and double rouding
(array([ 50248.,  50276.,  50052.,  50209.,  49826.,  50076.,  49742.,
        50010.,  49611.,  49698.,  49760.,  49748.,  49729.,  50136.,
        50004.,  50099.,  49640.,  49693.,  50177.,  49814.,  50266.,
        49853.,  49989.,  50108.,  50354.,  49498.,  49575.,  49922.,
        50058.,  49857.,  50064.,  49761.,  50299.,  50236.,  49821.,
        50052.,  49823.,  49932.,  49656.,  49782.,  50087.,  49633.,
        49804.,  49932.,  50109.,  49970.,  50024.,  50082.,  49457.,
        50155.,  49945.,  49463.,  49896.,  49963.,  50004.,  49960.,
        50032.,  49767.,  50007.,  49805.,  49999.,  50136.,  49612.,
        50339.,  49794.,  50080.,  50349.,  49945.,  50455.,  50128.,
        50318.,  49705.,  49801.,  50167.,  50004.,  50580.,  50384.,
        49707.,  50384.,  50034.,  50324.,  50099.,  49732.,  50233.,
        49778.,  49879.,  50374.,  50527.,  50047.,  50418.,  49967.,
        49980.,  49888.,  49917.,  49789.,  49673.,  49713.,  50094.,
        49950.,  50150.,  49966.,  49876.,  49946.,  49574.,  49856.,
        50122.,  50310.,  50116.,  49873.,  49698.,  49835.,  49701.,
        49887.,  49955.,  49949.,  50203.,  49832.,  50321.,  49710.,
        50491.,  50066.,  50138.,  50359.,  50176.,  50116.,  50162.,
        49882.,  49908.,  49875.,  49960.,  50120.,  49828.,  50435.,
        49835.,  50247.,  49894.,  50119.,  49937.,  50014.,  49952.,
        49978.,  50058.,  49879.,  49696.,  49721.,  50157.,  50024.,
        49967.,  50149.,  50514.,  49990.,  49621.,  49934.,  50251.,
        49769.,  50109.,  50050.,  49662.,  50328.,  49715.,  49898.,
        49685.,  49803.,  50083.,  50052.,  50000.,  49461.,  49678.,
        49863.,  50338.,  49931.,  50159.,  49783.,  49541.,  50053.,
        49733.,  49747.,  50014.,  50056.,  49954.,  50023.,  49931.,
        50170.,  49992.,  50127.,  49998.,  50327.,  49635.,  49542.,
        50156.,  49670.,  49929.,  49896.,  49843.,  49557.,  50210.,
        49592.,  49739.,  50590.,  49749.,  49785.,  50260.,  49705.,
        50151.,  49970.,  49750.,  49979.,  50135.,  49704.,  50184.,
        49872.,  49985.,  49824.,  49955.,  50033.,  50260.,  50181.,
        50230.,  50412.,  50106.,  49627.,  49807.,  49947.,  50297.,
        49873.,  49778.,  49982.,  49805.,  49682.,  50067.,  50125.,
        50476.,  50279.,  50064.,  49726.,  49704.,  49845.,  50176.,
        49619.,  50111.,  49891.,  50465.,  50057.,  50230.,  49984.,
        49941.,  49681.,  49566.,  49916.,  49836.,  99999., 100455.,
       100355.,  99960.,  99963.,  98834., 100255., 100215., 100245.,
       100166., 100008., 100232.,  99860., 100209.,  99613.,  99782.,
       100016.,  99580., 100105., 100219., 100079.,  99604.,  99527.,
        99915.,  99970., 100235., 100198., 100127.,  99685., 100536.,
       100107., 100297.,  99753., 100469.,  99798., 100343., 100742.,
        99373., 100262., 100203.,  99421.,  99887.,  99927.,  99399.,
        99471., 100320., 100000.,  99974.,  99746., 100175., 100059.,
       100037., 100237.,  99993., 100006., 100112., 100255.,  99617.,
        99736., 100002., 100048.,  99824.,  99782.,  99889.,  99827.,
        99990., 100069.,  99645.,  99563., 100009.,  99679.,  99684.,
       100046., 100220., 100553.,  99998., 100351., 100036., 101059.,
       100142.,  99886.,  99202.,  99912.,  99986., 100461.,  99454.,
        99700.,  99458., 100121.,  99486., 100217., 100249.,  99792.,
        99879., 100343., 100172., 100196.,  99611.,  99848., 100106.,
       100013.,  99935., 100454., 100272.,  99598., 100205., 100133.,
        99558.,  99821., 100123., 100652.,  99797., 100416., 100235.,
        99855., 100111.,  99819.,  99661.,  99800.,  99628., 100265.,
       100116.,  99630.,  99998.,  99485., 149902., 150242., 150276.,
       149290., 150426., 150161., 149309., 150374., 150186., 150091.,
       150003., 150372., 149603., 150283., 150584., 150790., 150483.,
       150864., 150385., 149691., 149989., 150124., 149672., 149732.,
       150019., 149607., 150002., 150192., 149922., 150124., 149666.,
       150763., 150100., 149815., 149434., 150195., 149678., 149763.,
       149759., 149983., 150189., 150442., 150443., 149526., 150206.,
       150533., 150271., 150379., 149422., 149648., 150280., 149944.,
       150686., 150121., 150265., 149553., 150385., 149668., 149729.,
       150059., 150029., 150112., 175159., 199518., 199392., 199214.,
       199752., 199868., 200022., 200125., 199299., 200021., 200465.,
       199476., 199290., 199590., 199765., 199293., 200563., 200461.,
       199647., 200175., 200136., 199828., 199948., 200258., 200096.,
       199940., 199766., 200063., 199594., 199893., 200408., 213276.,
       250226., 250679., 249512., 249842., 250148., 249336., 250439.,
       249624., 249574., 250140., 250302., 250572., 249868., 250660.,
       248613., 281517., 299657., 300446., 300652., 300655., 299816.,
       300876., 299223., 339692., 349983., 350774., 350034., 394989.,
       399893., 446860., 547911., 548166., 448215., 400532., 395899.,
       349972., 349722., 350753., 341078., 299777., 300257., 300759.,
       300436., 300503., 300153., 299809., 281599., 249562., 249282.,
       249261., 250186., 249725., 249276., 249920., 250036., 249551.,
       250440., 250495., 248745., 249038., 249829., 250824., 212448.,
       200150., 199626., 200296., 199640., 200611., 200321., 200257.,
       200262., 199241., 199569., 199398., 200281., 199803., 201379.,
       200202., 200082., 200062., 199800., 199749., 200276., 200068.,
       199943., 200247., 200219., 200659., 199802., 199395., 199661.,
       199887., 199414., 174617., 149870., 150176., 150167., 149732.,
       150805., 150225., 150394., 149693., 149810., 150012., 149400.,
       150106., 149560., 150549., 149724., 149456., 149757., 150191.,
       150053., 149313., 149866., 149543., 149729., 150850., 149622.,
       150391., 150397., 150387., 150518., 149840., 149428., 149993.,
       149748., 150242., 150217., 149804., 151182., 150120., 149994.,
       150999., 150083., 149896., 150355., 150470., 149519., 150455.,
       149983., 150088., 149989., 150694., 149062., 150019., 150709.,
       149856., 150137., 150003., 149927., 150206., 150471., 149479.,
       149427., 150116.,  99609.,  99979., 100199.,  99970.,  99543.,
       100438.,  99921., 100624.,  99899.,  99799.,  99581., 100389.,
        99516., 100420., 100361., 100051., 100114.,  99505., 100268.,
       100150., 100162., 100419.,  99890.,  99499., 100234.,  99850.,
        99579.,  99740., 100256.,  99890., 100095., 100598., 100066.,
        99376.,  99650., 100216., 100262., 100379.,  99982., 100445.,
        99743., 100251., 100240., 100460., 100215.,  99735.,  99653.,
       100289.,  99652., 100071., 100011.,  99768.,  99631., 100307.,
       100061.,  99819.,  99628., 100121., 100274.,  99754.,  99799.,
        99894., 100269., 100235., 100107.,  99747.,  99993.,  99618.,
        99942., 100208.,  99815.,  99954.,  99890., 100077., 100071.,
       100406., 100341.,  99948., 100470.,  99973., 100180.,  99761.,
        99953.,  99703., 100542.,  99797., 100130., 100224.,  99823.,
        99497., 100303., 100167., 100547.,  99736.,  99710., 100169.,
       100130., 100078., 100021., 100036., 100407.,  99985.,  99599.,
        99831.,  99475.,  99799.,  99758.,  99862., 100018., 100182.,
       100126.,  99666., 100086.,  99690., 100339., 100040.,  99916.,
       100161.,  99667., 100083., 100233., 100533., 100033.,  99686.,
       100144.,  50127.,  49796.,  50092.,  49844.,  49795.,  49882.,
        49867.,  50088.,  49988.,  49995.,  50197.,  49735.,  49999.,
        49960.,  50092.,  50158.,  50152.,  49736.,  50200.,  49729.,
        50087.,  50202.,  50409.,  50318.,  50288.,  49753.,  49783.,
        49861.,  49835.,  50199.,  50216.,  50089.,  49778.,  49752.,
        50127.,  49625.,  49774.,  49924.,  49785.,  50129.,  49917.,
        50080.,  49795.,  50079.,  49956.,  49966.,  50084.,  50199.,
        50134.,  49850.,  50207.,  50015.,  50247.,  50306.,  50193.,
        50020.,  49967.,  50311.,  49706.,  50156.,  49581.,  50160.,
        50262.,  50039.,  50212.,  50358.,  49886.,  50198.,  49917.,
        50278.,  49900.,  49896.,  50347.,  49558.,  49812.,  49699.,
        49702.,  50093.,  49215.,  50354.,  49962.,  49834.,  50417.,
        50365.,  50081.,  49911.,  50175.,  49918.,  50288.,  50194.,
        49885.,  49825.,  49920.,  50176.,  50273.,  50073.,  49756.,
        49977.,  50081.,  50212.,  49799.,  49732.,  50036.,  50170.,
        50497.,  50086.,  50043.,  49744.,  50179.,  49986.,  49862.,
        49869.,  49841.,  50030.,  50252.,  49971.,  49862.,  49924.,
        49622.,  50125.,  50248.,  50026.,  50078.,  49878.,  50227.,
        49976.,  49672.,  50230.,  50030.,  50440.,  50080.,  49890.,
        49995.,  50268.,  50276.,  49433.,  50177.,  49935.,  50172.,
        50240.,  49911.,  49760.,  49659.,  50047.,  50050.,  49879.,
        50251.,  49962.,  50252.,  50120.,  50165.,  49797.,  49800.,
        49993.,  49992.,  49815.,  49979.,  50174.,  49762.,  50170.,
        49885.,  50028.,  49936.,  50124.,  50129.,  49843.,  50030.,
        49649.,  50035.,  50079.,  50381.,  50023.,  49969.,  50119.,
        50124.,  50163.,  49823.,  50186.,  49671.,  49799.,  49778.,
        50154.,  50008.,  50107.,  50002.,  49957.,  49947.,  50072.,
        49942.,  50257.,  50476.,  49900.,  50128.,  49933.,  49942.,
        49883.,  50007.,  50153.,  49782.,  50322.,  49822.,  50053.,
        50008.,  49932.,  50020.,  50336.,  49989.,  49879.,  50132.,
        50149.,  50116.,  50199.,  49880.,  50256.,  50384.,  49641.,
        50004.,  50068.,  50115.,  50254.,  49575.,  50579.,  50388.,
        49726.,  50036.,  49996.,  49896.,  49545.,  49998.,  49900.,
        49874.,  49699.,  50072.,  50034.,  50266.,  50332.,  50062.,
        49920.,  50015.,  50073.,  50322.,  50217.,  49799.,  49808.,
        49946.,  49792.,  49679.,  49781.,  49980.,  50099.]), array([-2.98023213e-08, -2.97427166e-08, -2.96831120e-08, ...,
        2.96831115e-08,  2.97427162e-08,  2.98023208e-08]), <BarContainer object of 1000 artists>)

We finally check that double operations between float numpers remain floats.

for i in range(100000):
    i, j = random.randint(0, len(rnd32) - 1), random.randint(0, len(rnd32) - 1)
    d32 = numpy.float64(rnd32[i] * rnd32[j])
    d64 = numpy.float64(rnd32[i]) * numpy.float64(rnd32[j])
    if d32 != d64:
        raise AssertionError(
            "Issue with somme={0} = {1} + {2}".format(
                rnd32[i] + rnd32[j], rnd32[i], rnd32[j]
            )
        )

Interval length distribution

Let’s imagine now we want to define an intervalle in which a double is converted to the same float. Let’s find out about it length.

def find_interval(x):
    dx = numpy.abs(x - numpy.float32(x))  # usually not zero
    dx /= 100
    f = numpy.float32(x)
    x1 = x
    while numpy.float32(x1) == f:
        x1 -= dx
    x2 = x
    while numpy.float32(x2) == f:
        x2 += dx
    return x1 + dx, x2 - dx


length = numpy.zeros((2000,))
for i in range(length.shape[0]):
    x = rnd[i]
    x1, x2 = find_interval(x)
    length[i] = x2 - x1

min(length), max(length)

(1.441327634500722e-11, 5.9604616797770404e-08)

plt.hist(length, bins=50)
plot float and double rouding
(array([  45.,   27.,    0.,   75.,    0.,    0.,  111.,    0.,    0.,
          0.,    0.,    0.,  220.,    0.,    0.,    0.,    0.,    0.,
          0.,    0.,    0.,    0.,    0.,    0.,  511.,    0.,    0.,
          0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,
          0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,
          0.,    0.,    0.,    0., 1011.]), array([1.44132763e-11, 1.20621735e-09, 2.39802142e-09, 3.58982549e-09,
       4.78162956e-09, 5.97343363e-09, 7.16523770e-09, 8.35704177e-09,
       9.54884584e-09, 1.07406499e-08, 1.19324540e-08, 1.31242581e-08,
       1.43160621e-08, 1.55078662e-08, 1.66996703e-08, 1.78914743e-08,
       1.90832784e-08, 2.02750825e-08, 2.14668865e-08, 2.26586906e-08,
       2.38504947e-08, 2.50422988e-08, 2.62341028e-08, 2.74259069e-08,
       2.86177110e-08, 2.98095150e-08, 3.10013191e-08, 3.21931232e-08,
       3.33849272e-08, 3.45767313e-08, 3.57685354e-08, 3.69603395e-08,
       3.81521435e-08, 3.93439476e-08, 4.05357517e-08, 4.17275557e-08,
       4.29193598e-08, 4.41111639e-08, 4.53029680e-08, 4.64947720e-08,
       4.76865761e-08, 4.88783802e-08, 5.00701842e-08, 5.12619883e-08,
       5.24537924e-08, 5.36455964e-08, 5.48374005e-08, 5.60292046e-08,
       5.72210087e-08, 5.84128127e-08, 5.96046168e-08]), <BarContainer object of 50 artists>)

So we can approximate this interval by something like this:

ql = numpy.sort(length)[int(length.shape[0] * 0.8)]
ql

5.951835202466782e-08

An answer to the initial question

Let’s estimate \mathbb{P}\left(x_{64} < y_{64} \Longrightarrow x_{32} < y_{32} \; | \; |x-y| \leqslant d\right) ?

def inf_strict(x, y):
    f1 = x < y
    f2 = numpy.float32(x) < numpy.float32(y)
    return f1, f2


def count_events(fct):
    rows = []
    for di in range(1, 1001):
        d = di * ql / 100
        total = 0
        ok = 0
        rnd = numpy.random.random((2000 * 3,))
        for i in range(0, rnd.shape[0], 3):
            s = -1 if rnd[i + 2] < 0.5 else 1
            x, y = rnd[i], rnd[i] + rnd[i + 1] * d * s
            f1, f2 = fct(x, y)
            if f1:
                total += 1
                if f2:
                    ok += 1
        if (di + 10) % 100 == 0:
            print(di, d, ":", ok, total)
        rows.append(dict(d=d, ratio=ok * 1.0 / total, total=total))

    return pandas.DataFrame(rows)


df = count_events(inf_strict)
df.head()

90 5.3566516822201035e-08 : 651 1015
190 1.1308486884686886e-07 : 790 975
290 1.726032208715367e-07 : 863 977
390 2.321215728962045e-07 : 875 959
490 2.9163992492087233e-07 : 942 1015
590 3.5115827694554015e-07 : 935 993
690 4.1067662897020797e-07 : 966 1014
790 4.701949809948758e-07 : 925 965
890 5.297133330195436e-07 : 966 1000
990 5.892316850442114e-07 : 993 1030
d ratio total
0 5.951835e-10 0.023692 1013
1 1.190367e-09 0.054822 985
2 1.785551e-09 0.066084 1029
3 2.380734e-09 0.073242 1024
4 2.975918e-09 0.104288 1026 


df.plot(x="d", y="ratio")
plot float and double rouding
<Axes: xlabel='d'>

df.plot(x="d", y="ratio", logx=True)
plot float and double rouding
<Axes: xlabel='d'>

An answer to a similar question: what about not strict comparison?

Let’s estimate \mathbb{P}\left(x_{64} \leqslant y_{64} \Longrightarrow x_{32} \leqslant y_{32} \; | \; |x-y| \leqslant d\right) ?

def inf_equal(x, y):
    f1 = x <= y
    f2 = numpy.float32(x) <= numpy.float32(y)
    return f1, f2


df2 = count_events(inf_equal)
df2.head()

90 5.3566516822201035e-08 : 963 963
190 1.1308486884686886e-07 : 977 977
290 1.726032208715367e-07 : 979 979
390 2.321215728962045e-07 : 1009 1009
490 2.9163992492087233e-07 : 925 925
590 3.5115827694554015e-07 : 1022 1022
690 4.1067662897020797e-07 : 990 990
790 4.701949809948758e-07 : 980 980
890 5.297133330195436e-07 : 1018 1018
990 5.892316850442114e-07 : 988 988
d ratio total
0 5.951835e-10 1.0 1002
1 1.190367e-09 1.0 968
2 1.785551e-09 1.0 989
3 2.380734e-09 1.0 1006
4 2.975918e-09 1.0 968 


ax = df.plot(x="d", y="ratio", logx=True, label="<")
df2.plot(x="d", y="ratio", logx=True, label="<=", ax=ax)
plot float and double rouding
<Axes: xlabel='d'>

def sup_strict(x, y):
    f1 = x > y
    f2 = numpy.float32(x) > numpy.float32(y)
    return f1, f2


df3 = count_events(sup_strict)
df3.head()

90 5.3566516822201035e-08 : 631 1010
190 1.1308486884686886e-07 : 825 996
290 1.726032208715367e-07 : 851 972
390 2.321215728962045e-07 : 916 985
490 2.9163992492087233e-07 : 884 952
590 3.5115827694554015e-07 : 956 1027
690 4.1067662897020797e-07 : 894 947
790 4.701949809948758e-07 : 941 990
890 5.297133330195436e-07 : 973 1007
990 5.892316850442114e-07 : 984 1014
d ratio total
0 5.951835e-10 0.016244 985
1 1.190367e-09 0.040236 1019
2 1.785551e-09 0.052419 992
3 2.380734e-09 0.061000 1000
4 2.975918e-09 0.082234 985 


ax = df.plot(x="d", y="ratio", logx=True, label="<")
df2.plot(x="d", y="ratio", logx=True, label="<=", ax=ax)
df3.plot(x="d", y="ratio", logx=True, label=">", ax=ax)
plot float and double rouding
<Axes: xlabel='d'>

def sup_equal(x, y):
    f1 = x >= y
    f2 = numpy.float32(x) >= numpy.float32(y)
    return f1, f2


df4 = count_events(sup_equal)
df4.head()

90 5.3566516822201035e-08 : 1004 1004
190 1.1308486884686886e-07 : 988 988
290 1.726032208715367e-07 : 971 971
390 2.321215728962045e-07 : 1003 1003
490 2.9163992492087233e-07 : 967 967
590 3.5115827694554015e-07 : 998 998
690 4.1067662897020797e-07 : 1019 1019
790 4.701949809948758e-07 : 1006 1006
890 5.297133330195436e-07 : 990 990
990 5.892316850442114e-07 : 978 978
d ratio total
0 5.951835e-10 1.0 1004
1 1.190367e-09 1.0 996
2 1.785551e-09 1.0 966
3 2.380734e-09 1.0 990
4 2.975918e-09 1.0 1027 


ax = df.plot(x="d", y="ratio", logx=True, label="<")
df2.plot(x="d", y="ratio", logx=True, label="<=", ax=ax)
df3.plot(x="d", y="ratio", logx=True, label=">", ax=ax)
df4.plot(x="d", y="ratio", logx=True, label=">=", ax=ax)
plot float and double rouding
<Axes: xlabel='d'>

def inf_strict_neg(x, y):
    f1 = (-x) >= (-y)
    f2 = (-numpy.float32(x)) >= (-numpy.float32(y))
    return f1, f2


dfn = count_events(inf_strict_neg)
dfn.head()

90 5.3566516822201035e-08 : 978 978
190 1.1308486884686886e-07 : 957 957
290 1.726032208715367e-07 : 1043 1043
390 2.321215728962045e-07 : 1003 1003
490 2.9163992492087233e-07 : 1037 1037
590 3.5115827694554015e-07 : 1015 1015
690 4.1067662897020797e-07 : 994 994
790 4.701949809948758e-07 : 961 961
890 5.297133330195436e-07 : 977 977
990 5.892316850442114e-07 : 996 996
d ratio total
0 5.951835e-10 1.0 997
1 1.190367e-09 1.0 1007
2 1.785551e-09 1.0 947
3 2.380734e-09 1.0 985
4 2.975918e-09 1.0 1023 


ax = df.plot(x="d", y="ratio", logx=True, label="<")
dfn.plot(x="d", y="ratio", logx=True, label="-1 x >=", ax=ax)
plot float and double rouding
<Axes: xlabel='d'>

Conclusion

The result is expected. As soon as two float are rounded to the same value, the strict inequality no longer holds. However, if you need to write a code which has to handle double and float (in a template for example), you should use not strict inequalities. It is easier to compare the results but you should read some article like Is < faster than <=?. According to Processing costs of non-strict versus strict comparison, < is 5-10% faster than <=.

Total running time of the script: (0 minutes 36.898 seconds)

Gallery generated by Sphinx-Gallery