I have the following system
$$\begin{bmatrix} \dot{u}_1 \\ \dot{u}_2 \\ \vdots \\ \dot{u}_{999} \end{bmatrix} = -\frac{1}{1000^2}\operatorname{tridiag}(-1,2,-1)\begin{bmatrix} u_1 \\ u_2 \\ \vdots \\ u_{999} \end{bmatrix}$$
where the tridiagonal matrix (call it $A$) has dimensions $999\times999$ and $u_1(0) = \dots u_{999}(0) = 1.$
I'm asked what time step I would choose when using the explicit Euler method. My initial thoughts were to consider stability and using the explicit euler method I obtain $$u_{n+1} = (I -\frac{h}{1000^2}A)u_n.$$ I know from before that the eigenvalues of $A$ is $$-4\sin^2\left(\frac{\pi m}{2(M+1)}\right), \quad m=1,2,\dots,M$$ So I obtain the eigenvalues $$\lambda_m = 1 + \frac{4h}{1000^2}\sin^2\left(\frac{\pi m}{2(M+1)}\right)$$ which should be less than $1$ for stability which means but this is not possible for nonzero $h$. So what am I doing wrong?
Consider the system of ordinary differential equations given by $$y'(t) = By(t)$$ where $B$ is a square matrix and $u$ is a compatible vector. Euler's method takes the form $$u_{n+1} = u_n + hBu_n = (I + hB)u_n$$ where $I$ is the identity matrix of proper size. It follows that $$ u_n = (I + hB)^n u_0.$$ If the matrix $B$ is asymptotically stable, then the solution $y$ will decay to zero as $t$ tends to infinity. We desire that the computed solution $u$ exhibits the same behavior as $n$ tends to infinity. This is assured if the absolute value of the eigenvalues of the driving matrix, i.e., $G_h = I + hB$ are strictly less than unity.
In the case under discussion, $B \in \mathbb{R}^{999 \times 999}$ is the matrix $$ B = -\frac{1}{1000^2} \operatorname{tridiag}(-1,2,-1) $$ This matrix has eigenvalues which are strictly negative.
In general, loop as much information as possible into your symbols and avoid carrying signs around unless you absolutely have to.