Velocity Reviews > Generating equally-spaced floats with least rounding error

# Generating equally-spaced floats with least rounding error

Mark Dickinson
Guest
Posts: n/a

 09-24-2011
On Sep 24, 6:01*pm, Steven D'Aprano <steve
(E-Mail Removed)> wrote:
> Ah, I knew it was too easy!
>
> >>> from fractions import Fraction as F
> >>> start, stop, n = 1, 3.1, 7
> >>> [float(F(start) + i*(F(stop)-F(start))/n) for i in range(n+1)]

>
> [1.0, 1.3, 1.6, 1.9000000000000001, 2.2, 2.5, 2.8000000000000003, 3.1]

I believe that's still giving correctly-rounded results. Note that
the stop value of 3.1 isn't exactly 3.1 here: it's

3.100000000000000088817841970012523233890533447265 625

So the 4th value above is the closest float to 4/7 * 1.0 + 3/7 *
3.100000000000000088817841970012523233890533447265 625.

--
Mark

Arnaud Delobelle
Guest
Posts: n/a

 09-24-2011
On 24 September 2011 18:01, Steven D'Aprano
<(E-Mail Removed)> wrote:
> Mark Dickinson wrote:
>
>> Using Fraction for intermediate calculations actually works perfectly
>> for this, since conversions from float to Fraction are exact, while
>> conversions from Fraction to float are correctly rounded. Â*So if
>> you're using Python, you're not too bothered about efficiency, and you
>> want provably correctly-rounded results, why not use Fraction?
>>
>>>>> from fractions import Fraction
>>>>> start, stop, n = 0.0, 2.1, 7
>>>>> [float(Fraction(start) + i * (Fraction(stop) - Fraction(start)) / n)
>>>>> [for i in range(n+1)]

>> [0.0, 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1]

>
>
> Ah, I knew it was too easy!
>
>>>> from fractions import Fraction as F
>>>> start, stop, n = 1, 3.1, 7
>>>> [float(F(start) + i*(F(stop)-F(start))/n) for i in range(n+1)]

> [1.0, 1.3, 1.6, 1.9000000000000001, 2.2, 2.5, 2.8000000000000003, 3.1]

>>> start, stop, n = 1, 3.1, 7
>>> [((n-i)*start + i*stop)/n for i in range(n+1)]

[1.0, 1.3, 1.5999999999999999, 1.9000000000000001, 2.2, 2.5,
2.8000000000000003, 3.1]

>>>> start, stop, n = -1, 1.1, 7
>>>> [float(F(start) + i*(F(stop)-F(start))/n) for i in range(n+1)]

> [-1.0, -0.7, -0.39999999999999997, -0.09999999999999996,
> 0.20000000000000004, 0.5000000000000001, 0.8, 1.1]

>>> start, stop, n = -1, 1.1, 7
>>> [((n-i)*start + i*stop)/n for i in range(n+1)]

[-1.0, -0.7000000000000001, -0.39999999999999997,
-0.09999999999999996, 0.20000000000000004, 0.5, 0.8, 1.1]

On these examples, using fractions is no better than what I suggested
in my previous post.

--
Arnaud

Stefan Krah
Guest
Posts: n/a

 09-24-2011
Arnaud Delobelle <(E-Mail Removed)> wrote:
> >>>> start, stop, n = -1, 1.1, 7
> >>>> [float(F(start) + i*(F(stop)-F(start))/n) for i in range(n+1)]

> > [-1.0, -0.7, -0.39999999999999997, -0.09999999999999996,
> > 0.20000000000000004, 0.5000000000000001, 0.8, 1.1]

>
> >>> start, stop, n = -1, 1.1, 7
> >>> [((n-i)*start + i*stop)/n for i in range(n+1)]

> [-1.0, -0.7000000000000001, -0.39999999999999997,
> -0.09999999999999996, 0.20000000000000004, 0.5, 0.8, 1.1]
>
> On these examples, using fractions is no better than what I suggested
> in my previous post.

Why not use Decimal if one needs exact endpoints? Something like:

def drange(start, stop, step=1):
if (step == 0):
raise ValueError("step must be != 0")
with localcontext() as ctx:
ctx.traps[Inexact] = True
x = start
cmp = -1 if step > 0 else 1
while x.compare(stop) == cmp:
yield float(x)
x += step

>>> list(drange(Decimal(1), Decimal("3.1"), Decimal("0.3")))

[1.0, 1.3, 1.6, 1.9, 2.2, 2.5, 2.8]
>>> list(drange(Decimal(-1), Decimal("1.1"), Decimal("0.3")))

[-1.0, -0.7, -0.4, -0.1, 0.2, 0.5, 0.8]
>>> list(drange(Decimal(-1), Decimal("1.1"), Decimal("0.1823612873")))

[-1.0, -0.8176387127, -0.6352774254, -0.4529161381, -0.2705548508, -0.0881935635, 0.0941677238, 0.2765290111, 0.4588902984, 0.6412515857, 0.823612873, 1.0059741603]
>>> list(drange(Decimal(-1), Decimal("1.1"), Decimal(1)/3))

Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "<stdin>", line 10, in drange
File "/usr/lib/python3.2/decimal.py", line 1178, in __add__
ans = ans._fix(context)
File "/usr/lib/python3.2/decimal.py", line 1652, in _fix
context._raise_error(Inexact)
File "/usr/lib/python3.2/decimal.py", line 3836, in _raise_error
raise error(explanation)
decimal.Inexact: None

Stefan Krah

f.derainville@gmail.com
Guest
Posts: n/a

 09-24-2011
>>> import numpy
>>> numpy.arange(0, 3, 0.3)

array([ 0. , 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1, 2.4, 2.7])

Terry Reedy
Guest
Posts: n/a

 09-24-2011
On 9/24/2011 9:53 AM, Steven D'Aprano wrote:
> I'm trying to generate a sequence of equally-spaced numbers between a lower
> and upper limit. Given arbitrary limits, what is the best way to generate a
> list of equally spaced floats with the least rounding error for each point?
>
> For example, suppose I want to divide the range 0 to 2.1 into 7 equal
> intervals, then the end-points of each interval are:
>
> (0.0)---(0.3)---(0.6)---(0.9)---(1.2)---(1.5)---(1.---(2.1)
>
> and I'd like to return the values in the brackets. Using Decimal or Fraction
> is not an option, I must use floats. If the exact value isn't representable
> as a float, I'm okay with returning the nearest possible float.
>
> The width of each interval is:
>
> width = (2.1 - 0.0)/7

Calculating width is the fundamental problem. .3 cannot be exactly
represented in binary, and higher multiples with multiply the error.

>
> The relevant points can be calculated in either of two methods:
>
> #1 add multiples of the width to the starting value, 0.
>
> #2 subtract multiples of width from the ending value, 2.1.
>
> (Repeated addition or subtraction should be avoided, as it increases the
> amount of rounding error.)
>
> Mathematically the two are equivalent, but due to rounding, they may not be.
> Here's a run using Python 3.2:
>
>>>> [0.0 + i*width for i in range(]

> [0.0, 0.3, 0.6, 0.8999999999999999, 1.2, 1.5, 1.7999999999999998, 2.1]
>
> The 4th and 7th values have rounding errors, the rest are exact

No they are not. Their errors are just smaller and not visible with 16
digits.

>>> width = (2.1 - 0.0)/7
>>> ['%20.18f' % x for x in [0.0 + i*width for i in range(]]

['0.000000000000000000', '0.299999999999999989', '0.599999999999999978',
'0.899999999999999911', '1.199999999999999956', '1.500000000000000000',
'1.799999999999999822', '2.100000000000000089']

On 9/24/2011 10:18 AM, Arnaud Delobelle wrote:
> ((n-i)*a + i*b)/n for i in range(n+1)

>>> ['%20.18f' % x for x in [((7-i)*0.0 + i*2.1)/7 for i in range(]]

['0.000000000000000000', '0.299999999999999989', '0.599999999999999978',
'0.900000000000000133', '1.199999999999999956', '1.500000000000000000',
'1.800000000000000266', '2.100000000000000089']

In the two places where this disagrees with the previous result, I
believe it is worse. The *0.0 adds nothing. You are still multiplying an
inexactly represented 2.1 by increasingly large numbers. Even 7*2.1/7 is
not exact!

My answer is the same as Mark's. Do exact calculation with fractions
(whether using Fraction or not) and convert to float. I believe the
least common multiple of a.as_integer_ratio[1], b.as_integer_ratio[1],
and n is the common denominator needed to convert the numerators to a
range() output. For the current problem, that is lcm(10,7) = 70.

The best you can do for this example is
>>> ['%20.18f' % (i/10 ) for i in range(0, 22, 3)]

['0.000000000000000000', '0.299999999999999989', '0.599999999999999978',
'0.900000000000000022', '1.199999999999999956', '1.500000000000000000',
'1.800000000000000044', '2.100000000000000089']

Notice that the last place errors are all less than 50, so printing to
16 places will make them appear 'exact'.

Without being fancy (to see than you can cancel 7), you get the same
output with
>>> ['%20.18f' % (i/70 ) for i in range(0, 148, 21)]

['0.000000000000000000', '0.299999999999999989', '0.599999999999999978',
'0.900000000000000022', '1.199999999999999956', '1.500000000000000000',
'1.800000000000000044', '2.100000000000000089']

In the two places where these disagree with the first (your original)
they are *better*, with absolute error in the last places of 22 versus
89 and 44 versus 178 for .9 and 1.8. The last place errors greater than
50 gave you the 'inexact' answers in your post, and indeed, they are not
the best possible.

--
Terry Jan Reedy

Steven D'Aprano
Guest
Posts: n/a

 09-25-2011
Arnaud Delobelle wrote:

> On these examples, using fractions is no better than what I suggested
> in my previous post.

I'm sorry, I can't see your previous post. Perhaps it isn't available on my
news provider?

--
Steven

Steven D'Aprano
Guest
Posts: n/a

 09-25-2011
Terry Reedy wrote:

> On 9/24/2011 9:53 AM, Steven D'Aprano wrote:

[...]
>>>>> [0.0 + i*width for i in range(]

>> [0.0, 0.3, 0.6, 0.8999999999999999, 1.2, 1.5, 1.7999999999999998, 2.1]
>>
>> The 4th and 7th values have rounding errors, the rest are exact

>
> No they are not. Their errors are just smaller and not visible with 16
> digits.

Pardon. I meant that they were as close to exact as is possible in binary
floats. With the exception of 0.8999... and 1.7999... I don't believe any
other float can be closer to the exact value.

I did mention that "If the exact value isn't representable as a float, I'm
okay with returning the nearest possible float."

--
Steven

Steven D'Aprano
Guest
Posts: n/a

 09-25-2011
Chris Angelico wrote:

> On Sun, Sep 25, 2011 at 3:01 AM, Steven D'Aprano
> <(E-Mail Removed)> wrote:
>> Mark Dickinson wrote:
>>
>>> Using Fraction for intermediate calculations actually works perfectly
>>> for this, since conversions from float to Fraction are exact, while
>>> conversions from Fraction to float are correctly rounded. Â*So if
>>> you're using Python, you're not too bothered about efficiency, and you
>>> want provably correctly-rounded results, why not use Fraction?
>>>

>> Ah, I knew it was too easy!

>
> Try using Fraction for the start and stop too:

If you look again at the code I used, I did. I just did it inside the list
comp.

>>>> from fractions import Fraction as F
>>>> start,stop,n = F(0),F(21,10),7
>>>> [float(start+i*(stop-start)/n) for i in range(n+1)]

> [0.0, 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1]
>>>> [float(start+i*(stop-start)/n) for i in range(n+1)]

> [-1.0, -0.7, -0.4, -0.1, 0.2, 0.5, 0.8, 1.1]

Something looks a tad suspicious here... the two list comps are identical,
but give completely different results, even though you don't re-assign
start and stop between calls. You're not editing your results are you?
<wink>

But seriously, you were freaking me out there for a bit. I couldn't see why
pulling the conversion to fraction outside of the list comp was giving
different results. And then it hit me...

>>> 2.1 == F(21, 10)

False

You're testing different numbers from me. Try again with F(2.1) as the stop
value.

--
Steven

Guest
Posts: n/a

 09-25-2011
On Sep 24, 6:53*am, Steven D'Aprano <steve
(E-Mail Removed)> wrote:
> I'm trying to generate a sequence of equally-spaced numbers between a lower
> and upper limit. Given arbitrary limits, what is the best way to generatea
> list of equally spaced floats with the least rounding error for each point?
>
> For example, suppose I want to divide the range 0 to 2.1 into 7 equal
> intervals, then the end-points of each interval are:
>
> (0.0)---(0.3)---(0.6)---(0.9)---(1.2)---(1.5)---(1.---(2.1)
>
> and I'd like to return the values in the brackets. Using Decimal or Fraction
> is not an option, I must use floats. If the exact value isn't representable
> as a float, I'm okay with returning the nearest possible float.

Use numpy.

>>> from numpy import mgrid
>>> seq = mgrid[0.0:2.1:8j]
>>> seq

array([ 0. , 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1])

>>> for x in seq:

print repr(x)

0.0
0.29999999999999999
0.59999999999999998
0.89999999999999991
1.2
1.5
1.7999999999999998
2.1000000000000001

Chris Angelico
Guest
Posts: n/a

 09-25-2011
On Sun, Sep 25, 2011 at 3:31 PM, Steven D'Aprano
<(E-Mail Removed)> wrote:>> If you look again at
the code I used, I did. I just did it inside the list> comp.
You did, but you created a Fraction from a float; my version used
Fraction(21,10) instead of (effectively) Fraction(2.1). My description
was a little sloppy, but what I meant was to use Fraction for the
actual arguments to the "function" that this is implementing.
> Something looks a tad suspicious here... the two list comps are identical,
> but give completely different results, even though you don't re-assign
> start and stop between calls. You're not editing your results are you?
> <wink>

Whooops! I had a whole lot of other commands in the scrollback
(displaying intermediate results and such), and elided the rather
crucial reassignment of parameters along with them!

Fortunately, thanks to my reiplophobia, I still have the original
session up in IDLE. This is even the most recent thing in it.

>>> sys.version

'3.2 (r32:88445, Feb 20 2011, 21:29:02) [MSC v.1500 32 bit (Intel)]'
>>> from fractions import Fraction as F
>>> start,stop,n = F(0),F(21,10),7
>>> [float(start+i*(stop-start)/n) for i in range(n+1)]

[0.0, 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.1]
>>> start, stop, n = F(-1), F(11,10), 7
>>> [float(start+i*(stop-start)/n) for i in range(n+1)]

[-1.0, -0.7, -0.4, -0.1, 0.2, 0.5, 0.8, 1.1]

There, now it's honest. Sorry about that!!

ChrisA