# Truncate an infinite continued fraction at order 2000

I want to solve an equation which contains an infinite continued fraction $F(n)$. Then I must (obviously) truncate this continued fraction at $n=2000$.

The problem here is that Mathematica does not display the $2000$ terms of this fraction on the screen. The screen closes directly.

Please, how do I display this fraction up to $n=2000$?

The code is as follows:

F[n_] := 1/(1 + I A x - (n + 1)^2/(4 (n + 1)^2 - 1)
x^2 A (1 - I 2 B x) With[{nplus1 = n + 1}, Hold[F[nplus1]]])

Fold[(#1 /. Hold[F[#2]] :> F[#2]) &, F[1], Range[1, 2000]]


$A$ and $B$ are real constants.

• What do you mean by "The screen closes directly"? Are you getting a front end crash? – m_goldberg Jan 1 '16 at 14:42
• yes, a front end crash. – Betatron Jan 1 '16 at 15:05
• Take a look at this when you do 100 instead of 2000 terms. The answer is already pretty unusable because it has too many terms. By 200, I get a display error that says MaxFormatDepthExceeded – bill s Jan 1 '16 at 22:04
• So I could'nt display the 2000 term of this fraction !! Is there a way to solve the equation with 2000 terms without display these 2000 terms ? Thank's – Betatron Jan 2 '16 at 12:00
• (1) You can always suppress the resulting display of evaluation using ; for example try 2+2; (2) As bill s already said, even for smaller n the answer becomes unwieldy; do you really need the answer for n=2000 or could you truncate F(n) at smaller n? – Sascha Jan 2 '16 at 17:39

A common way, for a long time, of denoting a (generalized) continued fraction is to list the partial numerators and denominators, sometimes with $+$ and fraction bars like this: $$F = b_0+ \frac{a_1}{b_1+}\, \frac{a_2}{b_2+}\, \frac{a_3}{b_3+}\cdots$$ It is also an efficient way to store a continued fraction (cf. ContinuedFraction). What is needed is a way to transform a (finite) list representation

{b0, {a1, b1}, {a2, b2},...}


to an actual fraction. For that there is a standard recurrence, \begin{align} A_{-1}& = 1& B_{-1}& = 0\\ A_0& = b_0& B_0& = 1\\ A_{n+1}& = b_{n+1} A_n + a_{n+1} A_{n-1}& B_{n+1}& = b_{n+1} B_n + a_{n+1} B_{n-1}\, \end{align} where $F = A_n/B_n$ denotes the $n$-th convergent of the continued fraction.

The following constructs the sequence of numerators and denominators in the form

{b0, {a1, b1}, {a2, b2},..., {an, bn}}

Clear[F2, iF2];
iF2[0] = 0;
iF2[1] = {1, 1 + I A x};
iF2[n_] := {-n^2/(4 n^2 - 1) x^2 A (1 - I 2 B x), 1 + I A x};
F2[n_] := Table[iF2[k], {k, 0, n}];


The following function constructs the fraction represented by a generalized continued fraction.

fromGeneralizedContinuedFraction[cf : {_, {_, _} ...}] :=
Divide @@ Last@ Fold[{{0, 1}, #2}.#1 &, {{1, 0}, {First@cf, 1}}, Rest@cf];


Check. Note that the OP's code produces an equivalent c.f. if we set the held F[] to zero.

me = fromGeneralizedContinuedFraction[F2[3]];

op = Fold[(#1 /. Hold[F[#2]] :> F[#2]) &, Hold@F[1], Range[1, 3]];
op = op /. Hold[F[_]] -> 0  // Together // Simplify;

me - op // Together // Simplify
(*  0  *)


Nonetheless the c.f. for F2[2000] makes quite a large fraction and a very long computation when the variables A, B, and x are symbolic. It's somewhat fast for numerical computation, though:

Block[{A = 3., B = 2., x = 0.1},
fromGeneralizedContinuedFraction[F2[2000]] // AbsoluteTiming
]
(*  {0.020405, 0.91971 - 0.28252 I}  *)


As Mr. Wizard supposed in a comment, one can indeed use ContinuedFractionK[] here:

With[{A = 3., B = 2., x = 0.1},
1/(1 + I A x + ContinuedFractionK[-n^2/(4 n^2 - 1) x^2 A (1 - I 2 B x),
1 + I A x, {n, 2, 2000}])]
0.9197103744410972 - 0.28251974414934944 I


However, if what you want is to approximate the value of the infinite CF and not just a truncation, one can certainly do better. Michael displays one approach, using forward evaluation, but a possible danger is that the convergent's numerators and denominators might grow large even if the actual CF's value itself is perfectly modest. (As an easier example, consider the CF for $\phi$, where the numerators and denominators of the convergent are the Fibonacci numbers, which have an exponential growth rate.)

One might thus consider a modification of the forward method where the numerators and denominators are suitably scaled. I am aware of two, but there certainly might be others. Here, for instance, is Steed's algorithm:

With[{A = 3., B = 2., x = 0.1},
b = 1 + I A x; f = h = d = 1/b; k = 2;
While[a = -x^2 A (1 - I 2 B x) k^2/(4 k^2 - 1);
d = 1/(b + a d);
h *= d b - 1; f += h; k++;
Abs[h] > $MachineEpsilon Abs[f]]; f] 0.9197103744410973 - 0.28251974414934944 I  and here is the Lentz-Thompson-Barnett algorithm: With[{A = 3., B = 2., x = 0.1}, b = 1 + I A x; f = c = b; d = 0; k = 2; While[a = -x^2 A (1 - I 2 B x) k^2/(4 k^2 - 1); d = 1/(b + a d); c = b + a/c; h = c d; f *= h; k++; Abs[h - 1] >$MachineEpsilon Abs[f]];
1/f]
0.9197103744410975 - 0.2825197441493497 I


Both are useful, in that Steed is often faster than Lentz-Thompson-Barnett, but the latter is more easily modified for "problematic" cases. Note that in both methods, the number of numerators and denominators generated is much less than the 2000 that was desired by the OP.

• I didn't (and still don't) find this form of ContinuedFraction[] in the docs...(+1) – Michael E2 Jul 18 '16 at 12:38
• @Michael, ContinuedFraction[] and ContinuedFractionK[] are different functions. (The confusion in function names is unfortunate, but then ContinuedFractionK[] came way after ContinuedFraction[].) – J. M. is away Jul 18 '16 at 12:41
• Apparently, I'm going blind. I didn't see the K, nor ContinuedFractionK at the bottom of the CF page. (I was surprised, when I thought this wasn't in M.) – Michael E2 Jul 18 '16 at 12:43
• Don't fret, @Michael; I still inadvertently swap the two from time to time... :) – J. M. is away Jul 18 '16 at 12:49

I don't know if this is really an answer as the question itself seems misguided. (Why would one want to "display" a massive expression even if this could be done?) Nevertheless I would like to comment on the code used by the Betatron. Hold was used for manual control of recursion, however it would be easier to avoid unwanted recursion in the first place by using manual rule application instead. Please consider:

ClearAll[F]

rule =
F[n_] :> 1/(1 + I A x - (n + 1)^2/(4 (n + 1)^2 - 1) x^2 A (1 - I 2 B x) F[n + 1]);

Nest[# /. rule &, F[1], 3] /. F[_] -> 0


$\frac{1}{-\frac{4 A x^2 (1-2 i B x)}{15 \left(-\frac{9 A x^2 (1-2 i B x)}{35 (1+i A x)}+i A x+1\right)}+i A x+1}$

If one provides floating point values for A, B, and x this can easily, though somewhat inefficiently, be extended to 2000 terms:

Block[{A = 3.7, B = -1.2, x = 2.4},
Nest[# /. rule &, F[1], 2000] /. F[_] -> 0
]

-0.0194968 - 0.0955964 I