Adaptive quadrature and error indicators

Let $a<b$, let $n\ge 2$ be even, and let $f$ be real-valued at the nodes of the uniform grids with $n$ and $2n$ subintervals. The coarse Simpson value is $S_0=Q_S(f;a,b,n)$ and the refined value is $S_1=Q_S(f;a,b,2n)$, using Composite Simpson rule on an even uniform grid. The signed refinement indicator is $(S_1-S_0)\mathrm{/}15$ and its magnitude is ${\eta }_n(f;a,b):=\frac{\mathrm{\mid }{S}_1-S_0\mathrm{\mid }}{15}\mathrm{.}$

For the local two-panel case $n=2$, $S_1$ is also the sum of the two one-block Simpson values on the halves of $[a,b]$.

This definition alone is not a bound on the true error.

For data supplied only at the grid nodes, the calls to the Simpson procedure mean calls on the extension $\tilde{f}:[a,b]\to \mathbb{R}$ that agrees with those data at the nodes and is zero elsewhere. This extension is evaluable everywhere as required by the procedure. Its finite sums use only the prescribed nodes, so any other evaluable extension gives the same $S_0$, $S_1$ and indicator; no off-grid values or regularity are required.

The factor $1\mathrm{/}15$ has asymptotic Richardson justification under the hypotheses of Conditional asymptotic validity of the Simpson /15 indicator; it remains a heuristic indicator at any finite grid.

Let $a<b$, let $k\ge 1$ be an integer, set $n_k=2^k$ and $h_k=(b-a)\mathrm{/}{n}_k$, and let $Q_k=Q_S(f;a,b,n_k)$ for $f\in C^4([a,b])$. With $I={\int }_a^{b}f(x)dx$, one has $I-Q_k=-\frac{h_k^4}{180}{\int }_a^b{f}^{(4)}(x)dx+o(h_k^4)(k\to \mathrm{\infty })\mathrm{.}$ Here $r_k=o(h_k^4)$ means $r_k\mathrm{/}{h}_k^4\to 0$ along this sequence. The integral is the Riemann integral.

Input $a<b$, an absolute indicator tolerance $\tau >0$, a remaining depth $D\in \mathbb{N}$ and a positive panel-call budget $P\in \mathbb{N}$, and a deterministic evaluation routine $\mathrm{Eval}$ that returns a finite real value at a requested node or reports failure. It reports failure when a required node cannot be evaluated, including when its midpoint is unrepresentable in the evaluator's arithmetic. Apply the exact-real Simpson formulas of Composite Simpson rule on an even uniform grid to the returned values; floating-point rounding accuracy is not asserted.

Define the recursive call $A(u,v,d,p,t)$ for $u<v$, $d\ge 0$, $p\ge 1$ and $t>0$ as follows.

1. Compute $m=(u+v)\mathrm{/}2$ and the five values of $\mathrm{Eval}$ at $u,(3u+v)\mathrm{/}4,m,(u+3v)\mathrm{/}4,v$. On failure, return failure with no estimate.
2. Compute $S_0=Q_S(f;u,v,2)$ and $S_1=Q_S(f;u,v,4)$ from the two nested grids of Embedded Simpson refinement indicator, and put $\eta ={\eta }_2(f;u,v)=\mathrm{\mid }{S}_1-S_0\mathrm{\mid }\mathrm{/}15$.
3. If $\eta \le t$, return (accepted, estimate $S_1$, indicator $\eta$).
4. If $\eta >t$ and either $d=0$ or $p<3$, return (cap-reached, estimate $S_1$, indicator $\eta$).
5. Otherwise put $p_L=\lfloor (p-1)\mathrm{/}2\rfloor$ and $p_R=\lceil (p-1)\mathrm{/}2\rceil$. Recursively call $A(u,m,d-1,p_L,t\mathrm{/}2)$ and $A(m,v,d-1,p_R,t\mathrm{/}2)$. If either returns failure, return failure. Otherwise return the sum of their estimates and indicators, with status cap-reached if either child has that status and accepted otherwise.

The public call is $A(a,b,D,P,\tau )$. The budget counts recursive panel calls; each call uses at most five point evaluations. accepted means that the returned leaf indicators sum to at most $\tau$. Neither that status nor the indicator is a certified true-error statement. A cap result is an estimate with explicit cap-reached status, not a successful tolerance claim.

Let $a<b$, $f\in C^4([a,b])$, $n_k=2^k$ for $k\ge 1$, $h_k=(b-a)\mathrm{/}{n}_k$, $Q_k=Q_S(f;a,b,n_k)$, and $c=-(1\mathrm{/}180){\int }_a^b{f}^{(4)}(x)dx$. Let $I={\int }_a^{b}f$. Then the refinement indicator ${\eta }_k=\frac{\mathrm{\mid }{Q}_{k+1}-Q_k\mathrm{\mid }}{15}$ satisfies ${\eta }_k=\mathrm{\mid }c\mathrm{\mid }{h}_k^4\mathrm{/}16+o(h_k^4)$. If $c\mathrm{\ne }0$, then ${\eta }_k\mathrm{/}\mathrm{\mid }I-Q_{k+1}\mathrm{\mid }\to 1$. The Richardson-corrected value $R_{k+1}=Q_{k+1}+(Q_{k+1}-Q_k)\mathrm{/}15$ satisfies $I-R_{k+1}=o(h_k^4)$. When $c=0$, these statements give no relative error estimate for ${\eta }_k$; at no finite level does ${\eta }_k$ alone certify the true error.

