2
$\begingroup$

I wasn't sure what to title this Question. I tried calculating different forms of a continued fraction, and I found that I get a correct result when I set d to 0 but a wrong result (in the final two terms) if I set d to 0.0

ClearAll["Global`*"];
t = Sqrt[17];
d = 0.0;
n = 11;
a = ConstantArray[0, n];
Do[a[[i]] = Floor[t + d]; t = 1/(t - Floor[t + d]), {i, 1, n}];
a

{4, 8, 8, 8, 8, 8, 8, 8, 8, 7, 4}

I only caught this because I happen to know that this series should end in infinite 8's.

When I try to calculate the form of continued fractions that allows negative terms, I get the correct series from d=1/2 but wrong terms if d=.5. So I gather that this problem has something to do with precision and using rational numbers instead of floats but why is this happening for zero? Is Mathematica adding 0.0 incorrectly?

$\endgroup$
  • 1
    $\begingroup$ The presence of any machine precision number (including 0.0) will result in the calculation being done with machine precision. $\endgroup$ – Bob Hanlon Feb 4 at 21:33
  • 1
    $\begingroup$ By using an extended precision version of 0 you can get a small speedup from using numbers instead of full symbolic expressions but without sacrificing accuracy. Try d=0``1. $\endgroup$ – b3m2a1 Feb 4 at 21:47
  • 2
    $\begingroup$ Also in terms of programming style this is probably faster: d + NestList[1/(# - Floor[# + d]) &, t, 200] // Floor // RepeatedTiming $\endgroup$ – b3m2a1 Feb 4 at 21:47
4
$\begingroup$

Your algorithm is probing the limits of machine precision numbers.

When you add Sqrt[17]+0 you get $\sqrt{17}$ with infinite precision:

Sqrt[17] + 0
(* Sqrt[17] *)

When you add Sqrt[17]+0.0 you get much less precision:

Sqrt[17] + 0.0 // FullForm
(* 4.123105625617661` *)

The number you're effectively dealing with in the d=0.0 case is not $t=\sqrt{17}$ but

SetPrecision[Sqrt[17] + 0.0, ∞]
(* 4642204239785223/1125899906842624 *)

which is the closest available machine-precision number to $\sqrt{17}$.

If you change your code to

t = 4642204239785223/1125899906842624;
d = 0;
n = 11;
a = ConstantArray[0, n];
Do[a[[i]] = Floor[t + d]; t = 1/(t - Floor[t + d]), {i, 1, n}];
a

you get a similar answer to the d=0.0 case:

(* {4, 8, 8, 8, 8, 8, 8, 8, 8, 7, 2} *)

Another way of seeing this is that you're doing a continued-fraction representation, which differs when you convert a number to its machine-precision representation:

ContinuedFraction[Sqrt[17]]
(* {4, {8}} *)
ContinuedFraction[Sqrt[17] // N]
(* {4, 8, 8, 8, 8, 8, 8, 8} *)

Notice that the builtin ContinuedFraction command stops before spitting out "wrong" digits due to machine-precision issues.

$\endgroup$
1
$\begingroup$

Two comments rolled into one.

Try using extended precision numbers like:

d = 0``1

Also try a different programming style to get things cleaner and a bit faster:

t = Sqrt[17];
d = 0``1;
d + NestList[1/(# - Floor[# + d]) &, t, 200] // Floor // RepeatedTiming

{0.0055, {4, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
  8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
   8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
  8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
   8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
  8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
   8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
  8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8,
   8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
  8}}
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.