Indeed (φ^n - ψ^n) / (φ - ψ) is what you'd get if you diagonalize the matrix. However the field you'd be working over is Q[√5], although really we're intereste

Indeed (φ^n - ψ^n) / (φ - ψ) is what you'd get if you diagonalize the matrix. However the field you'd be working over is Q[√5], although really we're interested in the subring Z[φ] here. And since φ is (trivially) an eigenvalue of the matrix, let's call it M, this is equivalent to calculating in Z[M] by the Cayley-Hamilton theorem. Of course some of these constructions might be easier to work with, and certainly something like Z[x]/(x^2-x-1) allows for fairly straightforward calculation, but you're still doing the same multiplications essentially. And if you're working in a convenient ring you don't really need the ψ^n term, since you can just use that φ^n = F_{n-1} + F_{n} φ.