For valid inputs, Finite-depth adaptive Simpson with an indicator terminates after at most $P$ recursive panel calls and at most $5P$ evaluation requests. A finite result is a sum of the refined Simpson values on the returned leaf intervals, and its reported indicator is the sum of their nonnegative indicators. An evaluation failure, including an unrepresentable required node reported by the evaluator, returns failure with no estimate. If the status is accepted, that indicator sum is at most the requested tolerance $\tau$. If no evaluation fails and a rejected leaf cannot split because of the depth or panel-call cap, the result has status cap-reached. These claims do not bound the true quadrature error.

Input finite dyadic endpoints $a<b$, a positive dyadic tolerance $\tau$, integers $D\ge 0$ and $P\ge 1$, and a function $f\in C^4([a,b])$. Also supply two validated enclosure routines:

• $V(x)$ returns a finite real interval containing the exact value $f(x)$;
• $Y([u,v])$ returns a finite real interval containing every $f^{(4)}(x)$ for $x\in [u,v]$.

If either routine does not meet its enclosure contract, this algorithm makes no certification claim; soundness is a mathematical caller precondition and cannot be detected from samples. An explicit oracle failure or malformed or non-finite interval result returns failure. The intervals and operations are those of Real intervals, complex rectangles and ball enclosures and Outward-rounded interval, rectangle and ball operations. The function $\mathrm{ru}$ is upward directed rounding from Dyadic floating-point format and directed rounding. A midpoint that is not exactly representable or a non-finite interval operation returns failure with no enclosure; malformed or non-finite oracle results are covered by the explicit failure rule above.

Define $A([u,v],d,p,t)$ recursively. Form the exact midpoint $m=(u+v)\mathrm{/}2$; if it is not representable in the format, return failure. Query $V(u),V(m),V(v)$ and let $Y([u,v])=[y_{-},y_{+}]$. Put $M=\max (\mathrm{\mid }{y}_{-}\mathrm{\mid },\mathrm{\mid }{y}_{+}\mathrm{\mid })$ and let $L=v-u$ be the exact difference of the dyadic endpoints. Form an outward-rounded interval $W=[v,v]-[u,u]$ by interval subtraction, and compute $S_P:=(W\mathrm{/}[6,6])(V(u)+4V(m)+V(v))$ using outward interval operations. Compute $B_P=\mathrm{ru}(M{L}^5\mathrm{/}2880)$ by upward rounding the exact real expression; if it overflows, return failure. Set $C_P=S_P+[-B_P,B_P]$ by outward interval addition. If any interval operation has a non-finite status, return failure. Treat the endpoints of $C_P$ as exact dyadic rationals and report their exact rational midpoint $q_P$ and half-width $r_P$; compare $r_P$ with $t$ exactly.

1. If $r_P\le t$, return (certified, enclosure $C_P$, estimate $q_P$, error bound $r_P$).
2. If $r_P>t$ and either $d=0$ or $p<3$, return (inconclusive, enclosure $C_P$, estimate $q_P$, error bound $r_P$).
3. Otherwise split at $m$, allocate $p_L=\lfloor (p-1)\mathrm{/}2\rfloor$ and $p_R=\lceil (p-1)\mathrm{/}2\rceil$, and recurse with depth $d-1$ and tolerance $t\mathrm{/}2$ on each half. A failed child or non-finite interval sum returns failure. Otherwise outward-add the child enclosures to obtain $C$, with exact midpoint $q$ and half-width $r$; return certified if $r\le t$, and inconclusive otherwise.

The public call is $A([a,b],D,P,\tau )$. certified asserts an error at most $r\le \tau$; inconclusive still returns a valid enclosure and bound but does not assert the requested tolerance was met. The derivative enclosure routine is an explicit mathematical precondition, not something inferred from sampled values or the Simpson indicator. Outward rounding of all interval operations, and the separate Peano remainder bound, are essential.

Assume the input hypotheses and the two sound enclosure-oracle contracts in Certified adaptive Simpson from derivative enclosures. Whenever the algorithm returns certified or inconclusive with a finite enclosure $C=[L,U]$, the exact Riemann integral $I={\int }_a^{b}f(x)dx$ belongs to $C$. Its reported center $q=(L+U)\mathrm{/}2$ and half-width $r=(U-L)\mathrm{/}2$ satisfy $\mathrm{\mid }I-q\mathrm{\mid }\le r$. A certified status additionally means $r\le \tau$. The algorithm terminates after at most $P$ recursive panel calls. If an oracle reports failure or returns a malformed or non-finite interval, or an interval operation is non-finite, the algorithm returns failure and asserts no finite enclosure. An inconclusive result does not assert that the requested tolerance was met.