Trapezoidal rule

A 2016 Science paper reports that the trapezoid rule was in use in Babylon before 50 BCE for integrating the velocity of Jupiter along the ecliptic.{{cite journal |last=Ossendrijver |first=Mathieu |date=Jan 29, 2016 |title=Ancient Babylonian astronomers calculated Jupiter's position from the area under a time-velocity graph |journal=Science |doi=10.1126/science.aad8085 |pmid=26823423 |volume=351 |issue=6272 |pages=482–484 |bibcode=2016Sci...351..482O |s2cid=206644971 |url=https://www.science.org/doi/full/10.1126/science.aad8085}}{{Cite news |date=2016-01-29 |title=Ancient Babylonians 'first to use geometry' |url=https://www.bbc.com/news/science-environment-35431974 |access-date=2025-02-13 |work=BBC News |language=en-GB}}

When the grid spacing is non-uniform, one can use the formula

$$\backslash int\_\{a\}^\{b\}\; f(x)\backslash ,\; dx\; \backslash approx\; \backslash sum\_\{k=1\}^N\; \backslash frac\{f(x\_\{k-1\})\; +\; f(x\_k)\}\{2\}\; \backslash Delta\; x\_k\; ,$$

wherein $\backslash Delta\; x\_k\; =\; x\_\{k\}\; -\; x\_\{k-1\}\; .$

For a domain partitioned by $N$ equally spaced points, considerable simplification may occur.

Let

$\backslash Delta\; x\; =\; \backslash frac\{b-a\}\{N\},$ and $x\_k=a+k\; \backslash Delta\; x$ for {{tmath|1=k=0,1,\ldots, N}}.

The approximation to the integral becomes

$$\backslash begin\{align\}$$

\int_{a}^{b} f(x)\, dx

\approx \frac{\Delta x}{2}& \sum_{k=1}^{N} \left( f(x_{k-1}) + f(x_{k}) \right) \\[1ex]

&= \Delta x \left( \frac{f(x_N) + f(x_0) }{2} + \sum_{k=1}^{N-1} f(x_k) \right) .

\end{align}

File:Trapezium2.gif

The error of the composite trapezoidal rule is the difference between the value of the integral and the numerical result:

$$\backslash text\{E\}\; =\; \backslash int\_a^b\; f(x)\backslash ,dx\; -\; \backslash frac\{b-a\}\{N\}\; \backslash left[\; \{f(a)\; +\; f(b)\; \backslash over\; 2\}\; +\; \backslash sum\_\{k=1\}^\{N-1\}\; f\; \backslash left(\; a+k\; \backslash frac\{b-a\}\{N\}\; \backslash right)\; \backslash right]$$

There exists a number ξ between a and b, such that{{harvtxt|Atkinson|1989|loc=equation (5.1.7)}}

$$\backslash text\{E\}\; =\; -\backslash frac\{(b-a)^3\}\{12N^2\}\; f\text{'}\text{'}(\backslash xi)$$

It follows that if the integrand is concave up (and thus has a positive second derivative), then the error is negative and the trapezoidal rule overestimates the true value. This can also be seen from the geometric picture: the trapezoids include all of the area under the curve and extend over it. Similarly, a concave-down function yields an underestimate because area is unaccounted for under the curve, but none is counted above. If the interval of the integral being approximated includes an inflection point, the sign of the error is harder to identify.

An asymptotic error estimate for N → ∞ is given by

$$\backslash text\{E\}\; =\; -\backslash frac\{(b-a)^2\}\{12N^2\}\; \backslash big[\; f\text{'}(b)-f\text{'}(a)\; \backslash big]\; +\; O(N^\{-3\}).$$

Further terms in this error estimate are given by the Euler–Maclaurin summation formula.

Several techniques can be used to analyze the error, including:{{Harv|Weideman|2002|loc=p. 23, section 2}}

1. Fourier series
2. Residue calculus
3. Euler–Maclaurin summation formula{{harvtxt|Atkinson|1989|loc=equation (5.1.9)}}{{harvtxt|Atkinson|1989|loc=p. 285}}
4. Polynomial interpolation{{harvtxt|Burden|Faires|2011|p=194}}

It is argued that the speed of convergence of the trapezoidal rule reflects and can be used as a definition of classes of smoothness of the functions.

First suppose that $h=\backslash frac\{b-a\}\{N\}$ and $a\_k=a+(k-1)h$. Let $$g\_k(t)\; =\; \backslash frac\{1\}\{2\}\; t[f(a\_k)+f(a\_k+t)]\; -\; \backslash int\_\{a\_k\}^\{a\_k+t\}\; f(x)\; \backslash ,\; dx$$ be the function such that $|g\_k(h)|$ is the error of the trapezoidal rule on one of the intervals, $[a\_k,\; a\_k+h]$. Then

$$\{dg\_k\; \backslash over\; dt\}=\{1\; \backslash over\; 2\}[f(a\_k)+f(a\_k+t)]+\{1\backslash over2\}t\backslash cdot\; f\text{'}(a\_k+t)-f(a\_k+t),$$

and

$$\{d^2g\_k\; \backslash over\; dt^2\}=\{1\backslash over\; 2\}t\backslash cdot\; f\text{'}\text{'}(a\_k+t).$$

Now suppose that $\backslash left|\; f$(x) \right| \leq \left| f(\xi) \right|, which holds if $f$ is sufficiently smooth. It then follows that

$$\backslash left|\; f$$(a_k+t) \right| \leq f(\xi)

which is equivalent to

$-f$(\xi) \leq f(a_k+t) \leq f(\xi), or $-\backslash frac\{f$(\xi)t}{2} \leq g_k(t) \leq \frac{f(\xi)t}{2}.

Since $g\_k\text{'}(0)=0$ and $g\_k(0)=0$,

$$\backslash int\_0^t\; g\_k\text{'}\text{'}(x)\; dx\; =\; g\_k\text{'}(t)$$ and $$\backslash int\_0^t\; g\_k\text{'}(x)\; dx\; =\; g\_k(t).$$

Using these results, we find

$$-\backslash frac\{f$$(\xi)t^2}{4} \leq g_k'(t) \leq \frac{f(\xi)t^2}{4}

and

$$-\backslash frac\{f$$(\xi)t^3}{12} \leq g_k(t) \leq \frac{f(\xi)t^3}{12}

Letting $t\; =\; h$ we find

$$-\backslash frac\{f$$(\xi)h^3}{12} \leq g_k(h) \leq \frac{f(\xi)h^3}{12}.

Summing all of the local error terms we find

$$\backslash sum\_\{k=1\}^\{N\}\; g\_k(h)\; =\; \backslash frac\{b-a\}\{N\}\; \backslash left[\; \{f(a)\; +\; f(b)\; \backslash over\; 2\}\; +\; \backslash sum\_\{k=1\}^\{N-1\}\; f\; \backslash left(\; a+k\; \backslash frac\{b-a\}\{N\}\; \backslash right)\; \backslash right]\; -\; \backslash int\_a^b\; f(x)dx.$$

But we also have

$$-\; \backslash sum\_\{k=1\}^N\; \backslash frac\{f$$(\xi)h^3}{12} \leq \sum_{k=1}^N g_k(h) \leq \sum_{k=1}^N \frac{f(\xi)h^3}{12}

and

$$\backslash sum\_\{k=1\}^N\; \backslash frac\{f$$(\xi)h^3}{12}=\frac{f(\xi)h^3N}{12},

so that

$$-\backslash frac\{f$$(\xi)h^3N}{12} \leq \frac{b-a}{N} \left[ {f(a) + f(b) \over 2} + \sum_{k=1}^{N-1} f \left( a+k \frac{b-a}{N} \right) \right]-\int_a^bf(x)dx \leq \frac{f(\xi)h^3N}{12}.

Therefore the total error is bounded by

$$\backslash text\{error\}\; =\; \backslash int\_a^b\; f(x)\backslash ,dx\; -\; \backslash frac\{b-a\}\{N\}\; \backslash left[\; \{f(a)\; +\; f(b)\; \backslash over\; 2\}\; +\; \backslash sum\_\{k=1\}^\{N-1\}\; f\; \backslash left(\; a+k\; \backslash frac\{b-a\}\{N\}\; \backslash right)\; \backslash right]\; =\; \backslash frac\{f$$(\xi)h^3N}{12}=\frac{f(\xi)(b-a)^3}{12N^2}.

The trapezoidal rule converges rapidly for periodic functions. This is an easy consequence of the Euler-Maclaurin summation formula, which says that

if $f$ is $p$ times continuously differentiable with period $T$

$$\backslash sum\_\{k=0\}^\{N-1\}\; f(kh)h\; =$$

\int_0^T f(x)\,dx +

\sum_{k=1}^{\lfloor p/2\rfloor} \frac{B_{2k}}{(2k)!} (f^{(2k - 1)}(T) - f^{(2k - 1)}(0)) - (-1)^p h^p \int_0^T\tilde{B}_{p}(x/T)f^{(p)}(x) \, dx

where $h:=T/N$ and $\backslash tilde\{B\}\_\{p\}$ is the periodic extension of the $p$th Bernoulli polynomial.{{cite book|title=Numerical Analysis, volume 181 of Graduate Texts in Mathematics|first=Rainer|last=Kress |year=1998 |publisher=Springer-Verlag}} Due to the periodicity, the derivatives at the endpoint cancel and we see that the error is $O(h^p)$.

A similar effect is available for peak-like functions, such as Gaussian, Exponentially modified Gaussian and other functions with derivatives at integration limits that can be neglected.{{Cite journal |last=Goodwin|first=E. T. |date=1949 |title=The evaluation of integrals of the form $\backslash textstyle\backslash int\_\{-\backslash infty\}^\backslash infty\{f(x)e^\{-x^2\}dx\}$ | journal=Mathematical Proceedings of the Cambridge Philosophical Society |language=en |volume=45 |issue=2 |pages=241–245 |doi=10.1017/S0305004100024786 |bibcode=1949PCPS...45..241G |issn=1469-8064}} The evaluation of the full integral of a Gaussian function by trapezoidal rule with 1% accuracy can be made using just 4 points.{{Cite journal| last1=Kalambet|first1=Yuri |last2=Kozmin|first2=Yuri |last3=Samokhin|first3=Andrey |date=2018 |title=Comparison of integration rules in the case of very narrow chromatographic peaks |journal=Chemometrics and Intelligent Laboratory Systems|volume=179 |pages=22–30 |doi=10.1016/j.chemolab.2018.06.001|issn=0169-7439}} Simpson's rule requires 1.8 times more points to achieve the same accuracy.

For functions that are not in C², the error bound given above is not applicable. Still, error bounds for such rough functions can be derived, which typically show a slower convergence with the number of function evaluations $N$ than the $O(N^\{-2\})$ behaviour given above. Interestingly, in this case the trapezoidal rule often has sharper bounds than Simpson's rule for the same number of function evaluations.

The trapezoidal rule is one of a family of formulas for numerical integration called Newton–Cotes formulas, of which the midpoint rule is similar to the trapezoid rule. Simpson's rule is another member of the same family, and in general has faster convergence than the trapezoidal rule for functions which are twice continuously differentiable, though not in all specific cases. However, for various classes of rougher functions (ones with weaker smoothness conditions), the trapezoidal rule has faster convergence in general than Simpson's rule.{{Harv|Cruz-Uribe|Neugebauer|2002}}

Moreover, the trapezoidal rule tends to become extremely accurate when periodic functions are integrated over their periods, which can be analyzed in various ways.{{Harv|Rahman|Schmeisser|1990}}{{Harv|Weideman|2002}} A similar effect is available for peak functions.

For non-periodic functions, however, methods with unequally spaced points such as Gaussian quadrature and Clenshaw–Curtis quadrature are generally far more accurate; Clenshaw–Curtis quadrature can be viewed as a change of variables to express arbitrary integrals in terms of periodic integrals, at which point the trapezoidal rule can be applied accurately.

The following integral is given:

$$\backslash int\_\{0.1\}^\{1.3\}\{5xe^\{-\; 2x\}\{dx\}\}$$

{{ordered list| list-style-type = lower-alpha

| Use the composite trapezoidal rule to estimate the value of this integral. Use three segments.

| Find the true error $E\_\{t\}$ for part (a).

| Find the absolute relative true error $\backslash left|\; \backslash varepsilon\_\{t\}\; \backslash right|$ for part (a).

}}

Solution

{{ordered list| list-style-type = lower-alpha

| The solution using the composite trapezoidal rule with 3 segments is applied as follows.

$$\backslash int\_\{a\}^\{b\}\{f(x)\{dx\}\}\; \backslash approx\; \backslash frac\{b\; -\; a\}\{2n\}\backslash left\backslash lbrack\; f(a)\; +\; 2\backslash sum\_\{i\; =\; 1\}^\{n\; -\; 1\}\{f(a\; +\; \{ih\})\}\; +\; f(b)\; \backslash right\backslash rbrack$$

$$\backslash begin\{align\}$$

n &= 3 \\

a &= 0.1 \\

b &= 1.3 \\

h &= \frac{b - a}{n}

= \frac{1.3 - 0.1}{3}

= 0.4

\end{align}

Using the composite trapezoidal rule formula

$$\backslash begin\{align\}\; \backslash int\_a^b\; \{f(x)\{dx\}\}\; \backslash approx\; \backslash frac\{b\; -\; a\}\{2n\}\; \backslash left\backslash lbrack\; f(a)\; +\; 2\backslash left\backslash \{\; \backslash sum\_\{i\; =\; 1\}^\{n\; -\; 1\}\{f(a\; +\; \{ih\})\}\; \backslash right\backslash \}\; +\; f(b)\; \backslash right\backslash rbrack\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\backslash ;\; (3)\; \backslash end\{align\}$$

$$\backslash begin\{align\}$$

I &\approx \frac{1.3 - 0.1}{6}\left\lbrack f(0.1) + 2\sum_{i = 1}^{3 - 1}{f(0.1 + 0.4i)} + f(1.3) \right\rbrack\\

I &\approx \frac{1.3 - 0.1}{6}\left\lbrack f(0.1) + 2\sum_{i = 1}^{2}{f(0.1 + 0.4i)} + f(1.3) \right\rbrack\\

&= 0.2\lbrack f(0.1) + 2f(0.5) + 2f(0.9) + f(1.3)\rbrack\\

&= 0.2[5 \times 0.1 \times e^{- 2(0.1)}+2(5 \times 0.5 \times e^{- 2(0.5)})+2(5 \times 0.9 \times e^{- 2(0.9)}) + 5 \times 1.3 \times e^{- 2(1.3)}\rbrack\\

&= 0.84385

\end{align}

| The exact value of the above integral can be found by integration by parts and is

$$\backslash int\_\{0.1\}^\{1.3\}\; 5xe^\{-\; 2x\}\{dx\}\; =\; 0.89387$$

So the true error is

$$\backslash begin\{align\}$$

E_{t} &= \text{True Value} - \text{Approximate Value}\\

&= 0.89387 - 0.84385\\

&= 0.05002

\end{align}

| The absolute relative true error is

&= \left| \frac{0.05002}{0.89387} \right| \times 100\%\\

&= 5.5959\%

\end{align}

}}

• Gaussian quadrature
• Newton–Cotes formulas
• Rectangle method
• Romberg's method
• Simpson's rule
• Tai's model
• {{slink|Volterra integral equation#Numerical Solution using Trapezoidal Rule}}

{{notelist}}

{{refbegin}}

• {{citation |last=Atkinson |first=Kendall E. |year=1989 |title=An Introduction to Numerical Analysis |edition=2nd |publisher=John Wiley & Sons |location=New York |isbn=978-0-471-50023-0}}
• {{citation |last1=Rahman |first1=Qazi I. |last2=Schmeisser |first2=Gerhard |date=December 1990 |title=Characterization of the speed of convergence of the trapezoidal rule |journal=Numerische Mathematik |issn=0945-3245 |doi=10.1007/BF01386402 |volume=57 |issue=1 |pages=123–138|s2cid=122245944 }}
• {{citation |last1=Burden |first1=Richard L. |last2=Faires |first2=J. Douglas |year=2011|title=Numerical Analysis |edition=9th |publisher=Brooks/Cole }}
• {{citation |last=Weideman |first=J. A. C. |date=January 2002 |title=Numerical Integration of Periodic Functions: A Few Examples |journal=The American Mathematical Monthly |doi=10.2307/2695765 |jstor=2695765 |volume=109 |issue=1 |pages=21–36}}
• {{citation |last1=Cruz-Uribe |first1=D. |last2=Neugebauer |first2=C. J. |year=2002 |title=Sharp Error Bounds for the Trapezoidal Rule and Simpson's Rule |journal=Journal of Inequalities in Pure and Applied Mathematics |volume=3 |issue=4 |url=http://www.emis.de/journals/JIPAM/images/031_02_JIPAM/031_02.pdf }}

{{refend}}

{{wikibooks|A-level Mathematics|C2/Integration#Trapezium Rule|Trapezium Rule}}

• [http://www.encyclopediaofmath.org/index.php?title=Trapezium_formula&oldid=12696 Trapezium formula. I.P. Mysovskikh], Encyclopedia of Mathematics, ed. M. Hazewinkel
• [http://dedekind.mit.edu/~stevenj/trapezoidal.pdf Notes on the convergence of trapezoidal-rule quadrature]
• [http://www.boost.org/doc/libs/1_66_0/libs/math/doc/html/math_toolkit/trapezoidal.html An implementation of trapezoidal quadrature provided by Boost.Math]

{{Numerical integration}}

{{Calculus topics}}

Category:Numerical integration