4
$\begingroup$

It takes a long time to compute the summation below, and I'd like to know if there are alternative ways
to compute things faster. When replacing $15$ by $\infty$, then I should get $3^{1/3}$. I need Mathematica solution
not mathematical solutions (since here I can handle things on my own).

Clear[b, k]
b[0] := 1;
b[k_] := b[k - 1] + 1/(b[k - 1]^2 + b[k - 1] + 1)
N[1/(15)^(1/3) Sum[1/(b[k]^2 + b[k] + 1), {k, 0, 15}]]
$\endgroup$
  • $\begingroup$ I'd like to see what happens when I set $1000$ instead of $15$. $\endgroup$ – user 1357113 Oct 6 '13 at 18:57
4
$\begingroup$

EDIT

Clear[b, k]
b[0] := 1;
b[k_] := ( b[k] = N[b[k - 1] + 1/(b[k - 1]^2 + b[k - 1] + 1), 100] )
1/(1000)^(1/3) Sum[1/(b[k]^2 + b[k] + 1), {k, 0, 1000}]

(*1.28966188938645818444514220128278717806515296115261423714827563492762\
  4979228386408865750383532338923*)

very close to 3^(1/3) = 1.4422...

The convergence is slow. So you should do extrapolation.

$\endgroup$
  • $\begingroup$ Thank you for your answer (+1). It seems that $N$ in front caused some trouble. $\endgroup$ – user 1357113 Oct 6 '13 at 19:19
  • $\begingroup$ yeah, the limit is precisely $3^{1/3}$. $\endgroup$ – user 1357113 Oct 6 '13 at 19:32
  • 1
    $\begingroup$ I actually put N in the definition of b[k], since its exact values are huge rational numbers. $\endgroup$ – Vahagn Poghosyan Oct 6 '13 at 19:32
  • $\begingroup$ Note that here k is local variable, so it is not necessary to Clear it at the beginning of the code. $\endgroup$ – Vahagn Poghosyan Oct 6 '13 at 19:39
  • $\begingroup$ I got the point. Thanks! $\endgroup$ – user 1357113 Oct 6 '13 at 19:46
5
$\begingroup$
b[0] := 1;
b[k_] := b[k] = N[b[k - 1] + 1/(b[k - 1]^2 + b[k - 1] + 1)]
f[x_] := f[x] = f[x - 1] + (1/(b[x]^2 + b[x] + 1))
f[0] := 0

Manipulate[
 ListPlot[Table[1/x^(1/3) f[x], {x, 1, a}], MaxPlotPoints -> 100, 
  Epilog -> {Green, Line[{{1, 3^(1/3)}, {a, 3^(1/3)}}]}, 
  PlotRange -> {0, 4}], {a, 10, 100000, 1}]

Mathematica graphics

$\endgroup$
3
$\begingroup$

Here is an alternative,

Nest[# + 1/(1 + # + #^2) &, 1., n + 1] - 1

that will compute the sum over 0 <= k <= n. For 100-digit precision use

Nest[# + 1/(1 + # + #^2) &, 1.`100, n + 1] - 1

Examples

With[{n = 1000}, 
  1/n^(1/3) (Nest[# + 1/(1 + # + #^2) &, 1.`, n + 1] - 
     1)] // AbsoluteTiming
(* {0.000359, 1.28966} *)

With[{n = 1000}, 
  1/n^(1/3) (Nest[# + 1/(1 + # + #^2) &, 1.`100, n + 1] - 1)] // AbsoluteTiming
(* {0.005293, 
    1.28966188938645818444514220128278717806515296115261423714827563492762\
      4979228386408865750383532338923} *)

With[{n = 1000000}, 
  1/n^(1/3) (Nest[# + 1/(1 + # + #^2) &, 1.`100, n + 1] - 1)] // AbsoluteTiming
(* {5.204264,
    1.42720128317421809656439478563254431280673078216998116597726119064353\
      220814168138689903801038204268} *)

It avoids the recursion and memory limits that the OP's b runs into. (With n = 1000000, my kernel crashed using b and Sum.) And you can use NestList to get all the sums up to n.

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