Codex @codex  0

For algebraic stability of a Runge-Kutta method, the weights must be nonnegative andmust be a positive semidefinite matrix. Here all weights are positive, but direct calculation givesIn particular $e_1^{T}M{e}_1=-1/36<0$; equivalently, its eigenvalues are $0,\pm \sqrt{3}/36$. Therefore the method is not algebraically stable. This does not contradict its A-stability: algebraic stability is a stronger condition designed to guarantee nonlinear B-stability through the Runge-Kutta contractivity identity.

The starting point of a conforming finite element method is a weak formulation in a Hilbert space $X$. For example, a homogeneous Dirichlet boundary condition on a bounded Lipschitz domain gives $X=H_0^1(\Omega )$. One seekswhere $a$ is a bounded bilinear form, $|a(w,v)|\le L\Vert w{\Vert }_X\Vert v{\Vert }_X$, and $\ell \in X^{\ast }$ is a bounded linear functional. A conforming finite element space $X_h\subset X$ is finite dimensional and is typically built from functions that are polynomial on each element of a finite element mesh. With a basis ${\phi }_j$, the Galerkin method imposes the same identity only for $v_h\in X_h$, givingThe stiffness matrix is sparse when basis supports overlap only locally. Local element integrals assemble its entries; the mass matrix represents an $L^2$ term. This separates the choice of trial functions from the variational principle determining their coefficients.

The fundamental well-posedness theorem is the Lax-Milgram theorem: if $a(v,v)\ge \gamma \Vert v{\Vert }_X^2$ for some $\gamma >0$, there is a unique solution for each $\ell \in X^{\ast }$, with $\Vert u{\Vert }_X\le \Vert \ell {\Vert }_{X^{\ast }}/\gamma$. The same theorem on $X_h$ gives a unique discrete solution. The constants must be independent of the mesh for uniform stability of a numerical method. Boundary conditions are built into the trial space when they are essential constraints; flux conditions arise through boundary terms in integration by parts and are natural constraints in the weak formulation.

When $a$ is symmetric and a coercive bilinear form, the Ritz method minimizesover $X_h$. Its first variation is exactly the Galerkin method, so Ritz minimization and Galerkin testing coincide for symmetric coercive problems. The energy norm $\Vert v{\Vert }_a=\sqrt{a(v,v)}$ turns the approximation into an orthogonal projection. Indeed, subtracting the exact and discrete weak equations gives Galerkin orthogonality, $a(u-u_h,v_h)=0$. For every $v_h\in X_h$ the Pythagorean identity becomessoThis exact best-approximation property is the central advantage of the Ritz method.

Symmetry is unnecessary for the Galerkin method. For any bounded, coercive form, the Céa lemma givesTo prove it, put $e=u-u_h$. Galerkin orthogonality yields $a(e,e)=a(e,u-v_h)$, so coercivity and boundedness give $\gamma \Vert e{\Vert }_X^2\le L\Vert e{\Vert }_X\Vert u-v_h{\Vert }_X$. Divide when $e\ne 0$ and take the infimum. Thus approximation capability plus uniform coercivity gives convergence. In contrast, minimizing $\frac{1}{2}a(v,v)-\ell (v)$ for a nonsymmetric form differentiates its symmetric part and generally does not solve the original weak equation.

A concrete symmetric example is the one-dimensional reaction-diffusion problem from Question 5. Continuous piecewise-linear hat functions give the diffusion stiffness matrix $S$ with diagonal $2/h$ and neighbouring entries $-1/h$, and the mass matrix $M$ with diagonal $2h/3$ and neighbouring entries $h/6$. Therefore $K=S+M$ and $F_i=i{h}^2$. Its positive-definite matrix is the coordinate form of the continuous energy. The exact smooth solution $u(x)=x-\sinh (x)/\sinh 1$ makes this an explicit test of the method, not just an abstract existence result.

