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

2

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.9802321055427683e-08, 2.9802321721561498e-08)

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

(9.91019621920941e-09, 0.9999999740103884)

plt.hist(delta, bins=1000)
plot float and double rouding
(array([ 49933.,  50211.,  50124.,  50142.,  49965.,  50176.,  49854.,
        49792.,  49637.,  49803.,  50390.,  49874.,  49952.,  49763.,
        49934.,  49917.,  49852.,  49823.,  50244.,  50251.,  49714.,
        50292.,  49838.,  50156.,  50051.,  49874.,  49602.,  50190.,
        50168.,  49603.,  50060.,  49422.,  50202.,  50018.,  50208.,
        49919.,  50089.,  49876.,  49937.,  50208.,  49765.,  50164.,
        49935.,  50250.,  50073.,  49763.,  49961.,  49957.,  50238.,
        50052.,  49872.,  50178.,  49975.,  49901.,  50087.,  50132.,
        50534.,  49974.,  50212.,  50162.,  50001.,  50012.,  49924.,
        50143.,  49589.,  50291.,  50421.,  50204.,  50140.,  50070.,
        50077.,  49752.,  49769.,  49906.,  49941.,  50196.,  50447.,
        50228.,  49677.,  49700.,  49896.,  49729.,  49861.,  49916.,
        50222.,  49947.,  50407.,  50199.,  49961.,  49606.,  49923.,
        49608.,  49605.,  49833.,  49858.,  50042.,  49460.,  50093.,
        49861.,  49859.,  50027.,  50177.,  50414.,  50216.,  50217.,
        50451.,  50023.,  49773.,  50085.,  50331.,  50103.,  50218.,
        50115.,  50055.,  50215.,  49212.,  49979.,  50123.,  49769.,
        50164.,  49629.,  50229.,  50017.,  49724.,  49885.,  49504.,
        49924.,  49672.,  50246.,  50300.,  50135.,  50309.,  49773.,
        50258.,  49531.,  50094.,  49842.,  49952.,  50339.,  50133.,
        49820.,  50115.,  49872.,  50209.,  50264.,  49648.,  49812.,
        50089.,  49675.,  49804.,  50054.,  49976.,  49965.,  50201.,
        49714.,  49717.,  49757.,  49762.,  49965.,  49818.,  49608.,
        49826.,  50247.,  50037.,  49780.,  49844.,  49871.,  49934.,
        49935.,  49730.,  50171.,  49859.,  49852.,  50009.,  50205.,
        49996.,  49492.,  49888.,  49698.,  50011.,  50080.,  49959.,
        50088.,  50406.,  49802.,  49887.,  50188.,  49881.,  50441.,
        50334.,  50136.,  49781.,  49967.,  50177.,  49709.,  50169.,
        50160.,  49741.,  50044.,  50166.,  50229.,  50505.,  50119.,
        50176.,  50145.,  50322.,  50188.,  50091.,  49583.,  49840.,
        50235.,  49955.,  49862.,  50088.,  50080.,  49751.,  50064.,
        50114.,  50041.,  50062.,  49865.,  49821.,  50017.,  49886.,
        49683.,  50120.,  49622.,  50264.,  50064.,  50430.,  49694.,
        49816.,  50011.,  49672.,  50163.,  49868.,  49912.,  49982.,
        50110.,  50099.,  49982.,  49805.,  50119.,  50067.,  50356.,
        50065.,  49595.,  49931.,  49975.,  50233., 100017., 100095.,
        99719., 100001.,  99856.,  99882., 100271.,  99821.,  99861.,
       100335., 100261.,  99933., 100412.,  99757.,  99774., 100073.,
       100186., 100152., 100421.,  99387., 100394.,  99748.,  99789.,
       100323., 100208.,  99969.,  99710., 100245., 100056., 100055.,
        99572.,  99931.,  99888., 100133.,  99933.,  99314., 100031.,
       100083.,  99895., 100598., 100045., 100402., 100273.,  99627.,
       100165., 100125., 100183.,  99995., 100024., 100590., 100052.,
        99799.,  99856., 100354.,  99796.,  99707., 100075., 100373.,
       100179., 100206., 100190., 100080., 100069., 100151., 100133.,
       100119., 100093.,  99969., 100178., 100028.,  99774.,  99875.,
       100629., 100251.,  99605.,  99516., 100156., 100377., 100136.,
       100452.,  99568.,  99742., 100527.,  99929., 100080.,  99926.,
        99749.,  99888.,  99613., 100253.,  99969., 100308., 100041.,
       100323.,  99866.,  99384., 100466.,  99949.,  99784.,  99643.,
       100111., 100286., 100301.,  99711.,  99934., 100118., 100642.,
        99701.,  99512.,  99624., 100356.,  99803.,  99829., 100093.,
        99749., 100002.,  99794.,  99910.,  99911.,  99348., 100181.,
        99540., 100487., 100322.,  99972., 150480., 150304., 150431.,
       150513., 149820., 150113., 150867., 149769., 149829., 150206.,
       149709., 149843., 149940., 150536., 149616., 150406., 150081.,
       149314., 149190., 149867., 149532., 150664., 149127., 150455.,
       150383., 149923., 150148., 149780., 149438., 150602., 149840.,
       149912., 150232., 150353., 150403., 149796., 150108., 150318.,
       149658., 150176., 150041., 149664., 150404., 150544., 149937.,
       149753., 149338., 149691., 150170., 150511., 150224., 150202.,
       149102., 150186., 150036., 149908., 150246., 149879., 150024.,
       150433., 149805., 150276., 175660., 199597., 200557., 200545.,
       199562., 200452., 199952., 199484., 199814., 199956., 200014.,
       200444., 200015., 200399., 199806., 200585., 200201., 200183.,
       200093., 200285., 200318., 199533., 199660., 200249., 199567.,
       200052., 200014., 200393., 199838., 200248., 200887., 211886.,
       250242., 249519., 250284., 249680., 250398., 250632., 250077.,
       249789., 249622., 249599., 250217., 249764., 250530., 250694.,
       249060., 281422., 300561., 299899., 299114., 299654., 300692.,
       300667., 299642., 340225., 350652., 349699., 349983., 394727.,
       399944., 449050., 547320., 547204., 447194., 399792., 395037.,
       350239., 350128., 350125., 340207., 299647., 299787., 300104.,
       299121., 300420., 299964., 299784., 280898., 249849., 249750.,
       249125., 250254., 250816., 249811., 249858., 250340., 249795.,
       250335., 248952., 250597., 250064., 250020., 250558., 212938.,
       200182., 200668., 200037., 198918., 200097., 200152., 200113.,
       199732., 200067., 199304., 199824., 199326., 199630., 200091.,
       198549., 200027., 199177., 199080., 200152., 200064., 199703.,
       199877., 199787., 200239., 199667., 200728., 199703., 200099.,
       200223., 200108., 175418., 149525., 150058., 149313., 149717.,
       149877., 150227., 150624., 149566., 149852., 149739., 149930.,
       150646., 150056., 150751., 150256., 149810., 150198., 149927.,
       149666., 150235., 150179., 150484., 149809., 150985., 150583.,
       150205., 150003., 149722., 149711., 149136., 150083., 150134.,
       150445., 150382., 149502., 150547., 150270., 149586., 150354.,
       150217., 150508., 149911., 150508., 149965., 150282., 149700.,
       150468., 149912., 149288., 150814., 149879., 149906., 149837.,
       149268., 149745., 149763., 150008., 149539., 149638., 150243.,
       149720., 150065.,  99709., 100050.,  99804., 100403.,  99951.,
       100228.,  99925.,  99717., 100311.,  99933., 100606., 100941.,
        99660., 100224., 100098., 100280., 100126., 100344., 100144.,
        99787.,  99704.,  99689.,  99900.,  99671.,  99781., 100756.,
       100190.,  99864.,  99956., 100261., 100658., 100546.,  99691.,
        99954.,  99941.,  99521.,  99371., 100154.,  99969., 100374.,
       100428., 100532., 100156., 100087., 100037., 100168., 100391.,
        99969.,  99958.,  99707., 100271.,  99923., 100117., 100367.,
        99882., 100187., 100187., 100109.,  99161., 100074., 100447.,
        99634., 100459.,  99631., 100198.,  99715., 100291., 100049.,
        99956.,  99902., 100134., 100062., 100762.,  99852.,  99935.,
        99713., 100348.,  99534.,  99812.,  99583., 100295., 100434.,
       100376., 100074., 100072.,  99815., 100015., 100190.,  99633.,
        99684.,  99578.,  99744.,  99951., 100560., 100080.,  99419.,
        99942., 100044., 100209., 100338.,  99617., 100144.,  99658.,
        99571., 100066., 100050., 100128.,  99730., 100116.,  99668.,
       100256., 100700., 100081.,  99314., 100421.,  99602.,  99550.,
       100265., 100519.,  99780., 100039., 100319., 100364.,  99808.,
       100342.,  50334.,  49518.,  50134.,  50084.,  49605.,  49813.,
        50264.,  50027.,  50353.,  50083.,  49771.,  50150.,  50107.,
        50326.,  50056.,  50157.,  49854.,  50064.,  50016.,  50311.,
        49735.,  50175.,  50291.,  50026.,  50611.,  49832.,  50242.,
        49973.,  49612.,  49685.,  50183.,  50437.,  49764.,  50173.,
        50137.,  50085.,  50493.,  49858.,  49741.,  49812.,  49854.,
        49905.,  50109.,  49738.,  49873.,  50228.,  49923.,  49772.,
        50157.,  50287.,  49837.,  50168.,  49692.,  49916.,  49784.,
        49723.,  49604.,  50318.,  49552.,  49733.,  49805.,  49896.,
        50168.,  50096.,  50131.,  50403.,  50082.,  49760.,  49945.,
        49826.,  49846.,  49992.,  50104.,  50084.,  49807.,  50058.,
        49784.,  50089.,  50088.,  49969.,  49640.,  49969.,  49924.,
        49765.,  49847.,  50413.,  49963.,  49999.,  50078.,  49914.,
        50277.,  50109.,  49772.,  49940.,  49822.,  50111.,  50127.,
        50108.,  49892.,  49805.,  49629.,  49965.,  50191.,  50058.,
        49746.,  49879.,  49911.,  50122.,  50381.,  49707.,  49816.,
        49796.,  49936.,  50050.,  49886.,  49762.,  49998.,  49528.,
        49893.,  50268.,  50089.,  50056.,  49915.,  49906.,  49839.,
        49863.,  50229.,  49608.,  50276.,  50200.,  50814.,  50179.,
        50213.,  49913.,  50082.,  49895.,  50388.,  50194.,  49707.,
        50284.,  49935.,  49759.,  49826.,  50205.,  49984.,  49756.,
        49813.,  49873.,  50176.,  49845.,  50159.,  49739.,  50043.,
        49673.,  50112.,  49816.,  49902.,  49901.,  49941.,  49766.,
        49960.,  50187.,  49757.,  49888.,  50046.,  49975.,  50214.,
        49731.,  50004.,  49844.,  49674.,  50043.,  50038.,  49710.,
        50068.,  49795.,  50214.,  49570.,  49724.,  50172.,  49870.,
        50094.,  49766.,  50063.,  50257.,  50075.,  50145.,  50010.,
        49446.,  50264.,  50074.,  49815.,  49818.,  49753.,  50090.,
        50131.,  49989.,  50369.,  50013.,  50093.,  50036.,  50164.,
        49960.,  49835.,  49967.,  49836.,  50045.,  50331.,  50053.,
        50149.,  49874.,  49837.,  49903.,  49890.,  50066.,  50154.,
        50040.,  50075.,  50208.,  49680.,  49853.,  50110.,  49825.,
        49755.,  49622.,  50150.,  50246.,  49876.,  49951.,  50156.,
        49755.,  49871.,  50316.,  49868.,  49127.,  50075.,  50165.,
        49701.,  49947.,  50123.,  50321.,  49498.,  49773.,  50506.,
        49854.,  50044.,  49711.,  50126.,  49910.,  50418.]), array([-2.98023211e-08, -2.97427164e-08, -2.96831118e-08, ...,
        2.96831124e-08,  2.97427171e-08,  2.98023217e-08]), <BarContainer object of 1000 artists>)

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

