One-step method ($y_{n+1}=y_n+h\Phi (x_n,y_n,h)$)

A one-step method approximates an ordinary differential equation using the current approximate solution and evaluations within the next step. The function $\Phi$ is its increment function. A Runge-Kutta method may use several intermediate stages while remaining a one-step method: it does not require additional previous solution values as a linear multistep method does. Its consistency error is the exact-start update defect divided by $h$.

A Runge-Kutta method advances an initial value problem $y^{\prime }=f(t,y)$, $y(t_0)=y_0$, by sampling its vector field at several trial states within one time step. These samples are Runge-Kutta stages; a weighted combination of their slopes produces the next state. The construction improves on the single tangent line used by the Euler method while requiring evaluations of $f$, rather than explicit higher derivatives of the solution. It applies to scalar equations and systems of ordinary differential equations alike. The coefficient choices distinguish methods such as Heun's method, the classical fourth-order Runge-Kutta method and implicit Runge-Kutta methods.

For $s$ stages, write $A=(a_{ij})$ for the stage matrix, $c_i$ for relative time nodes and $b_i$ for output weights. A time step of length $h$ computesHere $k_i$ is a slope, with units of state per time; $h{k}_i$ is a state increment. Another common convention writes $K_i=h{k}_i$, so its stage states are $Y_i=y_n+{\sum }_j{a}_{ij}{K}_j$ and its output is $y_n+{\sum }_i{b}_i{K}_i$. These are equivalent conventions: including $h$ both inside $K_i$ and again in the output would count the time step twice. The textbook formulation uses the increment convention.^[ref]

Usually the Butcher tableau obeys $c_i={\sum }_j{a}_{ij}$, called internal consistency. Stage $i$ then combines slopes over an effective time $c_{i}h$. The second-order conditions with independent Runge-Kutta abscissae explain what changes when time nodes are specified independently. A Runge-Kutta stage is generally an auxiliary approximation, not a separate accepted solution value at its node. Several stages can even have the same time node but different state arguments, as in the classical fourth-order Runge-Kutta method.

The useful idea is to approximate the integral $y(t_n+h)-y(t_n)={\int }_{t_n}^{t_n+h}f(t,y(t))dt$ without knowing the intermediate exact trajectory. Trial states provide places to sample $f$, and Butcher order conditions make the resulting combination agree with the exact Taylor expansion. More stages alone do not guarantee a higher order of a numerical method: their coefficients must satisfy the appropriate conditions. Embedded Runge-Kutta pairs reuse these same samples for a second output and a practical local error estimator.

A Runge-Kutta stage consists of a trial state $Y_i$ at the time node $t_n+c_{i}h$ and the corresponding vector field evaluation $k_i$. In the slope convention, $Y_i=y_n+h{\sum }_j{a}_{ij}{k}_j$. Stages belong to the construction of one Runge-Kutta method step; they need not have the same accuracy as its final output. Stage states can leave the physically allowed region even when the endpoint lies inside it, so $f$ must be defined at every stage used.

An explicit Runge-Kutta method has $a_{ij}=0$ whenever $j\ge i$. Its stage matrix is strictly lower triangular: compute each Runge-Kutta stage using already known slopes, starting with $Y_1=y_n$. In contrast, an implicit Runge-Kutta method requires solving coupled stage equations, whose existence is addressed by stage solvability of an implicit Runge-Kutta method.

Here Heun's method means the explicit endpoint-averaging Runge-Kutta method: take an Euler predictor $Y_2=y_n+hf(t_n,y_n)$, sample $k_2=f(t_n+h,Y_2)$ there, and average it with $k_1=f(t_n,y_n)$ for the output. The second Runge-Kutta stage is evaluated at the predictor, not at the corrected endpoint. In particular it requires no implicit solve as the ODE trapezoidal rule does. Naming conventions vary between books; the formulas specify the method unambiguously.^[ref]

