BDF2 discrete energy identity Created 2026-10-06 Updated 2026-10-14
For three vectors in an inner product space, expansion of each squared norm proves
Therefore the second-order backward differentiation formula with a dissipative operator decreases the discrete energy . This quadratic form is equivalent to the product norm of the two time levels with constants independent of and . It proves uniform stability even when amplification roots of a polynomial coalesce inside the unit disk, where an eigenvector-separation proof may lose its bound.
Dissipative vector field 2026-10-05
A real vector field is dissipative in an inner product if
for all . Two solutions of the ordinary differential equation satisfy
Thus their distance never increases. The linear case recovers the dissipative operator condition. B-stability asks whether a Runge-Kutta method preserves this contraction for every positive step size.
The rotated global relations form a square system for the three missing boundary traces. Solve this system and substitute the resulting finite-time spectral boundary transforms into the integral representation. The pieces containing have zero contour contribution by analyticity, rotation of the relevant sectors and the Cauchy integral theorem. The remaining expression involves only . The finite-interval Airy spectral determinant has no nonremovable zeros in the contour domains, as established below, so this elimination introduces no extra spectral data there. Differentiating the representation, and using temporal inversion and compatibility at the corners, verifies the partial differential equation, initial value and the three prescribed boundary values.
For a direct Hadamard well-posedness check, the homogeneous difference satisfies
This proves uniqueness and continuous dependence for homogeneous boundary data. More generally subtract a smooth quadratic boundary lift satisfying the three prescribed endpoint conditions. The remaining problem has homogeneous boundary conditions and forcing . Its semigroup generator is on with . It is a dissipative operator; for every real the equation has a unique solution, since the homogeneous third-order boundary system is injective by the same energy identity and therefore invertible. Thus is maximal dissipative and generates a contraction semigroup. This gives existence and continuous dependence on the lifted data and forcing. Smooth compatible data therefore define a well-posed problem. This operator argument and the spectral elimination describe the same evolution.
A linear operator on a Banach space is called a dissipative operator when
On a Hilbert space, using the inner product linear in its first argument,
If , this proves the defining inequality. Conversely that inequality implies for every ; letting gives
For the next spectral assertion, dissipativity must be retained as a hypothesis. Let be closed, densely defined and dissipative. For , the dissipative bound makes injective. Its range is closed: if , the bound makes Cauchy, and then converges. Closedness produces a preimage of . If the range is all of , the same bound gives a bounded inverse, so is in the resolvent set. If the closed range is proper, there is a nonzero vector in its orthogonal complement, and the definition of the adjoint of a densely defined operator gives
This proves that the positive spectrum of a dissipative operator is adjoint point spectrum:
Without dissipativity the PDF's assertion is false. For example, multiplication by on is bounded, closed and densely defined. The value is spectral: normalized functions supported closer and closer to make tend to zero. But has no eigenvector for eigenvalue , since a square-integrable function supported on the singleton is zero.
For the given maximal-domain operator , finite-support sequences belong to , so the domain is dense. Coordinatewise limits also show it is closed: if and in , each defining coordinate equation passes to the limit, giving and hence .
Put and . For ,
Since by Cauchy-Schwarz inequality, has a finite limit . Its increments are square summable as well. Direct summation gives
The original inner-product series converges because . Taking therefore gives
Thus is dissipative, including its full maximal domain; the limiting boundary term has not been discarded.
To decide generation, solve for arbitrary and fixed . The first coordinate gives , and the rest give
Because all denominators are at least one,
Hence is bounded, so belongs to . The coordinate equations give , proving and surjectivity of . Together with density, closedness and dissipativity, the Lumer-Phillips theorem proves the contraction generator from an exponentially weighted difference operator conclusion:
The domain is crucial. The formal transposed matrix suggests a positive-eigenvalue solution
which is square summable because the products converge. It is not a true adjoint eigenvector: the admissible test vector has , so , while . These cannot agree for unless . This exemplifies how boundary terms exclude formal adjoint eigenvectors, and explains why a calculation using only finite-support test vectors would give the wrong generation answer.
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 . A time integrator then gives an amplification matrix , for example for a Runge-Kutta method. The useful stability requirement on a fixed physical time interval is
where is independent of the refining spatial and temporal meshes. The semidiscrete analogue is for . 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 is a normal matrix in the chosen discrete inner product, it has a unitary eigenbasis. Then
For dissipative problems it is sufficient that every 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 . Without that uniform bound the argument does not establish mesh stability. For finite-time estimates, an amplification modulus 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 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 , 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 on a unit interval with zero Dirichlet data, the centred second difference has a sine eigenbasis and eigenvalues
Explicit Euler method multiplies mode by . Requiring all factors to have modulus at most one yields on a given grid, and the sharp grid-independent sufficient limit is
These parabolic time-step restrictions express the growth of the largest-magnitude diffusion eigenvalue like . Backward Euler method instead has factor and the Crank-Nicolson method factor , both bounded by one for every . Their unconditional stability follows from A-stability, with much stronger damping of the shortest scales for backward Euler.
For , , explicit time stepping with a backward spatial difference gives
Thus the upwind finite difference scheme is contractive exactly for , 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 . Under a fixed positive hyperbolic Courant number, high-frequency data grow by a fixed factor in each of 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 , the excess per step is only , 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 symmetric negative definite, so a unitary eigenbasis still exists even without a closed eigenvalue formula. More generally, a dissipative operator in a discrete energy inner product obeys
which directly proves a uniform semigroup estimate. The dissipative Cayley-transform contraction also follows without normality: for ,
Taking 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. Let
Both eigenvalues are for every mesh, but applying the propagator to the second unit vector at gives norm at least . 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 for every , then for ,
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 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 and the error obeys , iteration gives
A residual tending uniformly to zero then gives convergence on . 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.
For a fully discrete linear partial differential equation evolution, write in a specified discrete norm. Finite-time stability of a numerical method means
with independent of the allowed meshes. For a multistep method, augment the state with its time-history values and include a stable starter. Bounds may grow as when the continuous problem has growth; demanding decay for all time would be a stronger assertion. The treatment of initial data, forcing, Dirichlet boundary conditions, and the norm is part of the hypothesis, not a detail supplied by an interior calculation.
For constant coefficients on an infinite or periodic uniform grid, von Neumann stability analysis substitutes the Fourier mode . The discrete Fourier transform and the Parseval identity turn a uniform multiplier estimate into a discrete L2 norm estimate. In a one-step scalar scheme, proves contractivity, while gives finite-time stability. In a multilevel scheme, every root of a polynomial in its amplification polynomial of a multilevel finite difference scheme matters, as does uniform control of the associated companion matrix.
For the Forward Euler diffusion scheme with ,
Thus for all frequencies exactly when
Sufficiency follows since ; necessity follows from the highest-frequency mode on even periodic grids. This gives the usual parabolic mesh restriction . For the Backward Euler diffusion scheme,
so every is stable. The Crank-Nicolson method has
again of modulus at most one for all , but stiff modes approach rather than zero. These are PDE manifestations of the A-stability and L-stability distinctions for time integration.
For advection , assume and set . The upwind finite difference scheme is , with
It is contractive when . This also has a direct maximum-norm proof: each newest value is a convex combination of two old values. The resulting Courant–Friedrichs–Lewy condition expresses that the numerical domain of dependence covers the physical one. For negative , the upwind direction must be reversed. By contrast, Forward Euler method time stepping with the centered first finite difference has and for nonzero . Under fixed nonzero refinement it is unstable, since a fixed nontrivial mode grows geometrically over steps. This is not a claim of instability under every imaginable mesh coupling: instead bounds its spurious finite-time growth, since .
The leapfrog advection scheme illustrates a multilevel subtlety. Its amplification polynomial of a multilevel finite difference scheme is
For , both roots of a polynomial have unit modulus and separation at least . A uniformly conditioned eigenbasis of the companion matrix then proves power boundedness of a two-level Fourier scheme for arbitrary history data. At and a grid admitting , the roots coincide on the unit circle, and a Jordan block produces growth. Thus the endpoint fails the ordinary arbitrary-history root condition for a multistep method; checking only that both moduli equal one misses the instability. Restricting the starter or filtering the parasitic mode is an additional hypothesis.
The energy method handles variable coefficients and finite boundaries for which Fourier stability analysis may be unavailable. For the Backward Euler diffusion scheme with homogeneous Dirichlet boundary conditions, define . The discrete summation by parts identity is
Take the inner product of with . The elementary identity gives
This proves unconditional contractivity including the actual boundary treatment. More generally, a dissipative operator gives of operator norm at most one: with , dissipativity implies , and the Cauchy-Schwarz inequality gives . In finite dimensions this also proves invertibility. For the second-order backward differentiation formula, the BDF2 discrete energy identity used in question 4 controls both time levels and proves unconditional diffusion stability. A repeated root strictly inside the disk is harmless here; the energy argument supplies the uniform bound without a singular eigenvector formula.
A second technique is eigenvalue stability analysis of a finite difference method. Under the method of lines, a time integrator advances by . If is a normal matrix, the scalar linear stability domain criterion on all controls its powers exactly in the corresponding L2 norm. For the centered Dirichlet discrete Laplacian, its eigenvalues lie between and zero; intersecting this interval with a time integrator's stability interval determines its mesh restriction. A uniform bound on the diagonalization of a matrix is needed for nonnormal diagonalizable systems. Eigenvalues alone can be misleading. For example,
has both eigenvalues inside the disk for , but the upper-right entry of is . Taking , , this grows like . Hence stable scalar eigenvalues do not give mesh-uniform stability. Boundary closures can introduce precisely the extra growth that an interior Fourier symbol overlooks, so boundary stability of a finite-difference method must be checked separately.
For nonlinear spatial discretizations, a useful replacement for Fourier stability analysis is a convexity argument. If the Forward Euler method map is nonexpansive in a chosen norm for , the second-order strong stability preserving Runge-Kutta method
is also nonexpansive under the same restriction. Indeed, for two inputs, the triangle inequality gives . Its Taylor expansion is , proving order two. Such strong stability preserving Runge-Kutta methods transfer suitable forward-step bounds without relying on a linear spectral calculation.
Finally, consistency of a numerical method explains what stability buys. For a linear Hadamard well-posed problem, the Lax equivalence theorem equates convergence of a consistent discretization with its stability, under the stated approximation-space and norm hypotheses. Directly, if the error satisfies , iteration gives the discrete Duhamel principle
Thus vanishing normalized local truncation error and initial error yield convergence. With a nonlinear Lipschitz step estimate, the discrete Gronwall inequality plays the same role. A stable but inconsistent stencil, such as the literal defective formulas in questions 3 and 4 under their usual refinement, is not rescued by any amplification bound. The decisive checks are a mesh-uniform evolution bound, a valid boundary closure, and consistency with the actual PDE.
For a closed densely defined dissipative operator, the estimate makes its range closed and its kernel trivial. Thus a positive spectral value means this range is proper, which is equivalent to a nonzero vector in . Dissipativity is essential: a real multiplication operator may have positive continuous spectrum with no adjoint eigenvector.
// Target: analysis.bigb