34
$\begingroup$

Today, I found the Euler Project. Problem #2 is

Each new term in the Fibonacci sequence is generated by adding the previous two terms. By starting with 1 and 2, the first 10 terms will be:

 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, ...

By considering the terms in the Fibonacci sequence whose values do not exceed four million, find the sum of the even-valued terms.

I have two solutions:

For[i = 1, i < 100, i++, If[Fibonacci[i] > 4000000, Break[]]];
Plus @@ Select[Fibonacci /@ Range[i - 1], Mod[#, 2] == 0 &]

4613732

Plus @@ Select[
Fibonacci /@ Select[Range@100, Fibonacci@# < 4000000 &], 
Mod[#1, 2] == 0 &]

4613732

However, I feel my solutions are not good and efficient. So my question is: can you show me more concise and efficient methods?

Also, how cam I test whether a number is a Fibonacci number?

$\endgroup$
  • 3
    $\begingroup$ Since you already have the answer, you may access the problem 2 forum on the Euler project. Here you will find the most efficient algorithms from other users. The third post there demonstrates an algorithm that can be easily calculated without a program. $\endgroup$ – travisbartley Nov 19 '13 at 5:38
45
$\begingroup$

The other methods described work well, but I am partial to exact answers. So, here it is. First, note that the Fibbonacci sequence has a closed form

$$F_n = \frac{\Phi^n - (-\Phi)^{-n}}{\sqrt{5}}$$

where $\Phi$ is the GoldenRatio. Also, as observed elsewhere, we are looking only at every third term, which gives the sum

$$\sum^k_{n = 1} \frac{\Phi^{3n} - (-\Phi)^{-3n}}{\sqrt{5}}$$

which is just the sum of two geometric series. So, after a few steps, this simplifies to

$$\frac{1}{\sqrt{5}} \left[ \Phi^3\frac{1 - (\Phi^3)^{k}}{1-\Phi^3} + \Phi^{-3}\frac{1 - (-\Phi^{-3})^{k}}{1 + \Phi^{-3}} \right]$$

where $k$ is the index of the highest even Fibbonacci number. To find $k$, we can reverse the sequence,

n[F_] := Floor[Log[F Sqrt[5]]/Log[GoldenRatio] + 1/2]
n[4000000]
(* 33 *)

So, $k = 11$. In code,

Clear[evenSum]; 
evenSum[(k_Integer)?Positive] := 
  Round @ N[
    With[
      {phi = GoldenRatio^3, psi = (-GoldenRatio)^(-3)}, 
      (1/Sqrt[5])*(phi*((1 - phi^k)/(1 - phi)) - psi*((1 - psi^k)/(1 - psi)))
    ],
    Max[$MachinePrecision, 3 k] (* needed for accuracy *)
]

evenSum[11]
(* 4613732 *)

which is confirmed by

Total @ Fibonacci @ Range[3, 33, 3]
(* 4613732 *)

Edit: The above implementation of n suffers round off problems for large numbers, for example n[9000] gives 21, but Fibonacci[21] is 10946. But, wikipedia gives a better choice

Clear[n2];
n2[fn_Integer?Positive] := 
 Floor@N@Log[GoldenRatio, (fn Sqrt[5] + Sqrt[5 fn^2 - 4])/2]

n2[10945]
(* 20 *)
$\endgroup$
  • $\begingroup$ Had to fix the formulas, I dropped quite a few minus signs in the original version. $\endgroup$ – rcollyer Nov 18 '13 at 18:04
  • $\begingroup$ Related to your error, of course facebook.com/physic.department/posts/360869444018997 $\endgroup$ – Dr. belisarius Nov 18 '13 at 18:25
  • $\begingroup$ @belisarius :)) ingenious. $\endgroup$ – Stefan Nov 18 '13 at 18:27
  • $\begingroup$ @belisarius according to my adviser, his adviser was infamous for making multiple canceling errors. His results were correct, the derivation not usually. :) $\endgroup$ – rcollyer Nov 18 '13 at 18:31
  • $\begingroup$ Could you explain how these functions n and evenSum should be used for big numbers? I tried to test your solution unsucessfully. While with mine approach I could simply find sums up to 10^10000. $\endgroup$ – Artes Nov 18 '13 at 18:36
20
$\begingroup$

Straight iteration over the even valued Fibonacci numbers is fast.

fibSum[max_] := Module[
  {tot, n, j},
  tot = 0;
  n = 0;
  j = 3;
  While[n = Fibonacci[j]; n <= max, j += 3; tot += n];
  tot]

Or one can use matrix products. This seems to be about the same speed. It has the advantage of not requiring a built in Fibonacci function.

fibSum2[max_] := Module[
  {tot, n, fp, mat, mat3},
  mat = {{0, 1}, {1, 1}};
  tot = 0;
  n = 0;
  fp = mat.mat;
  mat3 = mat.mat.mat;
  While[n = fp[[2, 2]]; n <= max, fp = fp.mat3; tot += n];
  tot]

These handle 10^10000 in slightly under a second on my desktop running Mathematica 9.0.1.

Here is a test for whether an integer is a Fibonacci number (a Fib detector?)

fibQ[n_] := 
 With[{k = 
    Round[Log[GoldenRatio, N[n, 20 + Length[IntegerDigits[n]]]] + 
      Log[GoldenRatio, Sqrt[5]]]},
  n == Fibonacci[k]]
$\endgroup$
20
$\begingroup$
Sum[Fibonacci[n], 
 {n, 3, InverseFunction[Fibonacci][4000000], 3}]

(* 4613732 *)
$\endgroup$
  • 1
    $\begingroup$ That's tyte, + goes w/o saying. $\endgroup$ – alancalvitti Oct 30 '14 at 18:18
14
$\begingroup$

Here I'd like to show my solution.

Since A014445 we know that the generating function of even fibonacci numbers is Fibonacci[3 n], we can write something like this:

Plus @@ Select[Table[Fibonacci[3 n], {n, 0, 30}], # < 4000000 &] // AbsoluteTiming

=> {0., 4613732}

An imperative style solution could look like this:

res = 0; n = 1; While[res < 4000000, res += Fibonacci[3 n]; n++]

Result will be in res and timings are negligible, zero.

Edit: To bring that into a cogent form. Here the full definition of SumEvenFibonacci:

ClearAll[SumEvenFibonacci];
SumEvenFibonacci::usage = 
  "SumEvenFibonacci[n] calculates the sum of even Fibonacci numbers \
   up to the upper-bound n.";
SumEvenFibonacci[n_Integer] /; n >= 0 := Module[{res = 0, i = 1},
    While[res < n, res += Fibonacci[3 i]; i++];
    res
]
$\endgroup$
14
$\begingroup$

The Fibonacci function, is Listable and fast increasing. We can observe a simple pattern, every third element is even (see e.g. Fibonacci[Range[12]]), moreover the first number exceeding 4000000 is: NestWhile[(# + 1) &, 1, Fibonacci[#] <= 4 10^6 &], thus we have:

Total[ Fibonacci[ Range[3, 33, 3]]]
4613732

and most likely this is the best approach:

sum[n_] := Total[ Fibonacci[ Range[3, NestWhile[(# + 1) &, 1, Fibonacci[#] <= n &] - 1, 3]]]

e.g.

sum[10^600]; // AbsoluteTiming
{0.096000, Null}

Concerning testing if given numbers are Fibonacci ones there is a simple test based on checking if consecutive elements generate Pythagorean triples (see e.g. this answer)

And @@ (#1^2 + #2^2 == #3^2 & @@@ 
  Array[{  Fibonacci[#] Fibonacci[# + 3], 
         2 Fibonacci[# + 1] Fibonacci[# + 2], 
           Fibonacci[# + 1]^2 + Fibonacci[# + 2]^2} &, 1000])
True
$\endgroup$
  • $\begingroup$ You already got my +1, but, yeah, it looks a little cretinous. $\endgroup$ – rcollyer Nov 18 '13 at 21:14
  • $\begingroup$ @rcollyer I like your approach, but I still cannot use it successfully. I guess Daniel Lichtblau's solution is the most efficient for big numbers, however I can see that another reliable (mine) solution is the least appreciated. $\endgroup$ – Artes Nov 18 '13 at 21:20
  • $\begingroup$ A little poking around indicates that n gives the nearest Fibbonacci number, not the one immediately below it, e.g. n[9000] returns 21, not 20. Likely a similar rounding issue is occurring in evenSum, but it tracks Daniel's exactly for a long time. So, I'm not sure of the specific issue. $\endgroup$ – rcollyer Nov 18 '13 at 21:30
12
$\begingroup$

Your solution looks fine. You can take advantage of the fact that every third Fibonacci number is even, which makes it a little faster. Fibonaccis are cheap to compute and they quickly exceed 4 million. Here's a comparison:

Select[Fibonacci[3 Range@33], # <= 4*^6 &] // Total // AbsoluteTiming

{0.000214, 4613732}

Select[Fibonacci@Range@100, # <= 4*^6 && EvenQ@# &] // Total // AbsoluteTiming

{0.000566, 4613732}

Plus @@ Select[Fibonacci /@ Select[Range@100, Fibonacci@# < 4000000 &], Mod[#1, 2] == 0 &] // AbsoluteTiming

{0.000717, 4613732}

AbsoluteTiming[ For[i = 1, i < 100, i++, If[Fibonacci[i] > 4000000, Break[]]]; Plus @@ Select[Fibonacci /@ Range[i - 1], Mod[#, 2] == 0 &]]

{0.000526, 4613732}

$\endgroup$
4
$\begingroup$

I will make use of the formula:

sum[x_] := (Fibonacci[3 Quotient[Floor[InverseFunction[Fibonacci][x]], 3]+2]-1)/2

For $x = 4 \times 10^6$ this gives:

sum[4 × 10^6] //AbsoluteTiming

{7.82619, 4613732}

This is pretty slow. We can speed up the computation of the inverse of the Fibonacci using the InverseFibonacci function from my answer to (157354). Below I give a version of that inverse function that works for arguments greater than 100:

if[x_?(GreaterThan[100])]:=Root[
    {
    Fibonacci[#]- x&, SetPrecision[Log[GoldenRatio, x Sqrt[5]],
    Min[2 + Log10[x], 10]]
    }
]

Now, using if above instead of InverseFunction[Fibonacci]:

fast[x_] := (Fibonacci[3 Quotient[Floor @ if[x], 3] + 2] - 1)/2

and applying this to the problem:

fast[4 10^6] //AbsoluteTiming

{0.000239, 4613732}

We can use this function on much larger arguments as well:

fast[
    92837410293847120384712309847213048971230498127340912874091273401298374012938741029
] //AbsoluteTiming

{0.002153, 88011840322506983234113472696205625385192191652246095943362996448287672522108009837}

$\endgroup$
1
$\begingroup$

Actually this sum can be closed-form analytically:

(* sum of even-valued Fibonacci terms that do not exceed z *)
sum[z_] := Module[{g, n, f},
  g = (1 + Sqrt[5])/2;
  n = Floor[Log[Sqrt[5] z]/Log[g]];
  f = Round[g^(n + 2)/Sqrt[5]];
  (f - 1)/2]

sum[4*^6] // RepeatedTiming    
(* {0.0000448, 4613732} *)

sum[1*^100] // RepeatedTiming
(* {0.000045, 12065...3520} *)

I used these identities on mathworld: (8), (22), although one should first recognize that only every third term is even.

$\endgroup$

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