For smooth $f$, a Taylor expansion about $(t_n,y_n)$ gives $k_2=f+h(f_t+f_{y}f)+O(h^2)$, where $f_y$ is the state Jacobian matrix. Consequently the output is $y_n+hf+h^2(f_t+f_{y}f)/2+O(h^3)$, matching the exact solution through degree two by the chain rule. This gives local truncation error $O(h^3)$ and, with the usual bounded smoothness and Lipschitz continuity hypotheses on a fixed finite interval, global error $O(h^2)$. The Euler-Heun embedded pair compares this output with its Euler predictor.

For $y^{\prime }=y$, $y_n=1$ and $h=0.2$, the sampled slopes are $1$ and $1.2$. Transporting the second slope to the starting point would end at $1.24$; the average slope $1.1$ instead ends at $1.22$. The exact exponential function reaches $e^{0.2}\approx 1.22140275816$. The magnified endpoint comparison separates the pair difference $0.02$ from the actual Heun local defect $0.00140275816$.

An $s$-stage Runge-Kutta method has stability function $R(z)=1+z{b}^T(I-zA{)}^{-1}1$, a rational function with numerator and denominator degrees at most $s$. Order $2s$ forces $R(z)-e^z=O(z^{2s+1})$, so uniqueness of the diagonal Padé approximant identifies $R$ with $[s/s{]}_{e^z}$. That approximant has no poles in the closed left half-plane and modulus at most one there. Hence every such maximal-order method is A-stable. Since the Padé approximant numerator and denominator are coprime and both have degree $s$, the original stage determinant cannot have hidden canceled poles in that half-plane.

This one-step family interpolates between the Forward Euler method at $a=0$, the trapezoidal rule at $a=1/2$, and the Backward Euler method at $a=1$. Its stability function is $R(z)=[1+(1-a)z]/[1-az]$. It is first order except at one half, where it is second order. Applying it to spatial diffusion gives the theta diffusion stability criterion.

On an orthogonal diffusion eigenmode, the theta method has the displayed amplification with $s=-\Delta t\lambda$. The identity $(1+as{)}^2-[1-(1-a)s{]}^2=2s+(2a-1)s^2$ proves contraction at all nonnegative $s$ exactly when $a\ge 1/2$. For $a<1/2$, contraction holds up to $s\le 2/(1-2a)$, with a positive denominator in that range. Orthogonality is what transfers these modal bounds to a full L2 norm estimate.

At $a=1/2$, arbitrarily stiff heat modes approach amplification $-1$: the Crank-Nicolson method is contractive but does not rapidly damp unresolved rough components. The Backward Euler method instead has limit zero. Thus unconditional linear stability and strong stiff decay are different requirements. High smooth-data consistency order alone does not imply a uniform high-order estimate down to the initial instant for arbitrary L2 initial conditions.

A strong stability preserving Runge-Kutta method writes its stages as convex combinations of previously computed states and suitable Forward Euler method steps. It transfers any convex-functional nonincrease property of the forward step, such as a norm bound, under a proportionally scaled time-step restriction.

For example, if $E_k(U)=U+kF(U)$ is nonexpansive in a norm, then $Y=E_k(U)$ and $U_{\mathrm{new}}=(U+E_k(Y))/2$ is also nonexpansive: the triangle inequality bounds the new difference by one half of the initial difference plus one half of the twice-advanced difference. Both are at most the initial difference. A Taylor expansion yields $U_{\mathrm{new}}=U+kF(U)+k^2{F}^{\prime }(U)F(U)/2+O(k^3)$, so this is a second-order Runge-Kutta method. The forward-step hypothesis must hold on the stage states; preservation of one chosen convex bound is not the same as unconditional B-stability.

The classical four-stage explicit Runge-Kutta method usesIt satisfies the fourth-order conditions for a Runge-Kutta method. Its polynomial stability function excludes A-stability, while $|R(iq){|}^2=1-q^6/72+q^8/576$ gives the imaginary-axis stability interval $|q|\le 2\sqrt{2}$.

