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

np.int64(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

np.int64(100000000)

Let’s study the distribution of the difference.

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

(np.float64(-2.9802322165650708e-08), np.float64(2.980232005622696e-08))

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

(np.float64(3.9351747416205285e-09), np.float64(0.9999999963626003))

plt.hist(delta, bins=1000)
plot float and double rouding
(array([ 49762.,  50032.,  49564.,  49549.,  49934.,  50051.,  49870.,
        50139.,  50310.,  49896.,  50027.,  49990.,  49855.,  49715.,
        49855.,  49669.,  49855.,  49657.,  49891.,  50002.,  50009.,
        49893.,  49940.,  49784.,  49696.,  50086.,  49904.,  50125.,
        49957.,  50237.,  50067.,  50543.,  49763.,  50064.,  49779.,
        50002.,  49881.,  50026.,  49610.,  50146.,  49905.,  49784.,
        50018.,  50007.,  50189.,  50287.,  49711.,  49909.,  49668.,
        50101.,  49954.,  49985.,  49555.,  49851.,  49799.,  50106.,
        49726.,  50004.,  49993.,  49776.,  49890.,  50114.,  50198.,
        49899.,  49814.,  49945.,  49632.,  49748.,  50111.,  49965.,
        49581.,  50341.,  50206.,  49797.,  50081.,  50062.,  49923.,
        50172.,  49718.,  50220.,  49969.,  49996.,  49800.,  50085.,
        49694.,  50022.,  49763.,  50169.,  50280.,  50164.,  49938.,
        50189.,  50494.,  50010.,  49823.,  49795.,  50462.,  50189.,
        49971.,  49865.,  49635.,  50117.,  49847.,  50138.,  49912.,
        50258.,  50165.,  49792.,  49807.,  49818.,  49724.,  49894.,
        50410.,  49689.,  49854.,  50231.,  49759.,  49744.,  49982.,
        49956.,  50113.,  50258.,  49679.,  50024.,  49819.,  49710.,
        50186.,  50177.,  50185.,  49908.,  50269.,  49869.,  50099.,
        50184.,  50272.,  49902.,  50349.,  50157.,  50061.,  49979.,
        50118.,  50114.,  49947.,  50041.,  49832.,  49713.,  49953.,
        50535.,  49948.,  50077.,  50310.,  49662.,  50373.,  49961.,
        49742.,  50119.,  49780.,  50125.,  50166.,  50301.,  49836.,
        49886.,  50173.,  49864.,  49897.,  50294.,  50152.,  50695.,
        50093.,  49763.,  49837.,  49922.,  49924.,  50348.,  49950.,
        49991.,  49692.,  49438.,  50040.,  49743.,  49739.,  49902.,
        50033.,  49788.,  49787.,  50335.,  50133.,  49806.,  49983.,
        49498.,  49974.,  50033.,  49744.,  50021.,  49855.,  49798.,
        49791.,  49621.,  49935.,  49827.,  50200.,  50096.,  49681.,
        50304.,  49863.,  50400.,  50331.,  50049.,  49936.,  50059.,
        50287.,  50306.,  50065.,  49630.,  50227.,  49831.,  49957.,
        50049.,  50216.,  49986.,  49632.,  50501.,  50107.,  49984.,
        49869.,  49766.,  49691.,  50310.,  49939.,  49812.,  49667.,
        49928.,  49842.,  49882.,  50111.,  49791.,  50470.,  49890.,
        49693.,  50030.,  49889.,  49793.,  49880.,  50063.,  49862.,
        49911.,  50005.,  50131.,  50020.,  50193.,  99739.,  99491.,
       100195., 100380., 100189., 100400.,  99859.,  99319., 100196.,
        99974., 100256.,  99772., 100219., 100377.,  99406., 100198.,
        99346., 100636., 100208., 100171., 100219.,  99548.,  99749.,
       100158.,  99558., 100255.,  99782.,  99496.,  99696., 100151.,
       100096., 100173., 100026.,  99786.,  99930.,  99856.,  99593.,
       100084., 100312., 100026., 100249., 100156., 100427.,  99960.,
       100016.,  99944., 100352., 100191.,  99936.,  99952.,  99698.,
        99852.,  99812., 100119., 100384.,  99842., 100516.,  99829.,
       100134.,  99790., 100074., 100232., 100406.,  99874., 100067.,
       100730., 100073.,  99995., 100056., 100120., 100249.,  99731.,
       100047.,  99865.,  99398.,  99664.,  99421., 100257., 100253.,
       100095.,  99766., 100056.,  99949., 100089., 100197., 100170.,
        99898.,  99590.,  99232., 100433., 100417., 100309.,  99560.,
        99826.,  99629., 100249.,  99802.,  99232., 100154.,  99315.,
        99941., 100041.,  99985., 100164., 100208.,  99674., 100602.,
        99957., 100039., 100135., 100047., 100148., 100165.,  99696.,
       100386., 100349., 100461., 100067., 100045., 100112., 100158.,
        99759.,  99984., 100252., 100464., 150188., 150069., 149979.,
       149128., 150146., 150177., 150151., 149896., 150591., 150620.,
       149570., 150136., 149844., 150332., 149922., 150278., 149977.,
       149791., 149838., 149859., 149766., 149608., 150566., 149450.,
       150261., 149758., 149838., 150039., 150500., 149787., 149538.,
       150130., 149569., 150782., 150475., 149594., 150247., 150081.,
       149282., 149905., 149589., 150105., 149076., 150266., 150321.,
       149854., 149352., 150312., 149996., 150521., 149549., 150079.,
       150837., 150763., 149389., 149793., 149851., 149650., 149943.,
       150457., 149854., 149926., 174931., 199944., 199907., 200094.,
       200133., 199369., 199819., 200232., 199581., 200077., 200120.,
       199688., 199181., 200431., 200577., 200103., 200361., 199856.,
       199837., 200215., 199934., 199728., 199756., 199829., 200230.,
       200399., 200539., 199893., 200114., 200201., 199746., 213052.,
       250250., 249778., 249759., 249507., 249569., 250275., 249795.,
       250261., 250306., 250860., 248900., 249558., 249025., 250413.,
       250235., 281338., 298611., 300798., 299935., 300131., 299456.,
       300148., 299623., 341069., 350318., 350210., 350121., 394603.,
       400300., 447986., 547597., 547169., 449043., 399792., 394929.,
       350966., 349641., 349953., 340994., 300345., 300532., 299202.,
       299999., 299767., 299686., 299450., 281436., 249811., 250521.,
       249911., 249812., 250579., 250320., 250494., 251135., 250158.,
       250014., 249115., 249611., 249782., 250099., 249784., 212967.,
       199967., 200476., 199286., 200280., 199411., 199867., 199654.,
       200056., 199781., 199969., 199378., 199745., 200028., 200040.,
       199900., 201056., 200030., 199850., 200154., 199800., 199658.,
       199499., 200484., 200231., 200596., 199912., 200590., 201046.,
       200701., 199237., 175572., 149965., 149717., 151150., 150096.,
       150369., 149801., 150031., 149415., 150620., 149231., 150192.,
       150476., 149896., 149669., 149498., 150110., 149953., 150034.,
       149704., 149286., 150404., 149981., 149846., 149850., 150040.,
       150305., 149852., 148946., 149772., 150825., 150698., 149495.,
       150292., 149535., 149548., 149873., 149869., 150079., 149865.,
       150063., 149500., 149720., 150277., 149950., 150182., 149502.,
       149432., 150115., 149951., 149154., 149426., 149921., 149628.,
       150063., 149483., 149997., 150018., 150475., 149948., 149665.,
       149949., 149770.,  99278., 100136., 100485.,  99321., 100142.,
       100169., 100162.,  99992., 100310., 100082., 100291., 100122.,
        99915.,  99411.,  99342., 100130., 100043.,  99593.,  99940.,
       100244., 100484., 100235., 100062.,  99891.,  99803., 100682.,
        99862.,  99981.,  99707., 100429., 100450.,  99751.,  99852.,
       100680., 100334., 100091., 100232., 100393., 100143., 100182.,
       100035., 100248.,  99771., 100015.,  99717., 100085., 100785.,
        99815.,  99991.,  99672.,  99834.,  99902., 100512.,  99943.,
        99944., 100609.,  99982., 100009., 100244.,  99949., 100295.,
       100088.,  99736.,  99770., 100282., 100144., 100212.,  99450.,
       100406., 100105.,  99984., 100035.,  99925.,  99378.,  99837.,
        99639.,  99692., 100551., 100780.,  99847., 100208., 100075.,
        99935.,  99828., 100359.,  99534.,  99829.,  99528., 100725.,
        99747.,  99480., 100048.,  99814.,  99623., 100405.,  99759.,
       100035., 100332., 100535.,  99813., 100301., 100331., 100234.,
       100011.,  99689., 100857.,  99670., 100549., 100521., 100508.,
        99897., 100002.,  99857., 100093., 100011.,  99746.,  99479.,
        99888.,  99859., 100127., 100609., 100332.,  99705.,  99410.,
       100371.,  49960.,  50360.,  49815.,  50312.,  49802.,  50065.,
        50019.,  50415.,  49875.,  49856.,  50241.,  49992.,  50064.,
        50386.,  49676.,  50133.,  49965.,  50440.,  50028.,  49751.,
        50190.,  49984.,  50449.,  50255.,  50128.,  49947.,  50334.,
        50260.,  49802.,  49828.,  50247.,  49833.,  50330.,  49793.,
        49382.,  50209.,  50387.,  50228.,  49834.,  50248.,  49920.,
        49830.,  50617.,  49749.,  50210.,  50221.,  50081.,  49816.,
        49561.,  49896.,  49986.,  50157.,  49825.,  50336.,  49861.,
        49902.,  49676.,  49585.,  50150.,  49898.,  50105.,  50154.,
        50127.,  49745.,  49912.,  49751.,  50054.,  49554.,  49395.,
        49773.,  49856.,  50069.,  50054.,  50187.,  50051.,  49867.,
        49931.,  50203.,  49758.,  49930.,  49822.,  49930.,  49932.,
        50232.,  49736.,  50197.,  49384.,  49934.,  50204.,  50036.,
        50274.,  49827.,  50402.,  49950.,  50118.,  49971.,  49948.,
        49615.,  49993.,  49943.,  50203.,  50412.,  50115.,  50732.,
        49821.,  49685.,  50031.,  49892.,  49868.,  49675.,  50287.,
        49706.,  50234.,  50107.,  49982.,  50450.,  50281.,  49960.,
        50192.,  50007.,  50106.,  49884.,  49915.,  50482.,  50342.,
        50046.,  49832.,  49913.,  49581.,  50207.,  50290.,  50337.,
        49979.,  50239.,  49742.,  49628.,  50316.,  50290.,  50126.,
        50274.,  50038.,  50253.,  50417.,  50223.,  50270.,  50106.,
        49814.,  50042.,  50251.,  50265.,  49950.,  50038.,  49899.,
        50020.,  50141.,  50154.,  50048.,  50010.,  49761.,  49965.,
        49852.,  50255.,  50005.,  50124.,  50048.,  49820.,  50169.,
        50026.,  50151.,  49993.,  50144.,  50192.,  50016.,  50640.,
        49446.,  49674.,  50330.,  49728.,  49977.,  49757.,  49926.,
        50380.,  50196.,  49871.,  49869.,  50061.,  49943.,  49825.,
        50430.,  50099.,  49882.,  49694.,  50292.,  49825.,  49988.,
        50039.,  50255.,  49980.,  49963.,  49996.,  49931.,  50039.,
        49633.,  50011.,  50512.,  50405.,  49962.,  50012.,  49949.,
        50091.,  49782.,  49690.,  49817.,  49768.,  50078.,  49872.,
        50504.,  50033.,  49982.,  50244.,  49783.,  49871.,  50037.,
        50388.,  50091.,  49754.,  50300.,  50138.,  50029.,  49439.,
        49813.,  50093.,  49666.,  49788.,  50255.,  49844.,  49956.,
        49886.,  49998.,  49915.,  49653.,  49861.,  49783.,  50063.,
        50205.,  50473.,  49504.,  49982.,  50061.,  49701.]), array([-2.98023222e-08, -2.97427175e-08, -2.96831129e-08, ...,
        2.96831108e-08,  2.97427154e-08,  2.98023201e-08], shape=(1001,)), <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)

(np.float64(1.447638189877587e-11), np.float64(5.960463056453591e-08))

plt.hist(length, bins=50)
plot float and double rouding
(array([  40.,   29.,    0.,   59.,    0.,    0.,  116.,    0.,    0.,
          0.,    0.,    0.,  216.,    0.,    0.,    0.,    0.,    0.,
          0.,    0.,    0.,    0.,    0.,    0.,  480.,    0.,    0.,
          0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,
          0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,    0.,
          0.,    0.,    0.,    0., 1060.]), array([1.44763819e-11, 1.20627947e-09, 2.39808255e-09, 3.58988563e-09,
       4.78168872e-09, 5.97349180e-09, 7.16529488e-09, 8.35709797e-09,
       9.54890105e-09, 1.07407041e-08, 1.19325072e-08, 1.31243103e-08,
       1.43161134e-08, 1.55079165e-08, 1.66997196e-08, 1.78915226e-08,
       1.90833257e-08, 2.02751288e-08, 2.14669319e-08, 2.26587350e-08,
       2.38505381e-08, 2.50423411e-08, 2.62341442e-08, 2.74259473e-08,
       2.86177504e-08, 2.98095535e-08, 3.10013566e-08, 3.21931596e-08,
       3.33849627e-08, 3.45767658e-08, 3.57685689e-08, 3.69603720e-08,
       3.81521751e-08, 3.93439781e-08, 4.05357812e-08, 4.17275843e-08,
       4.29193874e-08, 4.41111905e-08, 4.53029936e-08, 4.64947966e-08,
       4.76865997e-08, 4.88784028e-08, 5.00702059e-08, 5.12620090e-08,
       5.24538121e-08, 5.36456151e-08, 5.48374182e-08, 5.60292213e-08,
       5.72210244e-08, 5.84128275e-08, 5.96046306e-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

np.float64(5.953439408123984e-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.3580954673115855e-08 : 655 987
190 1.131153487543557e-07 : 860 1038
290 1.7264974283559553e-07 : 894 1035
390 2.3218413691683537e-07 : 930 1021
490 2.917185309980752e-07 : 942 1014
590 3.5125292507931505e-07 : 946 1004
690 4.107873191605549e-07 : 949 990
790 4.703217132417947e-07 : 957 998
890 5.298561073230346e-07 : 961 1003
990 5.893905014042744e-07 : 914 952
d ratio total
0 5.953439e-10 0.019507 974
1 1.190688e-09 0.036000 1000
2 1.786032e-09 0.063572 991
3 2.381376e-09 0.083495 1030
4 2.976720e-09 0.069430 965 


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.3580954673115855e-08 : 1025 1025
190 1.131153487543557e-07 : 1021 1021
290 1.7264974283559553e-07 : 965 965
390 2.3218413691683537e-07 : 953 953
490 2.917185309980752e-07 : 1002 1002
590 3.5125292507931505e-07 : 994 994
690 4.107873191605549e-07 : 1000 1000
790 4.703217132417947e-07 : 995 995
890 5.298561073230346e-07 : 979 979
990 5.893905014042744e-07 : 967 967
d ratio total
0 5.953439e-10 1.0 1016
1 1.190688e-09 1.0 1025
2 1.786032e-09 1.0 1019
3 2.381376e-09 1.0 1018
4 2.976720e-09 1.0 956 


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.3580954673115855e-08 : 624 968
190 1.131153487543557e-07 : 882 1053
290 1.7264974283559553e-07 : 862 997
390 2.3218413691683537e-07 : 924 1016
490 2.917185309980752e-07 : 909 975
590 3.5125292507931505e-07 : 951 1008
690 4.107873191605549e-07 : 956 1004
790 4.703217132417947e-07 : 923 972
890 5.298561073230346e-07 : 941 982
990 5.893905014042744e-07 : 1006 1037
d ratio total
0 5.953439e-10 0.024096 996
1 1.190688e-09 0.039722 1007
2 1.786032e-09 0.055441 974
3 2.381376e-09 0.060451 976
4 2.976720e-09 0.093047 978 


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.3580954673115855e-08 : 998 998
190 1.131153487543557e-07 : 986 986
290 1.7264974283559553e-07 : 1003 1003
390 2.3218413691683537e-07 : 1012 1012
490 2.917185309980752e-07 : 993 993
590 3.5125292507931505e-07 : 966 966
690 4.107873191605549e-07 : 976 976
790 4.703217132417947e-07 : 930 930
890 5.298561073230346e-07 : 994 994
990 5.893905014042744e-07 : 1021 1021
d ratio total
0 5.953439e-10 1.0 1000
1 1.190688e-09 1.0 1025
2 1.786032e-09 1.0 983
3 2.381376e-09 1.0 1015
4 2.976720e-09 1.0 1017 


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.3580954673115855e-08 : 975 975
190 1.131153487543557e-07 : 1002 1002
290 1.7264974283559553e-07 : 1013 1013
390 2.3218413691683537e-07 : 997 997
490 2.917185309980752e-07 : 1010 1010
590 3.5125292507931505e-07 : 1007 1007
690 4.107873191605549e-07 : 969 969
790 4.703217132417947e-07 : 996 996
890 5.298561073230346e-07 : 982 982
990 5.893905014042744e-07 : 992 992
d ratio total
0 5.953439e-10 1.0 974
1 1.190688e-09 1.0 979
2 1.786032e-09 1.0 978
3 2.381376e-09 1.0 969
4 2.976720e-09 1.0 1007 


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 29.369 seconds)

Gallery generated by Sphinx-Gallery