Solution (source code)

= Solution

Put $h=\Delta x$, $k=\Delta t$, with fixed positive $\mu=k/h$. Insert a smooth exact solution and divide the residual by $2k$:
$$
\frac{u(t+k)-u(t-k)}{2k}
-\frac{u(x+h,y)-u(x-h,y)+u(x,y+h)-u(x,y-h)}{2h}.
$$
<Taylor expansion> gives
$$
u_t-u_x-u_y+\frac{k^2}{6}u_{ttt}-\frac{h^2}{6}(u_{xxx}+u_{yyy})+O(k^4+h^4).
$$
The <differential equation> cancels the leading part, leaving normalized <local truncation error> $O(k^2+h^2)$. Since $u_{ttt}=(\partial_x+\partial_y)^3u$ contains mixed <derivatives>, there is no fixed positive <Courant number> canceling this leading error for every smooth solution. Thus \b[the method is second order in time and space]. Global second-order convergence also requires <stability> and initial values at both time levels accurate enough to retain this order.