3
$\begingroup$

I am interested in calculating the following series

$$ \sum_{k=0}^{+\infty}{ \frac{ x^k }{\, a (a+1) \cdots (a+k) \,} } $$

using this continued fraction (I expect) equivalent form:

$$ \cfrac{1}{a + \cfrac{-ax}{a+1 + \cfrac{x}{a+2 + \cfrac{-(a+1)x}{a+3 + \cfrac{2x}{a+4 + \cfrac{-(a+2)x}{a+5 + \cfrac{3x}{a+6 + \cdots}}}}}}} \, \mathrm{,} $$

which can be expressed as

$$ \cfrac{f_0}{g_0 + \cfrac{f_1}{g_1 + \cfrac{f_2}{g_2 + \cdots}}} \, \mathrm{,} $$

where

$$ f_k = \left\{ \begin{aligned} 1, &&& \text{if $k=0$}\\ \frac{k}{2}x , &&& \text{if $k$ is even, $k \neq 0$} \\ -\left(a+\frac{k-1}{2} \right)x , &&& \text{if $k$ is odd} \end{aligned} \right. $$

[NOTE: The definition of $f_k$ was initially wrong, and has been fixed according to one of the answers.]

and

$$ g_k = a + k \, \text{,} $$

for $k \geq 0$.

I have tried to implement this continued fraction in Mathematica. This is what I have done:

f[k_, a_, x_] := If[EvenQ[k], (k/2)*x, -(a + (k - 1)/2)*x];
g[k_, a_] := a + k;
CONTFRAC[a_, x_] := 
  ContinuedFractionK[f[k, a, x] // Evaluate, 
   g[k, a] // Evaluate, {k, 0, +Infinity}];

Then I try

N[CONTFRAC[3/2, -100/3]]

(which should return something similar to $0.029542920686933995848485613587409381469661178405257$), but all I get is this:

This is what I get

So how can/should I ask Mathematica to calculate this continued sum?

NOTE: Why am I interested in doing that? Mathematica numerically calculates the value of the series with enough accuracy, but in other programming environments it is not like that for some values of $a$ and $x$ (for more details, see this math.SE thread).

$\endgroup$
5
$\begingroup$

There are a couple things wrong with your code. For one, ContinuedFractionK needs a cutoff maximal value of k to work. Second, you need to make sure that f only evaluates after a numerical value of k has been given. (Else, EvenQ[k] evaluates to False. Finally, f[0, a, x] evaluated to zero.

Here is some code that works:

f[k_Integer, a_, x_] := 
 If[EvenQ[k], If[k == 0, 1, (k/2)*x], -(a + (k - 1)/2)*x]; 
g[k_, a_] := a + k; 
CONTFRAC[a_, x_, kmax_] := 
 ContinuedFractionK[f[k, a, x], g[k, a], {k, 0, kmax}];

You now need to specify an kmax to evaluate the continued fraction, e.g.

N[CONTFRAC[3/2, -100/3, 100], 50]
(*0.029542920686933995848485613587409381469661178405257*)

You will have to vary kmax to check for convergence. (Either manually, or you have to write a loop.)

$\endgroup$
  • $\begingroup$ Thank you! I am going to change the definition of $f_k$ in the question. I thought that Mathematica was able to numerically compute infinite continued fractions (up to a given precision). $\endgroup$ – Vicent Oct 4 '17 at 15:08
3
$\begingroup$

First things first: we have the relationship

$$\sum_{k=0}^{\infty}\frac{z^k}{a(a+1)\cdots(a+k)}=\exp(z)z^{-a}\gamma(a,z)$$

where $\gamma(a,z)$ is the lower incomplete gamma function, which is expressed in Mathematica's notation as Gamma[a, 0, z].

Having gotten that out of the way, let me present a routine that uses the Lentz-Thompson-Barnett algorithm to evaluate a rearranged version of the CF in the OP:

SetAttributes[gamlow, Listable];
gamlow[a_?InexactNumberQ, z_?InexactNumberQ] := Module[{ak, c, d, di, f, h, k, tol, t2},
   tol = 10^(-Internal`PrecAccur[{a, z}]); t2 = tol^2;
   k = 1; f = c = 1; d = 0;
   While[ak = z If[OddQ[k], -(2 a + k - 1), k]/(2 (a + k - 1) (a + k));
         di = 1 + ak  d; If[di === 0, di = t2]; d = 1/di;
         c = 1 + ak/c; If[c === 0, c = t2];
         h = c d; f *= h; k++;
         Abs[h - 1] > tol];
   Exp[-z] z^a/(a f)]

(This routine can be compiled, if desired.)

Let's compare with the built-in:

{Plot3D[gamlow[a, z], {a, 1/2, 9/2}, {z, 0, 4}, PlotRange -> All], 
 Plot3D[Gamma[a, 0, z], {a, 1/2, 9/2}, {z, 0, 4}, PlotRange -> All]} // GraphicsRow

built-in versus CF

Plot3D[Abs[1 - gamlow[a, z]/Gamma[a, 0, z]], {a, 1/2, 9/2}, {z, 0, 4},
       PlotLabel -> "Relative Error", PlotRange -> All]

relative error

$\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.