Re: Uniform partition of an interval

[John Smith]

The question was: 
  
   why is length([1.8:0.05:1.9]) = 2 ?

On my pentium linux, and probably many other 64-bit floating point 
systems, octave floating point numbers are represented as:

   sign x 1.(52 binary places of decimals) x 2^(exponent)

The sum in question is  (1.9 - 1.8) / 0.05. We can attempt to tabulate
the hex representations of these decimal numbers:

Decimal number    Hex Representation         Hex meaning of representation
   1.8           1.CCCCCCCCCCCCD x 2^0         1.CCCCCCCCCCCCD000...
   1.9           1.E666666666666 x 2^0         1.E666666666666000...
 1.9 - 1.8       1.9999999999990 x 2^-4        0.1999999999999000...
   0.05          1.999999999999A x 2^-5        0.0CCCCCCCCCCCCD00...

We can now do the calculation, which I think gives

(1.9-1.8)/0.05 = 2 * 1.9999999999990 / 1.999999999999A 
               = 2 * 1/(1+0.000000000000A/1.9999999999990)
     approx    = 2 * (1 - 0.00000000000064)
               = 2 *      0.FFFFFFFFFFFF96
               =          1.FFFFFFFFFFFF2C x 2^0

which gives:

Decimal number    Hex Representation         Hex meaning of representation
 (1.9-1.8)/0.05   1.FFFFFFFFFFFF3 * 2^0         1.FFFFFFFFFFFF3000 

( My machine seems to be giving (1.FFFFFFFFFFFF4 x 2^0); Don't 
  understand! Never mind, we can press on.... )

However 1.FFFFFFFFFFFF3 x 2^0 is definitely less than 2 so 

[1.8:0.05:1.9] gives [ 1.8 1.85 ] rather than [ 1.8 1.85 1.9 ]

As JWE says, octave tries hard to do inteligent_floor(1.FFFFFFFFFFF3) 
hoping that the answer will come to 2. (see Range::nelem_internal and 
related functions in liboctave/Range.cc). But the inteligence is not
great enough.

But during the (1.9-1.8) calculation we have subtracted two large 
similar sized numbers, and lost too much precision.

Could a fix be generated using something like

  double ct = 3.0 * DBL_EPSILON * min(abs(rng_base), abs(rng_limit)) 

in Range::nelem_internal to allow the tollerance to 
grow as precision is lost; or is this just too much of a kludge?

Any thoughts?

John


On Wed, 17 May 2000, John W. Eaton wrote:

> On 17-May-2000, address@hidden <address@hidden> wrote:
> 
> | Hello! 
> | 
> | What is the MAIN reason that 1.8:0.05:1.9 produces [1.8000 1.8500]
> | and not [1.8000 1.8500 1.9000]? 
> | I am using 2.0.14 version of Octave. 
> | Thank you for your answer. 
> | Best regards,
> | Emil Zagar
> 
> I'd guess that the MAIN reason is that there is a bug in the way
> Octave is trying (very hard) to compute the correct number of elements
> for ranges.  If you're in a debugging mood, the code to look at is in
> the Range::nelem_internal and related functions in liboctave/Range.cc.
> 
> jwe
> 
> -- 
> www.che.wisc.edu/octave | Thanking you in advance.  This sounds as if the
> www.che.wisc.edu/~jwe   | writer meant, "It will not be worth my while to
>                         | write to you again."        -- Strunk and White
> 
> 
> 
> -----------------------------------------------------------------------
> Octave is freely available under the terms of the GNU GPL.
> 
> Octave's home on the web:  http://www.che.wisc.edu/octave/octave.html
> How to fund new projects:  http://www.che.wisc.edu/octave/funding.html
> Subscription information:  http://www.che.wisc.edu/octave/archive.html
> -----------------------------------------------------------------------
> 







-----------------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.che.wisc.edu/octave/octave.html
How to fund new projects:  http://www.che.wisc.edu/octave/funding.html
Subscription information:  http://www.che.wisc.edu/octave/archive.html
-----------------------------------------------------------------------