The Fibonacci sequence is a driven by the second order linear difference equation

Fn+2 = Fn+1 + Fn, with boundary conditions F1 = 1, F2 = 1,

and thus can be solved exactly. As we know from practice that Fn is roughly exponential, we try Fn = Aa^n for A and a constants. This gives the quadratic

a^2 = a + 1, which happens to be the equation for the golden ratio Φ, and its inverse which I'll denote Φ'
(i.e. Φ' = 1/Φ, Φ' = Φ - 1)

As the equation is second order then it is a linear combination of these two solutions and the boundary conditions define the constants involved, i.e.

Fn = AΦ^n  + BΦ'^n

F0 = 0 (easy if you follow backwards) so A + B = 0
F1 = 1 . Using Φ = (1 + r)/2  and Φ' = (1 - r)/2  where r is the positive square root of 5, you can find

A - B = 2/r  yielding A = 1/r, B = -1/r

So Fn = (Φ^n /r) - (Φ'^n /r) = (Φ^n - Φ'^n)/r  for all n.

As can be seen, the even terms are when n is a multiple of 3, so using this formula add F3 + F6 + ... until you get a term greater than one million. Thus a program for this could be only a handful of lines long. A slightly further simplification would be to work out Φ^3 and Φ'^3, call them b and b' respectively. Then

F3k = (b^k - b'^k)/r  for k = 1,2,3...