Let $A$ be an $n\times n$ matrix, $B$ be an $n\times m$ matrix, $C$ an $m \times m$ matrix, and consider the sum $$\sum_{k = 0}^{N1} A^k B C^k.$$ Is there any smart way to rewrite this sum in a way similar to the partial sums of the geometric series; namely for $a,b \in \mathbb{R},$ $$\sum_{k = 0}^{N1} b a^k = b \frac{1a^N}{1a}?$$

1$\begingroup$ take the simple case $B=I$, $n=3$, what "simplication" of $I+AC+A^2C^2$ would you have in mind? $\endgroup$– Carlo BeenakkerAug 18, 2022 at 14:17

$\begingroup$ No, but you can solve the equation $\sum_{k = 0}^{N1} A^k X C^k = E$ on a computer, for given $A, C, E$, in $O(m^3+n^3)$ operations. Would you be interested in that? $\endgroup$– Federico PoloniAug 18, 2022 at 18:20
2 Answers
This is a long comment. Write $S$ for this sum. If we just directly imitate the usual geometric series argument we are led to consider
$$ASC = \sum_{k=1}^N A^k B C^k = S  B + A^N B C^N$$
so $ASC  S = A^N B C^N  B$. And here we're a bit stuck: $ASC  S$ just isn't a multiplication of $S$ by any matrix on either the left or the right so we can't hope to perform a simple matrix division here. We just need to talk about the inverse, if it exists, of the linear map $S \mapsto ASC  S$. But at least this is a system of linear equations whose solution, which we can hope is unique, gives the coefficients of $S$ so that's not too bad.

$\begingroup$ Thanks. Indeed I tried a similar idea but I don't need a specific solution, I need something some regularity in the dependency on $N$ because then afterwards I have to plug this into other series. $\endgroup$– tomateAug 19, 2022 at 18:45
With meager loss of generality, let's first assume that $A,C$ are diagonalizable so that $A=RDR^{1},C=SES^{1}$ where $R,S$ are invertible and $D,E$ are diagonal.
Then $$\sum_{k=0}^{N1}A^kBC^k=\sum_{k=0}^{N1}RD^kR^{1}BSE^kS^{1}$$
$$=R(\sum_{k=0}^{N1}D^kR^{1}BSE^k)S^{1}.$$
Suppose now that $R^{1}BS=(\alpha_{i,j})_{i,j}$ and that the diagonal entries in $D$ are $(\delta_1,\dots,\delta_n)$ and the diagonal entries in $E$ are $(\epsilon_1,\dots,\epsilon_m)$.
Then $$\sum_{k=0}^{N1}A^kBC^k=R(\sum_{k=0}^{N1}\delta_i^k\alpha_{i,j}\epsilon_j^k)_{i,j}S^{1}$$
$$=R\cdot\Big(\frac{(\delta_i\epsilon_j)^N1}{\delta_i\epsilon_j1}\cdot\alpha_{i,j}\Big)_{i,j}\cdot S^{1}.$$
It looks like there is a not quite as nice expression for our sum when $A,C$ are not diagonalizable, and $D,E$ are in Jordan normal form.
Generalization:
Suppose that $f(x,y)$ is a polynomial (or an analytic function when we assume convergence), and $f(x,y)=\sum_{k,l}b_{k,l}x^ky^l$. Then by using the same reasoning and matrices as before, we have $$\sum_{k,l}b_{k,l}A^kBC^l=R\cdot(f(\delta_i,\epsilon_j)\cdot\alpha_{i,j})_{i,j}\cdot S^{1}.$$

$\begingroup$ Indeed $A$ and $C$ are diagonalizable. I hoped to be able to find a basisindependent formula, but this sounds like a very interesting way to proceed. $\endgroup$– tomateAug 19, 2022 at 18:51