Recurrence equation for sequence of vectors Consider recurrent formula for a sequence of numbers $(y_n)$ (either real or complex):
$$a_k y_{n+k}+a_{k-1}y_{n+k-1}+\cdots+a_0y_n=\sum_{i=0}^k a_i y_{n + i} = 0$$
It's known that the exact explicit formula for such sequence is expressed via the roots of the characteristic polynomial:
$$P_k(\lambda)=\sum_{i=0}^k a_i \lambda^i = 0$$
In case of $k$ distinct roots $\lambda_1, \lambda_2,\ldots, \lambda_k$ the expression is simply:
$$y_n=\sum_{i=0}^k c_i \lambda_i^n = 0$$
Where value of coefficients $c_i$ depend on the initial conditions. In other words we need to know $k$ first values of the sequence $y_1, y_2,\ldots, y_k$ for it to be uniquely defined.
Now, I am trying to generalize such approach on the case of vector sequences $v_m \in \mathbb{R}^m, m > 1$:
$$A_k v_{n+k}+A_{k-1}v_{n+k-1}+\cdots+A_0v_n=\sum_{i=0}^k A_i v_{n + i} = 0$$
In this case $A_i$ is a square matrix $(A_i \in \operatorname{Mat}(m \times m))$.
What I've got so far is that in the simplest case:
$$v_{n+1}=A_0v_n$$
the solution can be expressed via eigenvalues of the matrix $A_0$. Again, we have $m$ not multiple distinct eigenvalues:
$$v_n=\sum_{j=1}^mc_j\lambda_j^nh_j$$
where $(\lambda_j, h_j)$ is eigenpair: $A_0h_j=\lambda_jh_j$.
Is it possible to generalize such approach? Maybe someone has suggestions or can link to the appropriate literature? Any help will be appreciated.
 A: Yes, it is possible to generalize this. At least after a suitable reformulation.
I assume that $a_k=1$, and that we solved $y_{n+k}$ the recurrence relation:
$$
y_{n+k}=-a_{k-1}y_{n+k-1}-\cdots-a_0v_n.
$$
Let's begin with the sequence of numbers. We can think of a segment of the sequence $\vec{y}^{(i)}=(y_i,y_{i+1},\ldots,y_{i+k-1})^T$ as a column vector describing the state of this system at the moment of time $t=i$. Then the evolution of the system is described by the matrix product
$$
\vec{y}^{(i+1)}=P \vec{y}^{(i)},
$$
where the matrix $P$ is the companion matrix of the characteristic polynomial $P_k(\lambda)$
$$
P=\left(\begin{array}{cccccc}
0&1&0&0&\cdots&0\\
0&0&1&0&\cdots&0\\
0&0&0&1&\cdots&0\\
\vdots&\vdots&\vdots&\ddots&\ddots&\vdots\\
0&0&\cdots&\cdots&0&1\\
-a_0&-a_1&\cdots&\cdots&\cdots&-a_{k-1}\\
\end{array}\right).
$$
It is east to see that the eigenvalues of $P$ are the zeros of the characteristic polynomial.  Let $\vec{x}_i, i=1,2,\ldots,k,$ of $P$, be the eigenvectors of $P\vec{x}_i=\lambda_i\vec{x_i}$. 
By writing the initial vector $\vec{y}^{(0)}$ as a linear combination of the eigenbasis
$$
\vec{y}^{(0)}=\sum_{i=1}^kc_i\vec{x}_i,
$$
we then see that the future is completely determined as
$$
\vec{y}^{(n)}=P^n\vec{y}^{(0)}=\sum_{i=1}^kc_i\lambda_i^n\vec{x}_i
$$
for all natural numbers $n$. Alternatively we can diagonalize $P$ and the find a formula for its power $P^n$.
The above works, when there is an eigenbasis of $P$. Should there be roots of $P_k(\lambda)$ with multiplicities $>1$ and not enough eigenvectors, we can no longer diagonalize $P$, but we can still put in Jordan canonical form (over $\Bbb{C}$ at least), and use that in the same way. This has the well known effect of introducing polynomial factors into the solutions.
On with your question. We can still form a "vector" describing the state of your sequence. Only this time the state is described by a vector in $\Bbb{R}^{mk}$ (or $\Bbb{C}^{mk}$). Say,
$$
\vec{y}^{(i)}=(\vec{v}_i,\vec{v}_{i+1},\ldots,\vec{v}_{i+k})^T
$$
with $k$ components. The time evolution is still described by a matrix equation
(I assume that $A_k$ is invertible and again solve for $\vec{v}_{n+k}$)
$$
\vec{y}^{(i+1)}=P\vec{y}^{(i)},
$$
but this time $P$ is an $mk\times mk$ matrix with a block structure
$$
P=\left(\begin{array}{cccccc}
0&I_m&0&0&\cdots&0\\
0&0&I_m&0&\cdots&0\\
0&0&0&I_m&\cdots&0\\
\vdots&\vdots&\vdots&\ddots&\ddots&\vdots\\
0&0&\cdots&\cdots&0&I_m\\
-A_0&-A_1&\cdots&\cdots&\cdots&-A_{k-1}\\
\end{array}\right)
$$
with $m\times m$ blocks according to the familiar scheme.
To solve such a recurrence relation you again want to find an eigenbasis of $P$
(or, failing that, a basis that puts it into a Jordan canonical form). The same techniques should apply.