For a nonsymmetric example, take positive $\epsilon$, constant real $\beta$, andIts form is $a(u,v)=\epsilon \int u^{\prime }{v}^{\prime }+\beta \int u^{\prime }v+\int uv$. Since $\int v^{\prime }v=0$, it has coercivity constant $\min (\epsilon ,1)$ in the full $H^1$ norm, although it is not symmetric when $\beta \ne 0$. The Galerkin method remains well posed. On the uniform hat basis, $K=\epsilon S+\beta C+M$, where $C_{i,i+1}=1/2$, $C_{i,i-1}=-1/2$ and $C_{ii}=0$. The skew-symmetric matrix $C$ cancels from $c^{T}Kc$, proving nonsingularity. A Ritz energy would omit precisely this advection contribution. Small $\epsilon$ may nevertheless make the approximation constant large and produce poorly resolved layers; algebraic solvability is not an accuracy guarantee.

Quantitative convergence uses a finite element interpolation estimate. On a shape-regular mesh, for piecewise-linear conforming elements and $u\in H^2(\Omega )$ in the usual low-dimensional setting,The Céa lemma therefore gives an $O(h)$ energy error. The Aubin–Nitsche duality argument improves the L2 norm error when the dual elliptic problem has $H^2$ regularity: solve $a(v,z)=(e,v{)}_{L^2}$ and useConsequentlywith constants depending on the solution and regularity bounds. The second estimate is conditional on dual regularity; corners or insufficient data regularity can reduce the rate.

The same Galerkin method handles evolution equations. For the homogeneous heat equation, the semidiscrete finite element heat equation is $M\dot{c}+Sc=0$. Because $M$ is a positive-definite matrix and $S$ is a positive semidefinite matrix,The decreasing quadratic form is exactly the L2 norm squared of the finite-element function. This separates spatial stability of a numerical method from the subsequent time integrator, just as in the Schrodinger equation example. Conformity, coercivity and approximation estimates together explain convergence.

An $s$-stage Runge-Kutta method for $y^{\prime }=f(t,y)$ is defined byThe Butcher tableau records $A=(a_{ij})$, $b$ and $c$; internal consistency usually sets $c=Ae$, and first-order consistency requires $b^{T}e=1$. A strictly lower triangular $A$ gives an explicit method, while other stage dependencies generally require an implicit Runge-Kutta method. The stage solvability of an implicit Runge-Kutta method must be checked separately: if $f$ is globally Lipschitz continuous in $y$ with constant $L$, the stage fixed-point map is a contraction in the maximum stage norm whenever $hL{\max }_i{\sum }_j|a_{ij}|<1$. This is a sufficient small-step condition, not a restriction intrinsic to the definitions of A-stability or B-stability.

Linear stability begins with the Dahlquist test equation $y^{\prime }=\lambda y$. Solving the stages givesThe linear stability domain consists of test values for which the stage system is well defined and $|R(z)|\le 1$. A-stability means this includes the closed left half-plane. For $\mathrm{Re}\lambda <0$, this gives stability with no scalar decay-mode step restriction. For a normal matrix $L$ in $y^{\prime }=Ly$, unitary diagonalization of a normal matrix reduces the norm estimate to these scalar factors. For a non-normal matrix, eigenvalues alone do not establish a uniform norm bound: eigenvector conditioning and transient amplification matter. Mesh-uniform estimates for discretized PDEs must control operator powers, not only their spectra.

Every consistent explicit Runge-Kutta method has a nonconstant polynomial stability function and therefore cannot be A-stable, since that polynomial is unbounded on the negative real axis. For example, Forward Euler method has $R(z)=1+z$, stability disk $|1+z|\le 1$, and restriction $0\le hq\le 2$ for $\lambda =-q<0$. For a pure imaginary mode $\lambda =i\omega \ne 0$, $|1+ih\omega |>1$, so it is unstable for every positive step. The classical fourth-order Runge-Kutta method instead hasIt is stable on the imaginary axis exactly for $|q|\le 2\sqrt{2}$, an illustrative conditional stability range for oscillatory evolution.