A Butcher tableau records a Runge-Kutta method's stage matrix $A$, time nodes $c$, and output weights $b$ in the block arrayRead row $i$ above the line as the recipe for Runge-Kutta stage $i$: the left entry gives $t_n+c_{i}h$, and the right entries give the coefficients in $Y_i=y_n+h{\sum }_j{a}_{ij}{k}_j$. Below the line, $b_i$ weights slope $k_i$ in the final update. Zeros on and above the diagonal identify an Explicit Runge-Kutta method. For the internally consistent convention, row sums of $A$ equal $c$; first-order consistency of a numerical method also requires ${\sum }_i{b}_i=1$. An embedded Runge-Kutta pair adds another weight row while retaining $A,c$.

A Butcher order condition equates the coefficient assigned by a Runge-Kutta method to a rooted tree with the corresponding exact Taylor coefficient. Satisfying every condition through tree order $p$ is equivalent to Runge–Kutta order at least $p$.

For a general nonautonomous Runge-Kutta method, Taylor expansion of stage $i$ gives $k_i=f+h[c_i{f}_t+(Ae{)}_i{f}_{y}f]+O(h^2)$. Matching the exact increment $hf+h^2(f_t+f_{y}f)/2$ gives the three displayed conditions. When internal consistency $c=Ae$ is imposed, the last two coincide. Internal consistency is a usual tableau convention, but is not logically necessary for order two of a generalized tableau.

With $c=Ae$ and $C=\mathrm{diag}(c)$, a Runge-Kutta method has order at least four exactly when the eight Butcher order conditionshold; powers of $c$ are componentwise. These conditions match the numerical and exact coefficients indexed by rooted trees through order four. Matching the scalar Dahlquist test equation alone is insufficient to check the nonlinear Butcher order conditions.

An implicit Runge-Kutta method defines one or more stages through equations involving those same stages.

For stages $k_1=f(y)$ and $k_2=f(y+(1-a)h{k}_1+ah{k}_2)$ with update $y+h(k_1+k_2)/2$, the method has order two for every fixed real $a$. The third-order coefficient of $f^{\prime \prime }{f}^2$ is always $1/4$, not $1/6$, so no parameter gives third order on general nonlinear equations. The displayed stability function is bounded throughout the left half-plane only when $a=1/2$, when it is the A-stable trapezoidal function $(1+z/2)/(1-z/2)$. It is not L-stable.

For a vector field uniformly Lipschitz continuous in the state with constant $L$, the stage fixed-point map of an implicit Runge-Kutta method is a contraction in the maximum stage norm if $hL{\max }_i{\sum }_j|a_{ij}|<1$. The contraction mapping theorem then gives unique stages. This is a sufficient small-step condition; a dissipative vector field may permit solvability for larger steps. An algebraic stability or B-stability estimate compares existing stage solutions and should not be mistaken for an unqualified stage-existence theorem.

The implicit midpoint rule advances an autonomous ordinary differential equation byIt is a second-order, A-stable, time-symmetric numerical method.

The implicit midpoint rule is the midpoint method with $y_{n+1}=y_n+hf(t_n+h/2,(y_n+y_{n+1})/2)$. It is a one-stage implicit Runge-Kutta method.

For a uniformly Lipschitz continuous state field $f$ with constant $L$, an implicit midpoint rule step is uniquely solvable when $hL/2<1$, by the contraction mapping theorem. Write its step map as ${\Psi }_h(u)=u+h{\Phi }_h(u)$, where ${\Phi }_h(u)=f(t+h/2,u+h{\Phi }_h(u)/2)$. Subtracting the increment equations givesTherefore $\mathrm{Lip}({\Phi }_h)\le L/(1-hL/2)$, and the triangle inequality gives the bound for ${\Psi }_h$. Combining it with the implicit midpoint consistency bound yieldsfor $L>0$. This follows by summing the geometric error recurrence with forcing $K{h}^3/(1-hL/2)$. For $hL\le 1$ and $nh\le T$, $R_h^n\le e^{2LT}$, proving second-order global error with exact starting data. If $L=0$, the direct bound is $\Vert e_n\Vert \le \Vert e_0\Vert +KT{h}^2$.