for i in range(0, 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.455067403371979e-11, 5.960464410925681e-08)

plt.hist(length, bins=50)
plot float and double rouding
(array([ 36.,  21.,   0.,  72.,   0.,   0., 117.,   0.,   0.,   0.,   0.,
         0., 279.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,
         0.,   0., 519.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,
         0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,   0.,
         0.,   0.,   0.,   0.,   0., 956.]), array([1.45506740e-11, 1.20635254e-09, 2.39815441e-09, 3.58995628e-09,
       4.78175815e-09, 5.97356002e-09, 7.16536189e-09, 8.35716375e-09,
       9.54896562e-09, 1.07407675e-08, 1.19325694e-08, 1.31243712e-08,
       1.43161731e-08, 1.55079750e-08, 1.66997768e-08, 1.78915787e-08,
       1.90833806e-08, 2.02751824e-08, 2.14669843e-08, 2.26587862e-08,
       2.38505880e-08, 2.50423899e-08, 2.62341918e-08, 2.74259937e-08,
       2.86177955e-08, 2.98095974e-08, 3.10013993e-08, 3.21932011e-08,
       3.33850030e-08, 3.45768049e-08, 3.57686067e-08, 3.69604086e-08,
       3.81522105e-08, 3.93440123e-08, 4.05358142e-08, 4.17276161e-08,
       4.29194179e-08, 4.41112198e-08, 4.53030217e-08, 4.64948236e-08,
       4.76866254e-08, 4.88784273e-08, 5.00702292e-08, 5.12620310e-08,
       5.24538329e-08, 5.36456348e-08, 5.48374366e-08, 5.60292385e-08,
       5.72210404e-08, 5.84128422e-08, 5.96046441e-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.95118410107176e-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.356065690964584e-08 : 665 1025
190 1.1307249792036345e-07 : 842 1005
290 1.7258433893108105e-07 : 882 990
390 2.3209617994179865e-07 : 910 993
490 2.9160802095251626e-07 : 949 1012
590 3.5111986196323386e-07 : 946 1002
690 4.1063170297395146e-07 : 924 991
790 4.7014354398466907e-07 : 942 978
890 5.296553849953867e-07 : 923 957
990 5.891672260061043e-07 : 958 992
d ratio total
0 5.951184e-10 0.028398 986
1 1.190237e-09 0.033138 1026
2 1.785355e-09 0.043605 1032
3 2.380474e-09 0.068702 1048
4 2.975592e-09 0.065025 1015 


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.356065690964584e-08 : 999 999
190 1.1307249792036345e-07 : 984 984
290 1.7258433893108105e-07 : 954 954
390 2.3209617994179865e-07 : 999 999
490 2.9160802095251626e-07 : 1014 1014
590 3.5111986196323386e-07 : 981 981
690 4.1063170297395146e-07 : 1028 1028
790 4.7014354398466907e-07 : 1018 1018
890 5.296553849953867e-07 : 993 993
990 5.891672260061043e-07 : 996 996
d ratio total
0 5.951184e-10 1.0 983
1 1.190237e-09 1.0 966
2 1.785355e-09 1.0 1034
3 2.380474e-09 1.0 989
4 2.975592e-09 1.0 1022 


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.356065690964584e-08 : 635 1020
190 1.1307249792036345e-07 : 788 973
290 1.7258433893108105e-07 : 879 1018
390 2.3209617994179865e-07 : 899 983
490 2.9160802095251626e-07 : 930 998
590 3.5111986196323386e-07 : 940 993
690 4.1063170297395146e-07 : 964 1007
790 4.7014354398466907e-07 : 943 984
890 5.296553849953867e-07 : 943 982
990 5.891672260061043e-07 : 980 1015
d ratio total
0 5.951184e-10 0.020937 1003
1 1.190237e-09 0.040583 961
2 1.785355e-09 0.042084 998
3 2.380474e-09 0.081755 1003
4 2.975592e-09 0.094303 1018 


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.356065690964584e-08 : 1032 1032
190 1.1307249792036345e-07 : 1007 1007
290 1.7258433893108105e-07 : 1009 1009
390 2.3209617994179865e-07 : 1046 1046
490 2.9160802095251626e-07 : 972 972
590 3.5111986196323386e-07 : 967 967
690 4.1063170297395146e-07 : 989 989
790 4.7014354398466907e-07 : 1011 1011
890 5.296553849953867e-07 : 943 943
990 5.891672260061043e-07 : 985 985
d ratio total
0 5.951184e-10 1.0 987
1 1.190237e-09 1.0 1003
2 1.785355e-09 1.0 1004
3 2.380474e-09 1.0 975
4 2.975592e-09 1.0 1014 


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.356065690964584e-08 : 1025 1025
190 1.1307249792036345e-07 : 982 982
290 1.7258433893108105e-07 : 998 998
390 2.3209617994179865e-07 : 991 991
490 2.9160802095251626e-07 : 962 962
590 3.5111986196323386e-07 : 980 980
690 4.1063170297395146e-07 : 1025 1025
790 4.7014354398466907e-07 : 1013 1013
890 5.296553849953867e-07 : 994 994
990 5.891672260061043e-07 : 968 968
d ratio total
0 5.951184e-10 1.0 981
1 1.190237e-09 1.0 1011
2 1.785355e-09 1.0 991
3 2.380474e-09 1.0 995
4 2.975592e-09 1.0 1001 


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: (1 minutes 1.116 seconds)

Gallery generated by Sphinx-Gallery