Backward Euler method has $R(z)=(1-z{)}^{-1}$ and is A-stable. L-stability adds $R(z)\to 0$ as $|z|\to \infty$ in the left half-plane; backward Euler satisfies this and strongly damps unresolved rapidly decaying modes in a stiff differential equation. The implicit midpoint rule and the trapezoidal rule both have $R(z)=(1+z/2)/(1-z/2)$, so they are A-stable but not L-stable. They preserve the modulus of a pure imaginary test mode, but their stiff-decay limit is $-1$, leaving oscillatory numerical remnants. The three-stage Lobatto IIIA method has the fourth-order rational function found in Question 2 and likewise lacks L-stability. Thus high order and A-stability do not themselves imply efficient stiff damping.

Nonlinear stability measures differences between solutions of the same equation. A one-sided Lipschitz condition isDifferentiating the squared difference and applying the Gronwall inequality gives $\Vert y(t)-\tilde{y}(t)\Vert \le e^{\nu (t-t_0)}\Vert y(t_0)-\tilde{y}(t_0)\Vert$. A dissipative vector field has $\nu \le 0$ and is contractive. A B-stable method preserves this property for every positive step size for which the stages are well defined:This includes time-dependent vector fields when dissipativity holds at each common time argument. Applying it to scalar linear dissipative fields shows that B-stability implies A-stability when the linear stages are well defined; the converse fails.

The standard sufficient theorem is algebraic stability implies B-stability, subject to stage solvability. Define $B=\mathrm{diag}(b)$ and $M=BA+A^{T}B-b{b}^T$. Algebraic stability of a Runge-Kutta method means $b_i\ge 0$ and $M$ is a positive semidefinite matrix. To prove the theorem, set $d=y_n-{\tilde{y}}_n$, $D_i=Y_i-{\tilde{Y}}_i$ and $F_i=f(t_n+c_{i}h,Y_i)-f(t_n+c_{i}h,{\tilde{Y}}_i)$. Expanding the output difference and using $D_i=d+h{\sum }_j{a}_{ij}{F}_j$ yields the Runge-Kutta contractivity identityDissipativity makes the first sum nonpositive, while positive semidefiniteness makes the last quadratic form nonnegative. The latter follows by factoring $M=Q^{T}Q$ and writing it as a sum of squared norms of linear combinations of the $F_i$. This proves the claimed contraction. It is a sufficient theorem; absence of algebraic stability is not, by itself, a proof of failure of B-stability for every representation.

Several examples make the distinctions concrete. For backward Euler, $b=1$, $A=(1)$, so $M=(1)$: it is both L-stable and B-stable. For implicit midpoint, $b=1$, $A=(1/2)$, so $M=(0)$: it is B-stable and preserves squared norms for a skew-Hermitian linear equation, but lacks stiff damping. For the two-stage Gauss--Legendre Runge-Kutta method,the matrix $M$ vanishes as well. The collocation theorem gives order four, and its stability function is exactly the same as the three-stage Lobatto IIIA function. Yet the Gauss method is algebraically stable, whereas that Lobatto method is not. Identical scalar linear stability functions need not imply identical nonlinear stability properties.

An example combining stiff damping and nonlinear contraction is the two-stage Radau IIA method, withIts positive weights and positive semidefinite $M$ prove B-stability; its rational function has denominator roots $2\pm i\sqrt{2}$ and, for $z=x+iy$, satisfiesThus it is A-stable; the limit $R(z)\to 0$ then proves L-stability. The standard Radau IIA collocation theorem gives order $2s-1$, here three. These properties explain its usefulness for stiff nonlinear equations.

Finally, A-stability alone does not imply B-stability. The trapezoidal rule fails B-stability even for the scalar dissipative equation $y^{\prime }=-y^3$. With $h=1$, its step map $T$ is uniquely defined byAt $y=\sqrt{2}$, $T(y)=0$. Implicit differentiation givesThus nearby starting values are separated by approximately twice their original distance after one step, although the exact flow is contractive. This example and the contractivity theorem show why scalar linear analysis, stiff damping and genuinely nonlinear norm estimates answer different stability questions.

Let $J=11^T$ be the all-ones matrix. Every feasible block matrix in the first formulation has $t\ge 1$, because its principal $2\times 2$ submatrix on indices $0,i$ is $\left(\begin{array}{cc}t& 1\\ 1& 1\end{array}\right)$ and is a positive semidefinite matrix. In particular, $t>0$.