For an exact solution $y^{\prime }=f(t,y)$, define the normalized local truncation error of the implicit midpoint rule byIf $\Vert y^{\prime \prime }\Vert \le M_2$, $\Vert y^{\prime \prime \prime }\Vert \le M_3$ and the state Jacobian matrix of $f$ has operator norm at most $L$ on the intervening state segments, then the displayed bound holds uniformly on a fixed finite interval. To prove it, put $t=t_n+h/2$. Symmetric Taylor expansions giveFor the first bound, integrate the second-order remainder of $y^{\prime }(t+s)$ over the symmetric step; for the second, average the second-order remainders of $y$ at the endpoints. The odd terms cancel. The Lipschitz condition for $f$ and $y^{\prime }(t)=f(t,y(t))$ finish the proof by the triangle inequality. The unnormalized step residual is $h{\tau }_n=O(h^3)$.

The Crank--Nicolson method applies the trapezoidal rule to a semidiscrete evolution equation. For a linear system $u^{\prime }=Lu$ its update is

For symmetric negative semidefinite matrices $A,B$, each factor $C_A(h)=(I-hA/2{)}^{-1}(I+hA/2)$ has norm at most one, since its eigenvalues are $(1+h\lambda /2)/(1-h\lambda /2)$ with $\lambda \le 0$. The same holds for $B$, so their ordered product is a contraction for every $h\ge 0$, whether or not $A,B$ commute. This proves unconditional discrete L2 norm stability of the split Crank-Nicolson method for directional diffusion. Its splitting error remains first order in general, because the product's quadratic defect is $h^2[A,B]/2$.

Given distinct nodes $c_i$, a collocation Runge-Kutta method uses the Lagrange basis ${\ell }_j$ with $a_{ij}={\int }_0^{c_i}{\ell }_j(\tau )d\tau$ and $b_j={\int }_0^1{\ell }_j(\tau )d\tau$. Its stages enforce the differential equation at the collocation nodes.

This applies the collocation method to the time-step polynomial of an ordinary differential equation.

For distinct nodes $c_1,c_2\in [0,1]$, the two-stage collocation Runge-Kutta method has stability functionThe denominator has no zeros in the closed left half-plane. On the imaginary axis the denominator's squared modulus minus the numerator's is $(d^2-r^2)y^4$. Thus $d\ge r$, equivalently $c_1+c_2\ge 1$, is necessary and sufficient for A-stability, using the maximum modulus principle and the bound at infinity.

For $s$ distinct collocation nodes, let ${\pi }_s(t)={\prod }_{i=1}^s(t-c_i)$. If ${\int }_0^1{t}^j{\pi }_s(t)dt=0$ for $0\le j<k$, and the next moment is nonzero, the associated interpolatory quadrature has order $s+k$, and the collocation Runge-Kutta method has the same classical order for sufficiently smooth differential equations. Here $0\le k\le s$. Equivalently, if the quadrature is exact through degree $p-1$ but not degree $p$, the method has order $p$. The generic order is at least $s$, Gauss collocation methods reach $2s$, and Radau IIA methods reach $2s-1$.

For $s$ collocation nodes, the residual of the collocation polynomial vanishes at each node and contains the factor $\omega (\tau )={\prod }_i(\tau -c_i)$. Quadrature exactness through degree $q-1$ implies ${\int }_0^1\omega (\tau )p(\tau )d\tau =0$ for every polynomial of degree at most $q-s-1$. Expand the smooth vector-field interpolation error and the derivative of the exact flow in powers of the step. Every endpoint-error term of degree below $q$ has this vanishing weighted moment, so the endpoint local error is $O(h^{q+1})$, although the interior stage accuracy is generally only $s$. A step Lipschitz estimate and discrete Gronwall inequality then give global order $q$. Gaussian nodes attain $q=2s$, while the first nonexact quadrature moment bounds the order from above.

