Past exam of the mathematics course of the University of Cambridge / 2010 / iii / Paper 63 (Numerical Solution of Differential Equations)

For an autonomous ordinary differential equation, the chain rule gives $y^{\prime \prime }=f^{\prime }(y)f(y)$. The formula is therefore a multiderivative multistep method. Insert the exact solution and expand the two earlier values about $t=t_{n+2}$:The unscaled local truncation error isThe coefficients through $h^3$ vanish, but the fourth-order coefficient does not. Hence the method has order three, not four: the one-step residual is $O(h^{p+1})$ for an order-$p$ multistep formula. Its zero-stability polynomial is $(\zeta -1)(\zeta -1/7)$, with a simple unit root and the other root inside the unit disk. Thus the stated order is also the convergence order for smooth solutions, starting errors $O(h^3)$, and the nearby implicit solution branch, by convergence of a zero-stable multiderivative method.

Yes: this is an A-stable third-order method. Its second derivative is essential; the Second Dahlquist barrier applies to ordinary first-derivative linear multistep methods, not this multiderivative multistep method.

Apply the Dahlquist test equation $y^{\prime }=\lambda y$ and put $z=h\lambda$. The amplification polynomial isWe verify the root condition for a multistep method throughout $\mathrm{Re}z\le 0$. For a quadratic $a{\zeta }^2+b\zeta +c$ with $|a|>|c|$, the complex quadratic Schur criterion isfor both roots to be strictly inside the unit disk. Here is a brief justification of the criterion rather than just a root plot. Set $P^{\#}(\zeta )={\zeta }^2\overline{P(1/\overline{\zeta })}$. The polynomial $\overline{a}P-c{P}^{\#}$ is $\zeta (D\zeta +E)$, where $D=|a{|}^2-|c{|}^2$ and $E=\overline{a}b-c\overline{b}$. Its two roots are inside precisely when $|E|<D$. On the unit circle $|P^{\#}|=|P|$, and $|c|<|a|$, so the Rouche theorem equates its interior root count with that of $P$. A boundary zero of $P$ would also be a boundary zero of the transformed polynomial, which the strict inequality excludes. This proves the strict criterion; the non-strict form follows by continuity, with any unit root simple because $|c/a|<1$.

For the present polynomial it is enough to show $8|Q-1|<|Q{|}^2-1$, except at $z=0$. Write $z=-r+iy$, $r\ge 0$, and put $k=r^2+3r$, $s=y^2$. ThenDirect expansion gives the nonnegative factorizationEvery term is nonnegative, and the sum is strictly positive unless $r=y=0$. Also $|Q{|}^2=X+1>1$, so the criterion applies. For every nonzero $z$ in the closed left half-plane both amplification roots are strictly inside the unit disk. At $z=0$ they are $1$ and $1/7$, and the unit root is simple. ThereforeThis proof establishes stability of both the principal and parasitic modes, including on the imaginary axis.

In the variational convention needed here, a positive definite symmetric operator on a real Hilbert space satisfiesFor a bounded operator on the whole space, symmetry is self-adjointness. On a complex Hilbert space, the corresponding convention is Hermitian positivity; the quadratic expression is real. Strict positivity should be distinguished from uniformly positive definite symmetric operator coercivity, which requires a fixed $\gamma >0$ such that $⟨Lv,v⟩\ge \gamma \Vert v{\Vert }^2$.

Symmetry matters in a real space. If “positive definite” means only the inequality for the quadratic form, the next variational assertion is false: $L=\left(\begin{array}{cc}1& -1\\ 1& 1\end{array}\right)$ has $⟨Lv,v⟩=\Vert v{\Vert }^2>0$, but its skew part disappears from that form. Its Euler-Lagrange equation is $v=f$, not $Lv=f$. Thus the usual symmetric/Hermitian convention, or an explicit symmetry assumption, is required.

For a real Hilbert space, differentiate along an arbitrary direction $w$:The last equality uses the symmetry in the definition of a positive definite symmetric operator. Vanishing of this first variation for every $w$ is equivalent to $Lv=f$, givingIf $u$ solves it, expand the functional about $u$ to findThis is strictly positive for nonzero $w$, so $u$ is the unique minimizer. Conversely a minimizer has zero first variation and solves the equation. On a complex Hilbert space, use the real functional $⟨Lv,v⟩-2\mathrm{Re}⟨f,v⟩$ and real directional derivatives; varying $w$ and $iw$ yields the same result.

This proves the quadratic variational principle for a symmetric positive operator, conditional on existence. Strict positivity in infinite dimension does not by itself guarantee a minimizer for every $f$: on ${\ell }^2$, the operator $(Lv{)}_j=v_j/j$ is strictly positive, but $f_j=1/j$ would require the non-square-summable solution $v_j=1$. A coercive bilinear form supplies existence through the Lax-Milgram theorem under its usual boundedness assumptions.

Interpret the differential expression as $Lu=(p{u}^{\prime \prime }{)}^{\prime \prime }+qu$, with a suitable differential-operator domain and the stated clamped boundary conditions. For sufficiently smooth $u,v$ in that domain, twice integrating by parts givesThe boundary term vanishes because $v=v^{\prime }=0$ at both ends. For real coefficients this form is symmetric/Hermitian. In particularIndeed nonnegativity follows from $p>0$ and $q\ge 0$. If the integral were zero, $u^{\prime \prime }=0$ almost everywhere in the interval, so $u$ would be affine. Its endpoint values then force $u=0$. This proves strict positivity on the intended clamped domain and gives the weak formulation of a clamped fourth-order equation.

There is a genuine domain issue in the literal wording. The closure of clamped smooth functions in the stated L2 norm is all of $L^2(-1,1)$, since it contains the dense set $C_c^{\infty }(-1,1)$. Boundary traces are not retained by that closure, and a fourth-order expression does not define an operator on every element of $L^2$. For example $p=(1-x^2{)}^2$, $q=0$ satisfy the given sign conditions, but $u=1_{(-1/2,1/2)}$ belongs to this closure and makes $(p{u}^{\prime \prime }{)}^{\prime \prime }$ a distribution rather than an $L^2$ function. Even a general twice-smooth function need not have the four derivatives required by the classical expression.

A rigorous classical repair is to take smooth enough coefficients, for example $p\in C^2([-1,1])$, $q\in C([-1,1])$, and use a clamped domain such as $H^4\cap H_0^2$ in $L^2$, with the appropriate operator realization. A weak repair uses the clamped second-order Sobolev space $H_0^2$ and the form just displayed, with hypotheses ensuring its boundedness and closure; uniform positivity $p\ge p_0>0$ makes it coercive there. The integration-by-parts and strict-positivity proof does not need that stronger uniform bound, but existence and completion statements must not be inferred from the original $L^2$ closure alone. The printed $L[f]$ with an expression in $u$ is also a harmless variable-name mismatch; the calculation consistently uses $Lu$.

Put $h=\Delta x$. The spatial operator is the sum of a centred $x$ difference and a dissipative second-order forward advection stencil in $y$. For a smooth function, their Taylor expansions giveThus the leading spatial local truncation error is $h^2(u_{xxx}/6-u_{yyy}/3)$, generally nonzero. The semidiscretization is second order in space. Time has not yet been discretized, so there is no separate time order or time-step condition at this stage.

Use von Neumann stability analysis on the infinite grid. For the Fourier mode $u_{m,k}=\hat{u}(t)e^{i(m\xi +k\eta )}$, the semidiscrete eigenvalue isIts real part is nonpositive for every frequency. The mode evolves by $\hat{u}(t)=e^{t{\lambda }_h}\hat{u}(0)$, whose modulus does not grow. The discrete Fourier transform and Parseval identity therefore yieldwith a bound independent of $h$, where $\Vert u{\Vert }_{{\ell }_h^2}^2=h^2{\sum }_{m,k}|u_{m,k}{|}^2$. Hence the semidiscrete method is stable. Frequencies with $\eta =0$ are nondissipative rather than unstable.

The same sign can be seen directly by summation by parts: the centred $x$ operator is skew-adjoint, and for the unitary $y$ shift $S_y$,This is an energy version of the dissipative second-order forward advection stencil. A subsequent time integrator would need its own stability condition; semidiscrete stability does not assert that every fully discrete method is stable.

For $y^{\prime }=f(t,y)$ and step $h$, a collocation Runge-Kutta method constructs a polynomial $P$ of degree at most $s$ on $[t_n,t_n+h]$ satisfyingThe step value is $y_{n+1}=P(t_n+h)$. Thus the differential equation holds exactly at the prescribed collocation nodes, while the initial value fixes the integration constant. For a smooth vector field and small enough $h$, the implicit stage equations have the local solution branch continuing the zero-step initial state; arbitrary large steps need not have a unique nonlinear solution.

Let ${\ell }_j(\tau )={\prod }_{k\ne j}(\tau -c_k)/(c_j-c_k)$ be the Lagrange interpolation basis. Set $Y_i=P(t_n+c_{i}h)$ and $F_j=f(t_n+c_{j}h,Y_j)$. Since $P^{\prime }$ is a polynomial of degree at most $s-1$, the collocation conditions giveIntegrate from zero to $c_i$ and to one to obtainThese are the stage and update equations of a Runge-Kutta method with nodes $c_i$. Conversely, stages satisfying these equations reconstruct $P$ by integrating the displayed interpolation polynomial, so the equivalence works in both directions. Because ${\sum }_j{\ell }_j=1$, they also satisfy the collocation tableau row-sum identity ${\sum }_j{a}_{ij}=c_i$.

The collocation order theorem identifies the classical order with the order of its interpolatory quadrature. If ${\sum }_i{b}_{i}g(c_i)$ integrates every polynomial of degree at most $p-1$ exactly, but fails for some polynomial of degree $p$, then the collocation method has order $p$ for sufficiently smooth differential equations and the local implicit stage branch. With $s$ nodes, $s\le p\le 2s$.

Equivalently, put ${\pi }_s(\tau )={\prod }_i(\tau -c_i)$. If its first $k$ moments vanish,and the next does not, then $p=s+k$. Polynomial division explains the quadrature part: any polynomial through degree $s+k-1$ is a remainder of degree below $s$ plus ${\pi }_s$ times a polynomial of degree below $k$. The remainder is integrated exactly by interpolation, and the remaining integral vanishes. Conversely a nonzero next moment is an explicit failed quadrature test. The differential-equation order assertion is the collocation theorem being stated here. In particular Gauss collocation methods achieve $2s$, Radau IIA methods achieve $2s-1$, and endpoint-including Lobatto IIIA methods achieve $2s-2$ for $s\ge 2$.

The tableau as printed is not a collocation tableau. Its second stage row is $(3/4,3/4)$, whose sum is $3/2$, while the corresponding printed node is $1$. This violates the collocation tableau row-sum identity. The original PDF confirms this coefficient; it is not only a TeX transcription defect.

The required integrated basis at nodes $1/3,1$ isso direct integration givesThus replacing the last stage coefficient by $1/4$ gives the intended two-stage Radau IIA collocation method. Its quadrature moments areIt integrates quadratics exactly but not cubics, so the collocation order theorem gives order three. Independently, $b^{T}Ac=1/6$, confirming the remaining third-order Butcher order condition.

For completeness the literal printed Runge-Kutta method has only order one. For the autonomous equation $y^{\prime }=y$, the Taylor expansion of its stability function beginswhereas $e^z=1+z+z^2/2+O(z^3)$. Since $b^{T}1=1$, it is consistent, but this discrepancy already rules out order two. Therefore neither the requested equivalence nor the third-order conclusion is true without the stated coefficient repair.

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,\dots ,N-1$; neither affects the following stability argument when the interior indexing is used consistently.

Use the discrete L2 norm $\Vert u{\Vert }_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 givesAll 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. ConsequentlyThis bound is independent of the spatial mesh, proving semidiscrete stability. No time integrator or time-step restriction has yet been introduced.

Put $\delta =\Delta t$ and $L=A+B$, keeping the spatial grid fixed for the temporal error assertion. The matrix exponential expansions giveHence the Lie-Trotter splitting commutator error isThe energy estimate from part 1 makes $e^{tA}$, $e^{tB}$ and $e^{tL}$ contractions, so both $S_{\delta }$ and $E_{\delta }$ have norm at most one. Telescope their $n$th powers:For $n\delta \le T$, this givesThus the Lie-Trotter splitting has global temporal error $O(\Delta t)$. The constant here may depend on the fixed spatial operators; a mesh-uniform joint-limit error estimate needs additional regularity and commutator bounds. If $A$ and $B$ commute, this particular exponential splitting is exact, but variable diffusion coefficients generally remove that commutativity.

The $[1/1]$ Padé approximant of the exponential is $r(z)=(1+z/2)/(1-z/2)$. Replacing each exponential in its original order gives the split Crank-Nicolson methodwith the analogous definition for $B$. Operationally, first solve the $B$ substep and then the $A$ substep. The matrices $I-\delta A/2$, $I-\delta B/2$ are invertible for every $\delta \ge 0$, since $A,B$ are symmetric negative definite.

An orthonormal eigenbasis of $A$ diagonalizes $C_A$. Its eigenvalue for $\lambda \le 0$ isTherefore $\Vert C_A{\Vert }_h\le 1$, and similarly $\Vert C_B{\Vert }_h\le 1$. Submultiplicativity gives the contractivity of split Crank-Nicolson diffusion:No commutativity or simultaneous diagonalization is needed for this product bound. This proves unconditional discrete $L^2$ stability, not positivity of every stencil weight or unconditional maximum-norm monotonicity. Nor do the individually second-order rational substeps make the ordered splitting second order: when $[A,B]\ne 0$, its leading splitting defect remains ${\delta }^2[A,B]/2$.

A-stability removes the scalar linear decay restriction on the step size; it does not replace accuracy, nonlinear solvability or a suitable system norm. Its importance is clearest for a stiff differential equation, where rapidly decaying modes coexist with slowly varying quantities of interest. On the Dahlquist test equation $y^{\prime }=\lambda y$, the exact solution decays for $\mathrm{Re}\lambda <0$. A stable time integrator should not turn that decay into numerical growth.

For a one-step method write $y_{n+1}=R(z)y_n$, where $z=h\lambda$ and $R$ is its stability function. Its linear stability domain consists of parameters with $|R(z)|\le 1$ and a well-defined update. The method is A-stable when this domain contains the closed left half-plane. A consistent order-$p$ one-step method has $R(z)=e^z+O(z^{p+1})$ near zero, but that local agreement alone says nothing about stability at large negative $z$.

For a Runge-Kutta method with stage matrix $A$, weights $b$ and $1=(1,\dots ,1{)}^T$, eliminating the stages on the test equation givesAn explicit finite-stage method has a polynomial $R$. A nonconstant polynomial is unbounded along the negative real axis, so no consistent explicit Runge-Kutta method is A-stable. The explicit Euler method, with $R(z)=1+z$, is stable in the disk $|1+z|\le 1$ and on the negative real axis only for $-2\le z\le 0$. A stiff mode with a large negative eigenvalue therefore forces an unnecessarily small step even after its exact contribution has almost disappeared.

The Backward Euler method has $R(z)=1/(1-z)$. For $\mathrm{Re}z\le 0$, $|1-z{|}^2=1-2\mathrm{Re}z+|z{|}^2\ge 1$, so it is A-stable. The trapezoidal rule hasand is also A-stable. The implicit midpoint method has the same scalar stability function, although its nonlinear update differs from the trapezoidal rule. The first method has order one and the latter two order two. Their implicit solves trade greater work per step for freedom from the scalar stiff stability restriction.

Damping of very stiff modes separates A-stability from L-stability. For a one-step method, the latter additionally requires $R(z)\to 0$ as $z$ tends to infinity in the left half-plane. Backward Euler method satisfies this. The trapezoidal and midpoint rules instead have $R(z)\to -1$ along the negative real axis: a very stiff mode can alternate in sign with little damping even though its amplitude remains bounded. L-stability is often preferable for diffusion and rapid relaxation, while retaining oscillatory modes can be desirable in conservative dynamics.

For a linear multistep method, writeIts test-equation amplification roots satisfy $\rho (\zeta )-z\sigma (\zeta )=0$. Absolute stability requires every root to lie in the unit disk or on its boundary with all unit-modulus roots simple; all parasitic modes matter. At $z=0$ this is the zero-stability root condition. Consistency requires $\rho (1)=0$ and ${\rho }^{\prime }(1)=\sigma (1)$, and the Dahlquist equivalence theorem combines consistency and zero-stability to characterize convergence. Thus a claim of stability obtained by canceling a common factor is not a proof about the original recurrence and its arbitrary consistent starting data.

Explicit linear multistep methods cannot be A-stable either. If ${\beta }_k=0$ and ${\alpha }_k\ne 0$, some lower coefficient of $\rho -z\sigma$ becomes unbounded as $z\to -\infty$, while the leading coefficient stays fixed. If all roots remained in the unit disk, every coefficient-to-leading-coefficient ratio, an elementary symmetric function of those roots, would remain bounded. This contradiction forces an unstable root.

An important implicit example is BDF2 method, whose characteristic equation is $(3-2z){\zeta }^2-4\zeta +1=0$. A unit root $\zeta =e^{i\theta }$ givesThe leading coefficient vanishes only at $z=3/2$, outside the left half-plane. Near a small negative $z$, both roots are inside the disk; no unit-circle crossing can occur in the connected open left half-plane. On the imaginary axis the boundary equality can occur only at $z=0$, where the roots are $1,1/3$. Hence BDF2 is A-stable and its two roots tend to zero for very large negative $z$, giving strong stiff-mode damping. This is a complete boundary-locus argument, not just a boundary plot.

The Second Dahlquist barrier says that a convergent A-stable ordinary linear multistep method has order at most two. This theorem does not limit implicit Runge-Kutta methods or multiderivative multistep methods. In particular the A-stable third-order two-step multiderivative method in question 1 uses $y^{\prime \prime }$ and satisfies different order conditions. Among collocation Runge-Kutta methods, the standard family results give Gauss collocation methods order $2s$ with diagonal $[s/s]$ Padé approximants to $e^z$ as stability functions, and Radau IIA methods order $2s-1$ with subdiagonal $[s-1/s]$ approximants. Both families are A-stable; the Gauss functions have a nonzero limit at infinity, whereas Radau IIA is L-stable. The corrected two-stage right-Radau method, for example, haswhich matches the exponential through cubic order and tends to zero at infinity. Its poles $2\pm i\sqrt{2}$ lie in the right half-plane. For $z=-r+iy$, $r\ge 0$, direct subtraction givesproving its A-stability directly.

Sectorial A-alpha stability relaxes the requirement to a sector around the negative real axis, useful when the spectrum of a stiff problem remains in that sector; higher-order backward differentiation formulas illustrate the usefulness of this weaker requirement. For a normal linear system $y^{\prime }=Ly$, unitary diagonalization reduces the numerical update to its scalar factors $R(h{\lambda }_j)$. For a nonnormal system, a scalar eigenvalue test alone does not give a uniform bound in the chosen norm. Nonlinear stability is different again: B-stability asks for contractivity when $\mathrm{Re}⟨f(u)-f(v),u-v⟩\le 0$, and algebraic stability of a Runge-Kutta method is a standard sufficient criterion. Thus A-stability addresses a precise and valuable linear test, while step selection still has to meet accuracy requirements and the actual implicit equations must be solved on the relevant branch.

Eigenvalue analysis proves mesh-uniform stability when it is accompanied by control of eigenvectors or an energy norm; the signs of eigenvalues alone are insufficient. For an evolution equation, the method of lines first produces $u_h^{\prime }=L_h{u}_h$. A time integrator then gives an amplification matrix $S_{h,\delta }$, for example $S_{h,\delta }=R(\delta L_h)$ for a Runge-Kutta method. The useful stability requirement on a fixed physical time interval iswhere $C_T$ is independent of the refining spatial and temporal meshes. The semidiscrete analogue is $\Vert e^{t{L}_h}{\Vert }_h\le C_T$ for $0\le t\le T$. These are bounds on propagation of initial errors and residuals, not merely statements that each finite matrix has a bounded solution over its own finite time interval.

If $L_h$ is a normal matrix in the chosen discrete inner product, it has a unitary eigenbasis. ThenFor dissipative problems it is sufficient that every $\delta {\lambda }_j$ lie in the integrator's linear stability domain. A diagonalizable family also suffices if its eigenvector matrices have uniformly bounded condition numbers: the same maximum is multiplied by $\Vert V_h\Vert \Vert V_h^{-1}\Vert$. Without that uniform bound the argument does not establish mesh stability. For finite-time estimates, an amplification modulus $1+O(\delta )$ can also be acceptable; strict contraction at every step is a stronger property.

For a constant-coefficient stencil on the full lattice or a periodic grid, von Neumann stability analysis provides the eigenvectors explicitly. A Fourier mode diagonalizes each translation-invariant spatial difference, giving a scalar symbol ${\lambda }_h(\theta )$ or a small matrix symbol for a system. The discrete Fourier transform and Parseval identity turn a uniform modal bound into an actual discrete L2 norm bound. The two-dimensional stencil in question 3 has $\mathrm{Re}{\lambda }_h=-(1-\cos \eta {)}^2/h\le 0$, so its continuous-time amplification is contractive; a later time integrator must still be checked. With boundaries one must analyze the boundary rows as well: a whole-line Fourier symbol cannot detect an unstable boundary closure.

For the heat equation $u_t=\kappa u_{xx}$ on a unit interval with zero Dirichlet data, the centred second difference has a sine eigenbasis and eigenvaluesExplicit Euler method multiplies mode $j$ by $1+\delta {\lambda }_j$. Requiring all factors to have modulus at most one yields $\kappa \delta /h^2\le 1/[2{\sin }^2((N-1)\pi /(2N))]$ on a given grid, and the sharp grid-independent sufficient limit isThese parabolic time-step restrictions express the growth of the largest-magnitude diffusion eigenvalue like $h^{-2}$. Backward Euler method instead has factor $(1-\delta {\lambda }_j{)}^{-1}$ and the Crank-Nicolson method factor $(1+\delta {\lambda }_j/2)/(1-\delta {\lambda }_j/2)$, both bounded by one for every $\delta \ge 0$. Their unconditional stability follows from A-stability, with much stronger damping of the shortest scales for backward Euler.

For $u_t=-a{u}_x$, $a>0$, explicit time stepping with a backward spatial difference givesThus the upwind finite difference scheme is contractive exactly for $0\le \nu \le 1$, a Courant–Friedrichs–Lewy condition matching the direction and speed of information propagation. A centred spatial difference has purely imaginary eigenvalues, and forward Euler gives $|S{|}^2=1+{\nu }^2{\sin }^2\theta >1$. Under a fixed positive hyperbolic Courant number, high-frequency data grow by a fixed factor in each of $O(1/h)$ steps, proving instability as the mesh is refined. This is not a claim that no exceptionally small coupled step can ever give a finite-time bound: if $\delta =O(h^2)$, the excess per step is only $O(\delta )$, but the usual hyperbolic scaling is unstable. Higher-order explicit integrators can possess nontrivial imaginary-axis stability intervals and are then suitable for skew-adjoint advection or the energy-form first-order wave equation.

Variable coefficients usually prevent Fourier diagonalization, but may preserve a self-adjoint spatial operator. The positive-face diffusion discretization in question 5 is an example. Its discrete energy identity makes $A+B$ symmetric negative definite, so a unitary eigenbasis still exists even without a closed eigenvalue formula. More generally, a dissipative operator $L_h$ in a discrete energy inner product obeyswhich directly proves a uniform semigroup estimate. The dissipative Cayley-transform contraction also follows without normality: for $\alpha \ge 0$,Taking $v=(I-\alpha L_h{)}^{-1}u$ proves contraction of the Crank–Nicolson factor. The contractivity of split Crank-Nicolson diffusion then follows from multiplying two contraction bounds; simultaneous eigenvectors of the two directional matrices are unnecessary.

The negative spectra do not imply uniform semidiscrete stability example shows what can go wrong. LetBoth eigenvalues are $-1$ for every mesh, but applying the propagator to the second unit vector at $t=1$ gives norm at least $1/(eh)$. The evolution is therefore not mesh-uniformly stable. For a discrete matrix, a unit-modulus eigenvalue with a nontrivial Jordan block gives powers growing polynomially in the step number. Even strictly interior eigenvalues can allow large transient growth when eigenvectors are badly conditioned. A pseudospectrum or numerical range of an operator, and especially a direct energy estimate, can reveal behavior missed by eigenvalues.

A resolvent estimate provides another way to see the needed uniformity. If $\Vert S_h^n\Vert \le C$ for every $n$, then for $|z|>1$,This shows directly that a uniform power estimate imposes more than the location of the spectrum. Conversely, resolvent-to-power estimates in varying dimension must themselves be uniform; a bound with constants depending on the number of grid points is not enough.

For a linear multistep method, each spatial eigenvalue gives an amplification polynomial, and all roots must satisfy the root condition for a multistep method, including simplicity of any unit-modulus root. Merely following the root approximating $e^{\delta \lambda }$ misses parasitic modes. The temporal and spatial analyses must therefore be combined with the appropriate zero-stability and uniform modal bounds.

Finally stability connects to convergence. If the exact grid samples have normalized residual ${\tau }^n$ and the error obeys $e^{n+1}=S_h{e}^n+\delta {\tau }^n$, iteration givesA residual tending uniformly to zero then gives convergence on $[0,T]$. For a well-posed linear initial-value problem and a consistent linear finite-difference approximation, this is the forward implication of the Lax equivalence theorem. The practical conclusion is to analyze the full operator including boundaries, select a time method whose stability region contains the relevant scaled spectrum, and verify that eigenvectors or an energy estimate make the resulting bounds uniform as the grid is refined.