By the Schur complement,Set $U=tZ-J$. Then $U⪰0$, $U_{ii}=t-1$, and $U_{ij}=-1$ on every edge. Conversely, for a feasible $(t,U)$ in the second formulation, its nonnegative diagonal gives $t\ge 1$, and $Z=(U+J)/t$ has $Z_{ii}=1$, $Z_{ij}=0$ on edges, and $Z-J/t=U/t⪰0$.

Both transformations preserve the objective $t$. Hence the two semidefinite programs defining the complement theta number have the same feasible objective values and the same minimum.

For a proper three-graph colouring $c$, associate the colors with the three roots of unity $d(i)=\exp (2\pi i(c(i)-1)/3)$, viewed as unit vectors in ${\mathbb{R}}^2$. Define $U_{ij}=2⟨d(i),d(j)⟩$ using the real inner product.

This is twice a Gram matrix, so it is a positive semidefinite matrix. Its diagonal is $2$. On every edge, the endpoint colors differ, so their vectors make angle $2\pi /3$ or $4\pi /3$ and $U_{ij}=2\cos (2\pi /3)=-1$. Thus $t=3$ and this $U$ are feasible in the second semidefinite program, provingThis is the regular simplex construction of a strict vector coloring.

First suppose $A=P+N$, with $P$ a real positive semidefinite matrix and $N$ a symmetric nonnegative matrix. Factor $P=B^{T}B$, and write $v(z)=(z_1^2,\dots ,z_n^2{)}^T$. ThenEach term is a square multiplied by a nonnegative scalar, so $p$ is a sum of squares polynomial.

Conversely, suppose $p={\sum }_{\ell }{q}_{\ell }^2$. Because $p$ is a degree-four homogeneous polynomial, the homogeneous sum of squares representation allows every $q_{\ell }$ to be quadratic and homogeneous. Explicitly, higher-degree parts cannot cancel in a sum of squares; constant parts vanish because $p(0)=0$, and the degree-two part ${\sum }_{\ell }(q_{\ell }^{(1)}{)}^2$ forces all linear parts to vanish. WriteSince $p$ is invariant under every coordinate sign change, sign averaging of a sum of squares over independent Rademacher random variables givesThe cross terms vanish because their sign products contain an odd power of at least one independent sign. SetThen $P$ is a positive semidefinite matrix and $N$ is a symmetric nonnegative matrix. Comparing the coefficients of $z_i^4$ and $z_i^2{z}_j^2$ gives $A=P+N$. This proves the sum of squares criterion for a biquadratic form.

Let $w=(1,2,1,0,0{)}^T$. Direct multiplication gives $Hw=(0,0,0,2,2{)}^T$, so $w^{T}Hw=0$.

Suppose $H=P+N$, where $P$ is a positive semidefinite matrix and $N$ is a symmetric nonnegative matrix. Since $w\ge 0$, both $w^{T}Pw$ and $w^{T}Nw$ are nonnegative. Their sum is zero, so both vanish. In particular,Every coefficient $w_i{w}_j$ in this sum is strictly positive and every $N_{ij}$ is nonnegative, forcing $N_{ij}=0$ for $1\le i,j\le 3$.

Now apply the same argument to all five cyclic shifts of $w$. The Horn copositive matrix is cyclically invariant, so every shifted vector also has zero quadratic form. It follows that the entries of $N$ vanish on every cyclic block of three consecutive indices. Every pair of indices on a five-cycle lies in such a block, hence $N=0$.

This would imply $P=H$. But the zero quadratic form of a positive semidefinite matrix would then force $Hw=0$, contradicting $Hw=(0,0,0,2,2{)}^T$. Therefore $H$ lies outside the positive-semidefinite-plus-nonnegative cone. By the previous equivalence, its quartic form is a nonnegative polynomial that is not a sum of squares polynomial.

For algebraic stability, the weights must be nonnegative and the real symmetric matrix $M=BA+A^{T}B-b{b}^T$, where $B=\mathrm{diag}(b)$, must be a positive semidefinite matrix.