For a collocation Runge-Kutta method, $a_{ij}={\int }_0^{c_i}{\ell }_j(t)dt$, where the Lagrange interpolation basis obeys ${\sum }_j{\ell }_j=1$. Hence each row sum is its collocation node. Failure of this identity disproves a claimed collocation tableau with those nodes. The identity alone is necessary rather than sufficient: the individual integrated basis coefficients must also agree.

A Lobatto IIIA method collocates at the Lobatto quadrature nodes, including both endpoints of the time interval. Its three-stage version hasThe Butcher order conditions give order four. Its stability function iswhich is A-stable but not L-stable. It is not algebraically stable: the first diagonal entry of the defining matrix is $2b_1{a}_{11}-b_1^2=-1/36$. This distinguishes linear A-stability from nonlinear algebraic stability.

The $s$-stage Gauss--Legendre Runge--Kutta method collocates at the Gauss--Legendre quadrature nodes on $[0,1]$. It is symmetric, A-stable, and has order $2s$.

Choose $c_1,\dots ,c_s$ as the polynomial roots of $P_s(2c-1)$ in $(0,1)$, and form the Lagrange interpolation polynomials ${\ell }_j(t)={\prod }_{m\ne j}(t-c_m)/(c_j-c_m)$. The Gauss collocation method has $a_{ij}={\int }_0^{c_i}{\ell }_j(t)dt$ and $b_j={\int }_0^1{\ell }_j(t)dt$. This gives explicit coefficient formulas after finding polynomial roots: if ${\ell }_j(t)={\sum }_r{d}_{jr}{t}^r$, these integrals are ${\sum }_r{d}_{jr}{c}_i^{r+1}/(r+1)$ and ${\sum }_r{d}_{jr}/(r+1)$. Gaussian quadrature gives order $2s$.

The two-stage Gauss--Legendre Runge-Kutta method uses nodes $c_{1,2}=1/2\mp \sqrt{3}/6$, weights $b_1=b_2=1/2$ and coefficients $a_{11}=a_{22}=1/4$, $a_{12}=1/4-\sqrt{3}/6$, $a_{21}=1/4+\sqrt{3}/6$. The Butcher order conditions give order four. Its stability function is the displayed rational function, which is A-stable. The algebraic stability matrix is zero, so positive weights also give nonlinear contractivity for dissipative vector fields.

An $s$-stage Radau IIA method collocates at the right-endpoint Radau nodes. It has order $2s-1$ and is both algebraically stable and L-stable.

Suppose an $s$-stage Runge-Kutta method has order $2s-1$, an invertible matrix $A$, and $b^T{A}^{-1}1=1$. Its stability function is $R(z)=1+z{b}^T(I-zA{)}^{-1}1$. The adjugate matrix formula gives a rational function with denominator degree $s$ and numerator degree at most $s$; its zero limit at infinity lowers the latter degree to at most $s-1$. The order of a numerical method gives $R(z)-e^z=O(z^{2s})$, identifying the $[s-1/s]$ Padé approximant. Uniqueness follows because the cross-multiplied difference of two candidates has degree at most $2s-1$ but a zero of order at least $2s$, so is identically zero. The A-stability of near-diagonal exponential Padé approximants then proves A-stability, and the zero limit proves L-stability. These conditions identify the stability function; they do not alone claim uniqueness of all entries of the Runge-Kutta method tableau.

A Runge-Kutta method is algebraically stable when $b_i\ge 0$ and the matrix $M$ with $M_{ij}=b_i{a}_{ij}+b_j{a}_{ji}-b_i{b}_j$ is positive semidefinite. Algebraic stability implies contractivity for dissipative differential equations.

For stage states $Y_i=y_n+h{\sum }_j{a}_{ij}{k}_j$ and a vector field with $Y_i^T{k}_i=0$, expansion of the new squared norm gives $\Vert y_{n+1}{\Vert }^2-\Vert y_n{\Vert }^2=-h^2{\sum }_{i,j}{M}_{ij}{k}_i^T{k}_j$, where $M_{ij}=b_i{a}_{ij}+b_j{a}_{ji}-b_i{b}_j$. Hence $M=0$ preserves the quadratic form exactly for every well-defined step. The same argument with a fixed symmetric matrix in the inner product preserves any quadratic invariant of the vector field. In particular it preserves the Euclidean norm for a skew-symmetric matrix differential equation.

