Generating equally-spaced floats with least rounding error

[Terry Reedy]

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.8)---(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(8)]
> [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(8)]]
['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(8)]]
['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