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.
(-2.9802321277472288e-08, 2.9802320833383078e-08)
(7.923312783653103e-09, 0.9999999805091377)
(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)
(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
<Axes: xlabel='d'>
<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
<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
<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
<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
<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)