<s>[INST] Finding a closed form for a recurrence relation $a_n=3a_{n-1}+4a_{n-2}$ Consider the sequence defined by
$$
\begin{cases}
a_0=1\\
a_1=2\\
a_n=3a_{n-1}+4a_{n-2} & \text{if }n\ge 2
\end{cases}
.$$
Find a closed form for $a_n$.

I tried listing out examples, but I don't see any common pattern between them.  All solutions are greatly appreciated.
 [/INST] This recurrence relation can be described with the following matrix:
$\operatorname A = \left[\begin{array}{cc}
3 & 4 \\ 
1 & 0 \\ 
\end{array} \right]$
This means that the equation, $$ \operatorname A\left[\begin{array}{cc}
a_{n-1} \\ 
a_{n-2} \\ 
\end{array} \right] = \left[\begin{array}{cc}
3a_{n-1} + 4 a_{n-2} \\ 
a_{n-1} \\ 
\end{array} \right] = \left[\begin{array}{cc}
a_{n} \\ 
a_{n-1} \\ 
\end{array} \right] $$
holds. Why is this beautiful? First of all, by repeated multiplication we have that $$\operatorname A^k\left[\begin{array}{cc}
a_{1} \\ 
a_{0} \\ 
\end{array} \right] = \left[\begin{array}{cc}
a_{k+1} \\ 
a_{k} \\ 
\end{array} \right]$$
which gives you a way to calculate the $ \mathbb N \ni n$th term in $O(\log n)$ time, using binary (matrix) exponentiation. Since the size of the matrix is really small, it is a very fast way indeed.
Okay, that's good and well, but you've asked for a closed form. Well, you can get that from this as well. Let's say that you've diagonalized this matrix, and you've gotten $ \operatorname A = S \Lambda S^{-1}$. In this case,  $$ \operatorname \Lambda = \left[\begin{array}{cc}
4 & 0 \\ 
0 & -1 \\ 
\end{array} \right]\quad \text{and} \quad S = \left[\begin{array}{cc}
4 & 1 \\ 
1 & -1 \\ 
\end{array} \right]. $$
Using $\operatorname A^k = S\Lambda^k S^{-1}$, we get that 
$$ \left[\begin{array}{cc}
a_{k+1} \\ 
a_{k} \\ 
\end{array} \right] = \operatorname A^k\left[\begin{array}{cc}
a_{1} \\ 
a_{0} \\ 
\end{array} \right]= S\Lambda^k S^{-1} \left[\begin{array}{cc}
a_{1} \\ 
a_{0} \\ 
\end{array} \right] = \left[\begin{array}{cc}
4^{k+1} & (-1)^k \\ 
4^k & (-1)^{k-1} \\ 
\end{array} \right] \left[\begin{array}{cc}
\frac{3}{5} \\ 
-\frac{2}{5} \\ 
\end{array} \right] =\left[\begin{array}{cc}
\frac{3}{5} 4^{k+1} + \frac{2}{5} (-1)^{k+1} \\ 
\frac{3}{5} 4^{k} + \frac{2}{5} (-1)^{k} \\ 
\end{array} \right].$$ which implies $$a_k = \frac{3}{5} 4^{k} + \frac{2}{5} (-1)^{k}.$$
Okay, okay this is all very well, but can't I do it faster? Well, in this case, you can.
Start by trying to find $a_n$ in the form of $q^n$, that is 
$$ a_n = 3a_{n-1} + 4a_{n-2} \quad \text{becomes}\quad  q^n = 3q^{n-1} + 4q^{n-2}$$
divide by $q^{n-2}$ and you'll get the quadratic equation $ q^2 - 3q - 4 = 0$. The roots of this polynomial are $-1$ and $4$. Assume that $a_n$ must be a linear combination of the $n$th power of these, that is $a_n$ has the form of $a_n = a4^n + b(-1)^n$.
To find $a$ and $b$, we can use that $a_0$ and $a_1$: substituting $n=0,1$ we get \begin{align} a + b &= a_0 = 1 \\  4a - b &= a_1 = 2 \end{align}
Solving this gives you $a = \frac{3}{5}, b= \frac{2}{5}$, which implies $$a_k = \frac{3}{5} 4^{k} + \frac{2}{5} (-1)^{k}.$$
If you are familiar with the meaning of the word characteristic polynomial, then you might notice that the "two" methods I've described above are closely connected, since $k_A(x) = x^2 - 3x - 4$ is the characteristic polynomial of $\operatorname A$. Moreover, calculating its roots is the same as calculating the eigenvalues of $\operatorname A$. The only big difference is that you don't have to know anything about matrices to do the second method, which is also faster, and the reason it works is hidden in the 1st method.
A word about calculating big powers of $\operatorname A$. This is more relevant if the matrix has a big size and a similarly nice characteristic polynomial. Let's say you wanted to calculate $A^{10^{10000}}$, and let's say that $A$ is an $m\times m$ matrix with a characteristic polynomial of $k_A$, where $\deg(k_A)$ is small, compared to $10^{10000}$.
By the Cayley–Hamilton theorem, we know that $k_A(A) = 0$.
So the beautiful trick is to calculate $h(x) = x^{10^{10000}} \mod k_A(x)$ first, and then $h(A)$. 
Essentially, given any polynomial $f(x)$, in order to calculate $f(A)$, one can find the corresponding element $r(x)$ in the ring $R[x]/k_A(x)$, and then calculate $r(A)$ using again, binary exponentiation.
</s>