Re: ContinuedFraction broken?
- To: mathgroup at smc.vnet.net
- Subject: [mg21551] Re: ContinuedFraction broken?
- From: Mark Sofroniou <marks at wolfram.com>
- Date: Sat, 15 Jan 2000 02:03:57 -0500 (EST)
- Organization: Wolfram Research Inc
- References: <853u1c$8g3@smc.vnet.net>
- Sender: owner-wri-mathgroup at wolfram.com
The input number 0.9 is not representable exactly in binary (binary
is used to represent numbers internally in Mathematica). When a number
is represented in a different base it often contains a representation
error because of the finite nature of floating point arithmetic.
To see this the following function shows that 9/10 has an infinite
(periodic) representation.
In[1]:= RealDigits[9/10, 2]
Out[1]= {{1, {1, 1, 0, 0}}, 0}
In contrast, the floating point number 0.9 has a finite representation.
In[2]:= RealDigits[0.9, 2]
Out[2]= {{1, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1,
> 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1,
> 1, 0, 0, 1, 1, 1, 0}, 0}
You can see which rational number 0.9 actually corresponds to as follows:
In[3]:= SetPrecision[0.9, Infinity]
8106479329266893
Out[3]= ----------------
9007199254740992
Because of representation errors, it is not reasonable to expect that
the continued fraction of an exact number will always coincide with
the continued fraction of a floating point approximation. We want to
know how the error is propagated as the computation of a continued
fraction proceeds, so that incorrect partial quotients are not generated.
The idea behind the error control used in ContinuedFraction
is due to D. H. Lehmer and is described in Donald Knuth's
The Art of Computer Programming: Seminumerical Algorithms,
Volume 2, Addison Wesley. The idea is essentially as follows.
If we compute the continued fraction of two nearby numbers w, y
which bound an input x as w < x < y then we can be sure that
the partial quotients of x are correct as long as the partial
quotients of the neighbouring numbers w and y coincide.
In the following examples only the first two partial quotients of the
inputs coincide and this is the reason that ContinuedFraction only
generates two partial quotients for the input 0.9.
In[4]:= ContinuedFraction[SetPrecision[0.9, Infinity]]
Out[4]= {0, 1, 9, 450359962737049, 2}
This computes the continued fraction of the rational representation
of the nearest two floating point numbers such that w < x=0.9 < y.
In[5]:= ContinuedFraction[SetPrecision[0.9 + $MachineEpsilon/2, Infinity]]
Out[5]= {0, 1, 9, 75059993789508, 6}
In[6]:= ContinuedFraction[SetPrecision[0.9 - $MachineEpsilon/2, Infinity]]
Out[6]= {0, 1, 8, 1, 112589990684261, 2}
Mark Sofroniou,
Research and Development,
Wolfram Research Inc.
Doug Nelson wrote:
> I've been using Mathematica 4.0 for a couple of weeks now (meaning I'm a
> newbie) and have discovered what I think may be a bug in the ContinuedFraction
> built-in function.
>
> In[116]:=
> ContinuedFraction[0.9]
>
> Out[116]=
> {0, 1} * incorrect *
>
> In[120]:=
> ContinuedFraction[0.90000000001]
>
> Out[120]=
> {0, 1, 9} * correct *
>
> In[121]:=
> ContinuedFraction[9/10]
>
> Out[121]=
> {0, 1, 9} * correct *
>
> The docs state that 0.9 is expanded to cover the precision of the input
> (but it's not going far enough). The results are all the same if you
> specify a number of terms (say ",20]").
>
> Is this somehow pilot error?
> Doug