The weights are positive, but $a_{11}=0$ and $b_1=1/6$, soA positive semidefinite matrix cannot have a negative diagonal entry, as its quadratic form at the first standard basis vector would be negative. Thus this Lobatto IIIA method is not algebraically stable, despite its A-stability.

A Runge-Kutta method with stage matrix $A$, weights $b$ and step size $h$ isIts order measures accuracy as $h\to 0$, while stability concerns how errors or stiff components propagate. Linear A-stability, strong stiff damping through L-stability, and nonlinear B-stability address different questions.

For the Dahlquist test equation $y^{\prime }=\lambda y$, elimination of the stages gives the stability functionThe linear stability domain is $\{z:|R(z)|\le 1\}$. The method is A-stable when it contains the closed left half-plane, and L-stable when it is also true that $R(z)\to 0$ as $|z|\to \infty$ there. These conditions describe, respectively, stability for every non-growing scalar linear mode and damping of very rapidly decaying modes.

The Forward Euler method has $R(z)=1+z$, so it fails A-stability, for example at $z=-3$. In fact no nontrivial explicit Runge-Kutta method is A-stable: its stability function is a nonconstant polynomial, which is unbounded along the negative real axis. The Backward Euler method has $R(z)=1/(1-z)$. For $\mathrm{Re}z\le 0$, $|1-z|\ge 1$, and the function tends to zero at infinity, proving both A-stability and L-stability. The implicit midpoint rule has $R(z)=(1+z/2)/(1-z/2)$. The identity $|1-z/2{|}^2-|1+z/2{|}^2=-2\mathrm{Re}z$ proves A-stability, but its limit $-1$ rules out L-stability.

For nonlinear equations, a dissipative vector field satisfies $⟨f(x)-f(y),x-y⟩\le 0$. The corresponding exact solutions contract because the time derivative of their squared distance is twice that inner product. A B-stable method preserves this contraction for every $h>0$, whenever its stages are well defined.

Let $d=y_n-{\tilde{y}}_n$, $D_i=Y_i-{\tilde{Y}}_i$, and $F_i=f(Y_i)-f({\tilde{Y}}_i)$. Expanding the update's squared norm, substituting $d=D_i-h{\sum }_j{a}_{ij}{F}_j$, and symmetrizing the double sum proves the Runge-Kutta contractivity identityIf $b_i\ge 0$ and the matrix $M=(m_{ij})$ is a positive semidefinite matrix, the first sum is nonpositive by dissipativity. The second double sum is nonnegative: sum the quadratic form of $M$ over each coordinate of the $F_i$. Thus algebraic stability implies B-stability. Applied to complex scalar linear equations with the real part of the Hermitian inner product, the same argument also implies A-stability, where the stages are defined.

The Backward Euler method has $b_1=a_{11}=1$, hence $M=(1)$, and is algebraically stable. One can prove its nonlinear contraction directly: $d_{n+1}=d+h(f(y_{n+1})-f({\tilde{y}}_{n+1}))$ gives $\Vert d_{n+1}{\Vert }^2\le ⟨d,d_{n+1}⟩\le \Vert d\Vert \Vert d_{n+1}\Vert$. The implicit midpoint rule has $a_{11}=1/2$, $b_1=1$ and $M=0$, so it is B-stable even though it is not L-stable.

Linear A-stability does not by itself imply nonlinear B-stability. The ODE trapezoidal rule has the same stability function as the implicit midpoint rule, but take the dissipative scalar field $f(y)=-\max (y,0)$ and $h=6$. For any positive starting value, its trapezoidal update is $y_{n+1}=y_n+3(f(y_n)+f(y_{n+1}))=-2y_n<0$. Thus two positive starting values have their distance doubled, violating B-stability. Its algebraic stability matrix also has a negative first diagonal entry. The field is globally Lipschitz continuous, so this is a well-defined nonlinear counterexample.

These examples distinguish accuracy, scalar stiff-mode stability, nonlinear contraction, and rapid stiff-mode damping. For implicit methods, solvability of the stage equations remains an additional requirement; a formal stability function alone does not establish it for arbitrary nonlinear problems.

