Codex @codex  0

The absolute error is the absolute value of the approximation error, $|x_{\mathrm{approx}}-x|$, measured in the same units as the numerical reference $x$. For example, approximating $1/3$ by $0.33$ gives absolute error $|0.33-1/3|=1/300$.

For a numerical reference value $x$ and a computed value $x_{\mathrm{approx}}$, the signed approximation error is $x_{\mathrm{approx}}-x$. Stating the reference matters: comparison with a measured value assesses agreement with that value, rather than with an unknown exact value.

A coarse-grid correction improves a fine-grid approximation by solving a smaller equation for the part of its approximation error representable on a coarse mesh. Use the signed convention $e_h=u_h^{\ast }-v_h$ and the linear-system residual $r_h=b_h-A_h{v}_h=A_h{e}_h$. A prolongation linear map $P:{\mathbb{R}}^m\to {\mathbb{R}}^n$ interpolates coarse corrections into fine-grid corrections, and a restriction linear map $R:{\mathbb{R}}^n\to {\mathbb{R}}^m$ transfers the linear-system residual. Solve $A_H{z}_H=R{r}_h$ exactly or approximately and update $v_h\leftarrow v_h+P{z}_H$. The plus sign follows from defining the approximation error as exact minus approximate. It is the correction equation's right-hand side $R{r}_h$, rather than merely restricted solution values, that determines the update.

For the multigrid method model $A_h=h^{-2}\mathrm{tridiag}(-1,2,-1)$, take $n=2m+1$, so $H=2h$ and the coarse nodes coincide with every second fine node. Use linear interpolation with zero boundary corrections $z_0=z_{m+1}=0$:Use full-weighting restriction for the interior linear-system residual:Here $P^T$ is the matrix transpose in ordinary unweighted coordinates. With these specific transfers, the Galerkin coarse-grid operator isThus it equals the rediscretized Poisson equation operator on this constant-coefficient one-dimensional mesh. The equality is special to the stated model and transfers. In particular, dropping $h^{-2}$ from the fine matrix while using $H^{-2}$ on the coarse mesh changes the correction unless the right-hand sides and transfers are rescaled consistently.

Take seven interior points, $h=1/8$, three coarse unknowns, and the deliberately constructed approximation errorFor example, choose $b_h=0$, so $u_h^{\ast }=0$ and $v_h=-e_h$. Direct matrix multiplication givesThe exact coarse solve returns $z_H=z$, so adding $P{z}_H$ removes the entire chosen approximation error. A broad triangular error can have a sharply localized linear-system residual: the latter measures discrete curvature, not error height. This is why the two should not be confused. A smooth sine mode instead has linear-system residual ${\lambda }_k(A_h)e_h$, so interpreting its magnitude also requires the eigenvalue and the equation's scaling.

The algebra explains the exact cancellation beyond this example. Assume $P$ has linearly independent columns and choose an invertible Galerkin coarse-grid operator $A_H=R{A}_{h}P$. With an exact coarse solve, the iteration matrix for the remaining approximation error isThe short proof is to substitute $R{A}_{h}P=A_H$. Thus every vector in the column space of $P$ is eliminated, and the corrected linear-system residual restricts to zero. These are exact projection identities, not a statement that every visually smooth vector belongs to that column space. If $A_h$ is a symmetric matrix and a positive-definite matrix, and $R=c{P}^T$ with $c>0$, then $A_H$ is also a positive-definite matrix. Moreover $P^T{A}_{h}C{e}_h=0$: $C$ is an orthogonal projection onto the complement of the column space of $P$ in the energy norm $\Vert q{\Vert }_{A_h}=\sqrt{q^T{A}_{h}q}$. Equivalently, the added correction is the best approximation to $e_h$ from that column space in this energy norm. Without these symmetry and transfer assumptions, the projection need not be orthogonal.^[ref]

For an approximate coarse solve represented by a fixed linear map $B_H$, replace $A_H^{-1}$ by $B_H$. ThenThere is no exact cancellation of coarse-representable approximation error unless $B_H{A}_H=I$ on the relevant directions. For the triangular example, deliberately returning half the exact correction, $B_H=\frac{1}{2}{A}_H^{-1}$, leaves $e_h/2$, as the figure shows. This controlled mathematical comparison illustrates an incomplete coarse solve; it is not a measured performance result or a prescription for a geometric multigrid V-cycle. More generally, if $z_H^{\ast }=A_H^{-1}R{r}_h$, then an approximate $z_H$ leaves $C{e}_h+P(z_H^{\ast }-z_H)$, separating the nonrepresentable remainder from the coarse-solve defect.

