Sums of terms from Fibonacci

Find the sum of all the even-valued terms in the Fibonacci sequence which do not exceed four million.

I would like to use this puzzle as an excuse to explore the Fibonacci sequence. For instance, we know from the definition that each term is the sum of the previous two: $x_n=x_{n-1}+x_{n-2},$ and the first two terms of the sequence are given by $x_1=0$ and $x_2=1$. A simple python function that computes the n-th term of this sequence can be then written (in a very naïve way) as follows:

def fib(n):
counter,a,b=0,0,1
if n==0: return 0
elif n==1: return 1
else:
while counter&lt;n:
a,b,counter=b,a+b,counter+1
return a


And a python script using this function that performs the task at hand could be something like this:

contador,sum=0,0


[3, 2]
[6, 8]
[9, 34]
[12, 144]
[15, 610]
[18, 2584]
[21, 10946]
[24, 46368]
[27, 196418]
[30, 832040]
[33, 3524578]
>>> sum
4613732

A much faster way to compute the same values can be done by exploiting some ideas from linear algebra: Notice that each two successive terms of the sequence satisfy

$\begin{pmatrix} x_{n}\\ x_{n-1} \end{pmatrix} = \begin{pmatrix} x_{n-1}+x_{n-2}\\x_{n-1} \end{pmatrix} = \begin{pmatrix} 1 & 1\\ 1& 0\end{pmatrix}\, \begin{pmatrix} x_{n-1}\\x_{n-2} \end{pmatrix}$

And tracing back to the first two terms of the sequence ($x_1=0$, $x_2=1$), we have

$\begin{pmatrix} x_n\\x_{n-1} \end{pmatrix} = \begin{pmatrix} 1&1\\1&0\end{pmatrix}^{n-1} \begin{pmatrix}1\\ 0\end{pmatrix}$

We have turned the code into linear (meaning that we need to do as many computations as a multiple of the input value $n$), as opposed as the previous code, in which we made a multiple of $x_n$ computations to obtain the n-th term of the sequence, $x_n.$   This is a great advantage, but we still need to compute the powers of the matrix $\big( \begin{smallmatrix}1&1\\1&0\end{smallmatrix} \big)$ each time. Again, to aid in this task we employ some more machinery from linear algebra:

The matrix above is diagonalizable, which means that we can find an invertible matrix $P,$ and two real numbers $\lambda_1, \lambda_2 \in \mathbb{R}$ such that

$\begin{pmatrix} 1&1\\1&0\end{pmatrix} = P \begin{pmatrix}\lambda_1&0\\ 0&\lambda_2\end{pmatrix} P^{-1}$

In order to compute the matrix $P$ and the eigenvalues $\lambda_1, \lambda_2$, we use the characteristic polynomial of the matrix:

$\det \begin{pmatrix}1-\lambda&1\\1&-\lambda \end{pmatrix} = -(1-\lambda)\lambda-1 = \lambda^2-\lambda-1,$

which has two different roots $\lambda \in \{\frac{1}{2}(1\pm \sqrt{5})\},$ The corresponding eigenvectors give us the entries of the matrix $P$: I will leave to the reader as an exercise the verification that, if we denote $\lambda_1=\tfrac{1}{2}\big(1+\sqrt{5}\big)$ and $\lambda_2=\tfrac{1}{2}\big(1-\sqrt{5}\big),$ then $P=\big(\begin{smallmatrix}\lambda_1&\lambda_2\\1&1\end{smallmatrix}\big).$ Being this the case, we have

$\begin{pmatrix}1&1\\1&0\end{pmatrix}^{n-1} = \begin{pmatrix} \lambda_1&\lambda_2\\1&1\end{pmatrix}\, \begin{pmatrix}\lambda_1^{n-1}&0\\ 0&\lambda_2^{n-1}\end{pmatrix}\, \begin{pmatrix} \lambda_1&\lambda_2\\1&1\end{pmatrix}^{-1}.$

It only remains to put all this information together to be able to write a compact code that computes the n-th term of the Fibonacci sequence in constant time (that is, it takes exactly the same time for any value of $n$):

$\begin{pmatrix}x_n\\x_{n-1}\end{pmatrix}=\begin{pmatrix}\lambda_1&\lambda_2\\1&1\end{pmatrix}\,\begin{pmatrix}\lambda_1^{n-1}&0\\ 0&\lambda_2^{n-1}\end{pmatrix}\,\begin{pmatrix}\lambda_1&\lambda_2\\1&1\end{pmatrix}^{-1}\,\begin{pmatrix}1\\ 0\end{pmatrix}.$

We can then write an explicit formula for each term of the Fibonacci sequence:

$x_n = \dfrac{1}{\sqrt{5}}\, \bigg( \dfrac{1+\sqrt{5}}{2} \bigg)^{n}-\dfrac{1}{\sqrt{5}}\, \bigg( \dfrac{1-\sqrt{5}}{2} \bigg)^{n}$

We are now ready to tackle the puzzle with more sophisticated tools. Notice first that the parity of the Fibonacci numbers has the periodic pattern even, odd, odd, even, odd, odd, even, odd, odd, This helps us identify the even-valued terms of the Fibonacci sequence as $x_{3n}$ for all $n \geq 0.$ The answer to our puzzle can then be written as follows: Let $M$ such that $x_M \leq 4\cdot 10^6 \leq x_{M+1}$, then the sum of all even-valued terms of the Fibonacci sequence up to $x_M$ is given by

$\begin{array}{rcl} \displaystyle{\sum_{n=1}^{M/3}} x_{3n} &=& \dfrac{1}{\sqrt{5}}\, \displaystyle{\sum_{n=1}^{M/3}} \lambda_1^{3n} - \dfrac{1}{\sqrt{5}}\,\displaystyle{\sum_{n=1}^{M/3}} \lambda_2^{3n}\\ \\&=&\dfrac{1}{\sqrt{5}}\,\dfrac{\lambda_1^3 \big( 1-\lambda_1^M\big)}{1-\lambda_1^3} - \dfrac{1}{\sqrt{5}}\,\dfrac{\lambda_2^3\big(1-\lambda_2^M\big)}{1-\lambda_2^3}\end{array}$

A nice compact formula: The only missing piece of data is the value of $M.$   How would the reader find its value?