# 1st-order linear ODE with tridiagonal matrix. Efficient solutions?

I have a 1st-rder linear ODE system where the system is characterized by $A$. Given an initial state $x_0$, I want the state at some later time $t$, efficiently.

$A$ happens to be a symmetric tridiagonal matrix where the coefficients are the same in each diagonal line.

$$A = \begin{pmatrix} a & b & 0 & 0 & ... \\ b & a & b & 0 & ... \\ 0 & b & a & b & ... \\ 0 & 0 & b & a & ... \\ \vdots &&\ddots&\ddots&\ddots\end{pmatrix}$$

I need to solve this for many different $A$'s and also many different $x_0$'s. The general solution to that is of course

$$x_t = e^{A t}x_0$$

One way to do this is to eigen-decompose $A$ into $A = F D F^\top$. Then the solution becomes

$$x_t = (F e^{D t} F^\top) x_0$$

Which is a little better because now the multiplication with $x_0$ is $O(n^2)$, but the eigen-decomposition is $O(n^3)$ still. But this doesn't take advantage of the constant values of $a$ and $b$ or the symmetric tridiagonal nature of $A$. So it seems like I should be able to do better.

-

The good news is that there are plenty of ways to numerically solve your ODE without having to resort to the matrix exponential, and these are all almost certainly more efficient. As an example, consider Euler's method. You want the solution at time $t = t_f$ to $x'(t) = Ax(t)$ with initial condition $x(0) = x_0$. Pick $K + 1$ equally-spaced time-points $0 = t_0 < t_1 < \cdots < t_K = t_f$, and let $h = (t_f - t_0)/K$ be the spacing between them. Let $x_k$ denote the approximation to the solution at time $t_k$. Euler's method is just the iteration $x_{k + 1} = x_k + hAx_k$.
Each step requires one matrix-vector multiply and one addition of two vectors. Since your matrix is tridiagonal, the matrix-vector multiply will be $O(3n)$, and the vector addition will be $O(n)$. With $K + 1$ steps, you're looking at an overall computational complexity of $O(nK)$. As far as accuracy goes, Euler's method converges at a rate of $O(h) = O(K^{-1})$, so if you double the number of points you use, you cut the error in half.
At any rate, there's no need to do an $O(n^3)$ eigenvalue computation if you don't want to!