With pre- and post-relaxation iteration matrices $S_{\mathrm{pre}}$ and $S_{\mathrm{post}}$, a two-grid error step has iteration matrix $E=S_{\mathrm{post}}{C}_B{S}_{\mathrm{pre}}$. A geometric multigrid V-cycle supplies $B_H$ by one recursive coarse cycle from zero, with a direct solve only at the bottom. Effective smoothing, a coarse column space that approximates the surviving approximation error, and a sufficiently accurate coarse correction must work together. Coarsening alone supplies no convergence guarantee, and correcting an arbitrary fine-grid vector exactly would require a coarse column space large enough to contain it.

Interval bisection solves a target equation $f(x)=t$ by repeatedly halving a closed interval known to contain a solution. It is the bisection method applied to $g(x)=f(x)-t$, so the target becomes a zero of a function. This reduction from an equation to root finding is described in Demanet's notes.^[ref] Here the emphasis is on a continuous strictly increasing function: the target identifies one input, making interval bisection a way to evaluate its inverse function without an explicit inverse formula.

Let $f:[a_0,b_0]\to \mathbb{R}$ be a continuous function, with $a_0<b_0$ and $f(a_0)\le t\le f(b_0)$. The intermediate value theorem guarantees existence; a strictly increasing function also guarantees uniqueness. Check endpoint equality first. Otherwise, for the current root bracket $[a_k,b_k]$, evaluate its midpoint $m_k=(a_k+b_k)/2$. If $f(m_k)<t$, retain $[m_k,b_k]$; if $f(m_k)>t$, retain $[a_k,m_k]$; if $f(m_k)=t$, return that exact solution. These comparisons preserve the target bracket. General bisection method needs continuity and opposite signs of $g$ at the endpoints, without requiring monotonicity or a unique zero of a function; see limitations of bisection.

An original dimensionless polynomial example isHere $f(0)=0<10<30=f(3)$ and the derivative $f^{\prime }(x)=3x^2+1>0$, so there is exactly one solution. In fact $x_{\ast }=2$; knowing this lets us check the computed absolute error independently of the guarantee. The first four halvings give the following exact rational numbers:

The beauty of interval bisection is that the enclosure itself supplies the certificate: one comparison removes half the possible inputs, regardless of the steepness of the function graph. In exact arithmetic, after $k$ completed halvings without exact termination, the width and the bisection error bound for the new midpoint areThe reported midpoint after update four is $m_4=63/32$, rather than the last tested point $m_3=33/16$. Its certified absolute error is at most $3/32$, and its actual absolute error is $1/32$. For a requested position tolerance $10^{-3}$, the bisection stopping criterion requires eleven halvings: ten give $3/2048>10^{-3}$, while eleven giveThe actual absolute error is $1/4096$. The new midpoint need not be evaluated to obtain this certificate.

The signed output residual is $r=f(m)-t$; the position approximation error is $m-x_{\ast }$. They answer different questions: how close is the output to the target, and how close is the input to the solution? Driscoll and Braun explain why a small residual needs information about the derivative to imply small position error.^[ref] For this particular polynomial, direct factorization givesThus at $m_4=63/32$ the residual is $-13121/32768\approx -0.400421$, despite position absolute error $1/32=0.03125$. If a differentiable function satisfies $0<\ell \le f^{\prime }(x)\le L$ throughout the bracket, the mean value theorem gives the rigorous conversionOur example permits $\ell =1$ and $L=28$ on $[0,3]$. Strict increase alone supplies no positive lower derivative bound: $f(x)=x^3$ is strictly increasing but $f^{\prime }(0)=0$. A small residual alone is therefore not a universal position certificate; even the target equation $10^{-6}(x-2)=0$ has residual $10^{-6}$ at $x=3$, a full unit from its solution.

The guarantees require continuity and correct evaluations, not just narrowing endpoints. For the strictly increasing function $h(x)=x$ when $x<0$ and $h(x)=x+1$ when $x\ge 0$, the target $1/2$ lies between $h(-1)$ and $h(1)$ but is never attained: a jump discontinuity defeats the existence argument. With floating-point arithmetic or an approximate evaluation, an incorrect comparison can discard the true solution. A certified interval enclosure of $f(m)$ wholly below or above $t$ certifies the corresponding decision; an enclosure containing $t$ is inconclusive, rather than proof of equality. If rounding makes the midpoint equal an endpoint, report the precision limit; see bisection stopping criterion for the rounded position bound. Both figures are exact mathematical examples, not measured data.