An algebraically stable Runge-Kutta method is B-stable: for a dissipative vector field, a well-defined step cannot increase the distance between two starting values. This assumes that the compared stage equations have solutions, as required by stage solvability of an implicit Runge-Kutta method. The Runge-Kutta contractivity identity expresses the distance change as a nonpositive weighted term coming from the dissipative vector field minus a nonnegative quadratic form involving the matrix defining algebraic stability. Positive weights also rule out stage poles in the open left half-plane on the scalar test equation, giving A-stability there; zero-weight redundant stages need separate solvability attention.

For differences $D_i$ between corresponding stages and differences $F_i$ between their vector fields, a Runge--Kutta step satisfiesTo derive the identity, expand the squared norm of $d_{n+1}=d_n+h{\sum }_i{b}_i{F}_i$ and substitute $d_n=D_i-h{\sum }_j{a}_{ij}{F}_j$ in the linear term. Symmetrizing the double sum gives the displayed coefficients $m_{ij}$. For a dissipative vector field, the first sum is nonpositive when $b_i\ge 0$, and the second sum is nonnegative when the matrix $(m_{ij})$ is a positive semidefinite matrix. Thus algebraic stability implies B-stability.

As a quadrature rule, the trapezoidal rule replaces an integrand by the linear interpolation joining its two endpoint values:For nonnegative endpoint values, this is the area of a trapezoid; in general it is a signed area. The quadrature nodes are $a,b$, each with quadrature weight $(b-a)/2$. Unlike the midpoint quadrature rule, it uses endpoint samples, which can be reused between adjacent panels.

Every affine function $f(x)=\alpha x+\beta$ agrees with its linear interpolation, so the trapezoidal rule integrates it exactly. On $[0,2]$, $f(x)=2x+1$ gives $T(f)=(2/2)(1+5)=6=I(f)$. On $[0,1]$, the polynomial $f(x)=x^2$ instead gives $T(f)=1/2$, exceeding $I(f)=1/3$ by $1/6$. Its degree of precision of a quadrature rule is exactly one.

If $f\in C^2[a,b]$, its signed approximation error satisfiesA convex function lies below its endpoint chord, so the trapezoidal rule overestimates its integral; a concave function gives the reverse inequality. The error formula assumes that the second derivative exists and is a continuous function throughout the closed interval.^[ref] A short derivation is the identityobtained by two integrations by parts. The nonnegative Peano kernel has integral $(b-a{)}^3/12$, and the weighted mean value theorem for integrals gives the formula.

The trapezoidal rule for an ordinary differential equation averages vector fields at the old and new states and is A-stable.

This time-stepping use applies the same trapezoidal rule to the integral of the derivative in an ordinary differential equation. That connection explains its place under implicit Runge-Kutta methods; its quadrature behavior and the nonlinear issue in trapezoidal rule fails B-stability concern different uses of the same endpoint averaging.

For $f\in C^2[a,b]$ with $a<b$, the trapezoidal rule satisfies the displayed bound, and the constant $1/12$ is best possible. The polynomial interpolation error for the endpoint linear interpolant satisfies $|f(x)-p(x)|\le \frac{1}{2}M(x-a)(b-x)$, where $M=\max |f^{\prime \prime }|$. Integrating this inequality and using the triangle inequality for an integral proves the bound because ${\int }_a^b(x-a)(b-x)dx=(b-a{)}^3/6$. Equality holds for $f(x)=(x-a{)}^2$: the exact integral is $(b-a{)}^3/3$, its trapezoidal rule estimate is $(b-a{)}^3/2$, and $M=2$.

On an equal-panel partition of an interval, with $h=(b-a)/n$, sum the trapezoidal rule over adjacent panels:There are $n+1$ distinct quadrature nodes. Each interior sample belongs to two panels and therefore has quadrature weight $h$, while each endpoint has quadrature weight $h/2$. Equivalently, the composite trapezoidal rule integrates the piecewise linear function obtained by linear interpolation exactly.^[ref]

