Trying to understand a solution to Project Euler #2 I'm trying to understand the Haskell solution given here (The third of the four solutions given) to Problem 2 to project Euler. I think I've interpreted the code right but I do not think I understand the math involved. I've tried to reduce it to Mathematical terms as much as possible. $(a,b)[1]:=a, (a,b)[2]:=b$
Define $(a,b)\times (c,d)=(ad+bc-4ac,ac+bd)$
Consider the Fibonacci series starting with $1,1$ so that $F_{3n}$ is always even.
Fix arbitary natural $k$
Let $l=F_{3k}/2$
Let $m$ be the maximum natural such that $(1,4)^{2^m}[1]\le l$
Consider the sequence $A_i$
$A_1 = (0,1)$
$A_{i+1} = f((1,4)^{2^{m+1-i}},A_i)$
where $f(x,y)=x\times y$ if $(x\times y)[0] \le l$, else $y$
Then show that $A_m[1]+A_m[2]=F_{3k+2}$
I know the fibonacci Sum formula ($f_1+f_2+\cdots +f_n=f_{n+2}-1$) and the recurrence relation of the even only case $(f_{n+6}=4f_{n+3}+f_n)$
I'm not sure if I'm using the correct tags, feel free to change if required.
 A: For reference, here is the code you are referencing:
sumEvenFibsLessThan n = (a + b - 1) `div` 2
  where
    n2 = n `div` 2
    (a, b) = foldr f (0,1)
             . takeWhile ((<= n2) . fst)
             . iterate times2E $ (1, 4)
    f x y | fst z <= n2 = z
          | otherwise   = y
      where z = x `addE` y
addE (a, b) (c, d) = (a*d + b*c - 4*ac, ac + b*d)
  where ac=a*c

times2E (a, b) = addE (a, b) (a, b)

Let $H_n=(1/2)F_{3n}$ be half of the $n^{th}$ even Fibonacci number. Using $F_n=F_{n-1}+F_{n-1}$, you can show
$$
H_{n+1}=4H_n+H_{n-1}\tag1
$$
This lets you get a simple formula for the sum of the first $n$ even Fibonacci numbers:
\begin{align}
2(H_1+H_2+\dots+H_n)
&=2\big(\frac14(H_2-H_0)+\frac14(H_3-H_1)+\dots+\frac14(H_{n+1}-H_{n-1})\big)
\\&=\frac12\big(H_{n+1}+H_n-H_1-H_0\big)
\\&=\frac12\big(H_{n+1}+H_n-1\big)
\end{align}
Therefore, once you have the ordered pair $(a,b)=(H_n,H_{n+1})$, the answer you want is $(a+b-1)/2$ (notice that formula in the first line of the code).
How does the program compute $(H_{n},H_{n+1})$? They build this vector up using versions of that vector for smaller $n$. Specifically, you can prove the following Fibonacci identities:
\begin{align}
H_{n+m+1}&=H_{n+1}H_{m+1}+H_nH_m,\\\tag{*}
H_{n+m}&=H_{n+1}H_m+H_nH_{m+1}-4H_nH_m
\end{align}
Therefore, the vector $(H_{n+m},H_{n+m+1})$ can be recovered from the vectors $(H_n,H_{n+1})$ and $(H_m, H_{m+1})$ using the funny multiplication operation defined in that program:
$$
(H_{n+m},H_{n+m+1})=(H_n,H_{n+1})\star (H_m, H_{m+1}),\text{ where}\\
(a,b)\star (c,d) = (ad+bc-4ac,ac+bd)
$$
In particular,
$$
(H_{n},H_{n+1})=(H_1,H_2)\star (H_{n-1},H_n)=(1,4)\star (H_{n-1},H_n),
$$
which implies
$$
(H_n,H_{n+1})=(1,4)^n \star (H_0,H_{-1})=(1,4)^n \star(0,1)
$$
The program computes $(1,4)^n$ using exponentiation by squaring; this is what is happening with the sequence $A_i$, and with the iterate/takeWhile part of the program above. 

To prove $(*)$, note that $(1)$ implies the following matrix equation:
$$
\begin{bmatrix}4&1\\1&0\end{bmatrix}\begin{bmatrix}H_{n}\\H_{n-1}\end{bmatrix}
=
\begin{bmatrix}H_{n+1}\\H_{n}\end{bmatrix}\tag2
$$
Iterating this equation, and examining the base cases carefully, leads to the realization that
$$
\begin{bmatrix}4&1\\1&0\end{bmatrix}^n=\begin{bmatrix}H_{n+1}&H_n\\H_n&H_{n-1}\end{bmatrix}\tag3
$$
By examining the components of both matrices in the equation
$$
\begin{bmatrix}4&1\\1&0\end{bmatrix}^n\cdot \begin{bmatrix}4&1\\1&0\end{bmatrix}^m=\begin{bmatrix}4&1\\1&0\end{bmatrix}^{n+m},\tag4
$$
you can see that $(3)$ and $(4)$ together imply
\begin{align}
H_{n+m+1}&=H_{n+1}H_{m+1}+H_nH_m,\\
H_{n+m}&=H_{n+1}H_m+H_nH_{m-1}\\
&=H_{n+1}H_m+H_nH_{m+1}-4H_nH_m
\end{align}