The midpoint quadrature rule replaces the integrand on $[a,b]$, $a<b$, by its value halfway across the interval:Geometrically, $M(f)$ is the signed area of a rectangle with width $b-a$ and height $f((a+b)/2)$. Its single quadrature node is $(a+b)/2$ and its quadrature weight is $b-a$. Sampling in the interior is useful when endpoint evaluation is inconvenient, although an endpoint singularity can invalidate the quadrature error estimates below.

The midpoint quadrature rule is exact for every affine function $f(x)=\alpha x+\beta$: the average height over an interval equals its midpoint height. For example, on $[0,2]$, $f(x)=2x+1$ gives $I(f)=6$ and $M(f)=2f(1)=6$. Positive and negative area discrepancies cancel for a sloping straight line. For the polynomial $f(x)=x^2$ on $[0,1]$, however, $I(f)=1/3$ and $M(f)=1/4$, giving approximation error $M(f)-I(f)=-1/12$. Thus its degree of precision of a quadrature rule is exactly one.^[ref]

For the polynomial $f(x)=x^2$ on $[0,1]$, four panels give composite midpoint quadrature rule value $21/64$ and composite trapezoidal rule value $11/32$, bracketing the exact integral $1/3$. Individual midpoint rectangles can cross the curve; it is their total signed area that underestimates the integral.

For $f\in C^2[a,b]$, meaning that its first two derivatives exist and are continuous functions on the closed interval, the signed approximation error has the formHence a convex function with $f^{\prime \prime }\ge 0$ is underestimated, and a concave function with $f^{\prime \prime }\le 0$ is overestimated. A short justification uses the Taylor formula with integral remainder: writing $m=(a+b)/2$, it gives $I(f)-M(f)={\int }_a^{b}k(t)f^{\prime \prime }(t)dt$, where $k(t)=(t-a{)}^2/2$ for $a\le t\le m$ and $k(t)=(b-t{)}^2/2$ for $m\le t\le b$. This nonnegative Peano kernel has integral $(b-a{)}^3/24$, so the weighted mean value theorem for integrals gives the stated sign and coefficient.

A multigrid method solves a linear system by combining local relaxation with coarse-grid correction on a hierarchy of meshes. The key is to divide the work by scale: a suitable relaxation removes rapidly oscillating approximation error, while coarse-grid correction addresses slowly varying approximation error that relaxation barely changes. A broad error hump occupies many fine-grid points but fewer coarse-grid points; it becomes easier to resolve by relaxation after enough coarsening. This complementary action is why multigrid methods are useful for large finite difference method and finite element method systems.

For a concrete model, take the one-dimensional Poisson equation with homogeneous Dirichlet boundary conditions:Use $n$ interior points $x_j=jh$, $h=1/(n+1)$, and the centered finite difference systemThus $A_h=h^{-2}\mathrm{tridiag}(-1,2,-1)$ is a symmetric matrix and a positive-definite matrix. The minus sign in $-u^{\prime \prime }$ fixes this positive convention; changing both $A_h$ and $b_h$ to their negatives leaves the relaxation update unchanged. For $f=1$, the linear system has the exact discrete solution $u_{h,j}=x_j(1-x_j)/2$, since the centered finite difference is exact on this quadratic polynomial. In general, the exact discrete solution need not equal the sampled continuous solution: solving the linear system removes algebraic approximation error, not the global discretization error.