For $f\in C^2[a,b]$ and $K={\max }_{x\in [a,b]}|f^{\prime \prime }(x)|$, summing panel quadrature errors givesfor some $\xi \in (a,b)$. Thus panel quadrature errors of order $O(h^3)$ accumulate to order $O(h^2)$, using big O notation.^[ref] On a nonuniform partition of an interval with widths $h_j$ and largest width $H$, the bound is $K{\sum }_j{h}_j^3/12\le K(b-a)H^2/12$.

For the polynomial $f(x)=x^2$ on $[0,1]$, direct summation yieldsBoth quadrature errors shrink by exactly four when $n$ doubles, with opposite signs. This exact relation belongs to this polynomial example; the trapezoidal rule is not always twice as inaccurate as the midpoint quadrature rule.

The left plot uses exact polynomial sums; the right adds a common sample bias $\delta =0.002$ to the composite trapezoidal rule. That mathematical model converges to $1/3+\delta$, showing why smaller panels alone do not ensure accuracy for noisy samples. It illustrates noise sensitivity of a quadrature rule, rather than measured data.

On a partition of an interval symmetric about zero, the composite trapezoidal rule gives opposite quadrature nodes equal quadrature weights. Samples of an odd function cancel pairwise, and a zero quadrature node contributes zero. For a continuous odd function, its exact integral over $[-L,L]$ is also zero. Adding any affine function preserves exactness by linearity of the integral and of the quadrature rule. For example, $1+\sin x$ is nonpolynomial but both the exact integral and every uniform-grid composite trapezoidal rule value on $[-1,1]$ are two.

On a uniform partition of $[a,b]$ with width $h$, the right Riemann sum and the composite trapezoidal rule have equal weights at every interior mesh point. They differ only by adding half a final-endpoint weight and removing half an initial-endpoint weight. This proves the displayed identity. For continuous $f$, both converge to its definite integral; their finite endpoint correction need not be confused with the smaller trapezoidal discretization error for a smooth integrand.

If $f$ is twice continuously differentiable on $[a,b]$, take a uniform bound $|f^{\prime \prime }|\le M$. On one panel $[u,v]$, two integrations by parts give the signed quadrature errorThe nonnegative kernel has integral $(v-u{)}^3/6$, so the absolute error is at most $M(v-u{)}^3/12$. Summing $n$ equal panels proves the displayed bound. The maximum exists on the closed interval; a supremum or any known uniform bound also suffices. The estimate certifies a sufficient panel count, which need not be the least count achieving a particular actual error.

The trapezoidal rule is A-stable but is not B-stable. For the dissipative vector field $f(y)=-y^3$ at step size one, its uniquely defined step map $T$ obeys $T+T^3/2=y-y^3/2$. At $y=\sqrt{2}$, $T=0$ and implicit differentiation gives $T^{\prime }=-2$, so nearby states expand rather than contract. This is a concrete nonlinear obstruction, rather than merely failure of the sufficient algebraic stability criterion.

A one-step numerical method with step map ${\Phi }_h$ is time symmetric when ${\Phi }_{-h}={\Phi }_h^{-1}$. Its leading local-error power is odd, so every nontrivial consistent time-symmetric method has even order.

The order is the highest power through which the one-step expansion matches the exact Taylor expansion.

For an explicit Runge-Kutta method with coefficients $A,b$, put $c=A1$ and let $c^2$ mean componentwise squaring. The conditions $b^{T}1=1$, $b^{T}c=1/2$, $b^T{c}^2=1/3$ and $b^{T}Ac=1/6$ match the exact flow expansion through degree three: the two cubic differentials are $f^{\prime \prime }[f,f]$ and $(f^{\prime }{)}^{2}f$. For sufficiently smooth vector fields these give local error $O(h^4)$ and, under the usual finite-interval Lipschitz continuity assumptions, global error $O(h^3)$.