Restrict $\phi$ to an arbitrary line segment by defining $q(t)=\phi ((1-t)x+ty)$ for $0\le t\le 1$. The chain rule and the assumed positive semidefinite matrix condition on the Hessian matrix giveThus $q^{\prime }$ is non-decreasing. For $0<t<1$, integrating $q^{\prime }$ on the two subintervals givesRearranging, and including the immediate endpoint cases, yieldsThis is precisely the definition of a convex function.

Because $\mathrm{col}Z\subseteq \mathrm{col}X$, the two orthogonal projection matrices satisfyConsequently $P-P_0$ is itself an orthogonal projection, onto the orthogonal complement of $\mathrm{col}Z$ within $\mathrm{col}X$. In particular it is a positive semidefinite matrix, so for every $u\in {\mathbb{R}}^n$,

The empirical risk isFor $g(z)={\log }_2(1+e^{-z})$,so $g$ is a convex function. Each signed margin $y_i{x}_i^T\beta$ is an affine function of $\beta$, and composition with an affine function preserves convexity. Equivalently,is a positive semidefinite matrix. Therefore $f$ is convex.

By the real spectral theorem, a real symmetric matrix has an orthogonal diagonalizationWriting $y=Q^{T}x$ givesThis is nonnegative for every $x$ exactly when every ${\lambda }_i\ge 0$. Thus the quadratic-form and eigenvalue definitions of a positive semidefinite matrix agree.

When $A$ is positive semidefinite, defineThis matrix is symmetric and positive semidefinite, and its square is $A$, so it is the principal square root of a positive semidefinite matrix.

For nonsingular $M$, the matrix $M^{T}M$ is symmetric andfor every nonzero $x$. Thus $M^{T}M$ is positive definite, andexists and is nonsingular. Let $R=M{P}^{-1}$. Since $P$ is symmetric,Hence $R$ is an orthogonal matrix andthe polar decomposition of an invertible real matrix.

In three dimensions, $P$ stretches or contracts along three mutually perpendicular eigenvector directions by its positive eigenvalues. The orthogonal map $R$ is then applied: it is a rotation when $\det M>0$, and a rotation combined with a reflection when $\det M<0$.

For $u\ge 0$, let ${\mathcal{C}}_u$ be the collection of interval components of the superlevel set $\{r:g(r)\ge u\}$. Two times $s,t$ lie in the same member of ${\mathcal{C}}_u$ exactly when $u\le m_g(s,t)$. ThereforeFor arbitrary real ${\lambda }_1,\dots ,{\lambda }_n$, the Tonelli theorem now givesThus $(m_g(s_i,s_j){)}_{i,j}$ is a positive semidefinite matrix.

For an orthogonal projection matrix $P=P^{†}=P^2$,Thus $P$ is a positive semidefinite matrix.

No. Consider the two orthogonal projectionsTheir product isFor $x=(-1,2{)}^T$,Thus $PQ$ is not a positive semidefinite matrix. It is also not Hermitian; products of orthogonal projections need not remain orthogonal projections unless the factors commute.

For algebraic stability of a Runge-Kutta method, the weights must be nonnegative andmust be positive semidefinite. With the intended collocation weights,Positive semidefiniteness is therefore possible only at $\alpha =1/3$; substitution gives nonnegative weights and a positive-semidefinite $M$. Hence the intended family is algebraically stable exactly whenWith the sign printed in the paper, one instead obtainsEquality forces $\alpha =1/2$, where $M$ still has nonzero off-diagonal entries and is indefinite. The literal printed tableau is consequently algebraically stable for no value of $\alpha$.

A positive-semidefinite kernel on a nonempty set $\mathcal{X}$ is a symmetric function $k:\mathcal{X}\times \mathcal{X}\to \mathbb{R}$ such that, for every $n\ge 1$, every $x_1,\dots ,x_n\in \mathcal{X}$, and every $c_1,\dots ,c_n\in \mathbb{R}$,Equivalently, every finite kernel matrix $(k(x_r,x_s){)}_{r,s}$ is a positive semidefinite matrix.