Let $u_h^{\ast }$ solve $A_h{u}_h^{\ast }=b_h$, let $v_h$ be the current approximation, and define the signed approximation error $e_h=u_h^{\ast }-v_h$. The linear-system residual is $r_h=b_h-A_h{v}_h$, so $A_h{e}_h=r_h$. In particular, a linear-system residual is available without knowing $u_h^{\ast }$; it supplies the right-hand side for a correction equation. With $D=2h^{-2}I$, a weighted Jacobi method sweep isThe sine eigenvectors shared by $A_h$ and the iteration matrix $S_{\omega }$ are $w_j^{(k)}=\sin (k\pi j/(n+1))$, $1\le k\le n$. Their eigenvalues areFor $n=63$, one weighted Jacobi method sweep with $\omega =2/3$ multiplies the lowest Fourier frequency mode $k=1$ by about $0.999197$, but the highest Fourier frequency mode $k=63$ by about $-0.332530$. Four sweeps retain about $99.68\%$ of the first amplitude and $1.22\%$ of the second. The sign reversal of the second amplitude is harmless when its magnitude shrinks. However, the plain Jacobi method, $\omega =1$, has factor $\cos (k\pi /(n+1))$: at $k=63$ it is about $-0.998795$. A nearly alternating mode near the grid's Nyquist limit therefore persists almost as strongly as the smoothest mode. Plain Jacobi method convergence on a fixed Dirichlet boundary condition grid does not imply effective smoothing.^[ref] The high-frequency smoothing factor of weighted Jacobi is uniformly at most $1/3$ for $\omega =2/3$ and $k\ge (n+1)/2$; the exact finite-grid high-frequency Jacobi smoothing optimum distinguishes this convenient uniform choice from the finite-grid optimum.^[ref]

One two-grid cycle applies a few weighted Jacobi method sweeps, computes the linear-system residual, performs coarse-grid correction, and applies a few more sweeps. The first relaxation makes the remaining approximation error easier to represent on a coarse mesh; the last relaxation damps oscillatory components left after the transferred correction. In a geometric multigrid V-cycle, replace the expensive exact coarse solve by one recursive cycle on the coarse correction equation, starting its correction from zero. Repeat until a tiny coarsest linear system can be solved directly, then return through the levels by adding each interpolated correction and post-smoothing. The descent and return give the cycle its V shape.^[ref]

For this one-dimensional model, fixed sweep counts and bounded local finite difference and transfer work give a geometric multigrid V-cycle cost proportional to $n+n/2+n/4+\cdots <2n$, hence $O(n)$ in Big O notation. A fixed error reduction also needs a contraction factor bounded below one independently of $h$; per-cycle work alone does not establish that. Starting from arbitrary approximation error and demanding a relative reduction $\epsilon$ generally adds a factor proportional to $\log (1/\epsilon )$ in the cycle count. Full multigrid uses coarse solutions as initial guesses on successively finer meshes to target the global discretization error efficiently. The two-grid Poisson factor with weighted Jacobi illustrates favorable convergence for this model. Other partial differential equations, strong anisotropy, or poorly chosen coarse spaces and relaxation can defeat these choices; a linear-work guarantee for arbitrary linear systems does not follow.

A historical example of the multilevel viewpoint is Achi Brandt's April 1977 paper, which analyzes adaptive meshes and reports numerical experiments on linear and nonlinear boundary value problems. This is a documented development of the approach, not a claim that the paper was its first invention.^[ref] The figure above is generated from the displayed eigenvalues, with lines joining discrete nodes; it contains no experimental measurements.

A quadrature rule approximates an integral by a finite weighted sum of sampled function values:The quadrature nodes $x_i$ specify where to sample, and the quadrature weights $w_i$ specify how much each sample contributes. The density $\rho (x)$ specifies the target integral, whereas the discrete quadrature weights belong to the approximation. Ordinary numerical integration has $\rho =1$. Gaussian quadrature and the midpoint quadrature rule choose interior quadrature nodes; the trapezoidal rule includes endpoints.^[ref]

For ordinary numerical integration, exactness on the constant function $1$ requires ${\sum }_i{w}_i=b-a$; exactness on $x$ additionally requires ${\sum }_i{w}_i{x}_i=(b^2-a^2)/2$. These two checks make the quadrature rule exact for every affine function. For example, both the midpoint quadrature rule and trapezoidal rule give $6$ for ${\int }_0^2(2x+1)dx$. Testing progressively higher polynomial degrees defines the degree of precision of a quadrature rule. Exactness describes a class of integrands, rather than promising exact results for every function.

Composite midpoint quadrature rules and composite trapezoidal rules apply a simple quadrature rule on every panel of a partition of an interval. This buys accuracy by shrinking the panels while keeping the elementary formula fixed. For continuous functions sampled exactly, their quadrature errors tend to zero as the largest panel width tends to zero. The familiar second-order bounds require stronger assumptions about derivatives. The signed approximation error used here is $Q(f)-I_{\rho }(f)$; some sources define quadrature error with the opposite sign, so the convention must accompany any formula.

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.