= Solution
Write $h=\Delta x=1/N$ and collect only the interior values, so the actual vector dimension is $(N-1)^2$. The PDF's $N^2$ <matrix> size and its indices $k,l$ are inconsistent with the declared interior grid $m,k=1,\ldots,N-1$; neither affects the following stability argument when the interior indexing is used consistently.
Use the discrete <L2 norm> $\|u\|_h^2=h^2\sum_{m,k=1}^{N-1}|u_{m,k}|^2$, extending $u$ by zero to the boundary. Let $A$ be the horizontal flux operator and $B$ the vertical one. Each face coefficient is shared by its two adjacent grid values, so <summation by parts> gives
$$
\langle u,Au\rangle_h=-\sum_{k=1}^{N-1}\sum_{m=0}^{N-1}
a_{m+1/2,k}|u_{m+1,k}-u_{m,k}|^2,
$$
$$
\langle u,Bu\rangle_h=-\sum_{m=1}^{N-1}\sum_{k=0}^{N-1}
a_{m,k+1/2}|u_{m,k+1}-u_{m,k}|^2.
$$
All boundary-face terms are included. Positivity of $a$ makes $A,B$ symmetric negative definite on the interior grid: zero in either energy sum forces each corresponding grid line to be constant, and its zero boundary endpoint forces that constant to vanish. Consequently
$$
\frac d{dt}\|u(t)\|_h^2=2\operatorname{Re}\langle u,(A+B)u\rangle_h\le0,
\qquad
\boxed{\|u(t)\|_h\le\|u(0)\|_h.}
$$
This bound is independent of the spatial mesh, proving semidiscrete stability. No time integrator or time-step restriction has yet been introduced.
Back to article page
