Skip to main content
Mathematics LibreTexts

3.6: Numerical Integration

  • Page ID
    130088
    • Gilbert Strang & Edwin “Jed” Herman
    • OpenStax

    \( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)

    \( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)

    \( \newcommand{\id}{\mathrm{id}}\) \( \newcommand{\Span}{\mathrm{span}}\)

    ( \newcommand{\kernel}{\mathrm{null}\,}\) \( \newcommand{\range}{\mathrm{range}\,}\)

    \( \newcommand{\RealPart}{\mathrm{Re}}\) \( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)

    \( \newcommand{\Argument}{\mathrm{Arg}}\) \( \newcommand{\norm}[1]{\| #1 \|}\)

    \( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)

    \( \newcommand{\Span}{\mathrm{span}}\)

    \( \newcommand{\id}{\mathrm{id}}\)

    \( \newcommand{\Span}{\mathrm{span}}\)

    \( \newcommand{\kernel}{\mathrm{null}\,}\)

    \( \newcommand{\range}{\mathrm{range}\,}\)

    \( \newcommand{\RealPart}{\mathrm{Re}}\)

    \( \newcommand{\ImaginaryPart}{\mathrm{Im}}\)

    \( \newcommand{\Argument}{\mathrm{Arg}}\)

    \( \newcommand{\norm}[1]{\| #1 \|}\)

    \( \newcommand{\inner}[2]{\langle #1, #2 \rangle}\)

    \( \newcommand{\Span}{\mathrm{span}}\) \( \newcommand{\AA}{\unicode[.8,0]{x212B}}\)

    \( \newcommand{\vectorA}[1]{\vec{#1}}      % arrow\)

    \( \newcommand{\vectorAt}[1]{\vec{\text{#1}}}      % arrow\)

    \( \newcommand{\vectorB}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)

    \( \newcommand{\vectorC}[1]{\textbf{#1}} \)

    \( \newcommand{\vectorD}[1]{\overrightarrow{#1}} \)

    \( \newcommand{\vectorDt}[1]{\overrightarrow{\text{#1}}} \)

    \( \newcommand{\vectE}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash{\mathbf {#1}}}} \)

    \( \newcommand{\vecs}[1]{\overset { \scriptstyle \rightharpoonup} {\mathbf{#1}} } \)

    \( \newcommand{\vecd}[1]{\overset{-\!-\!\rightharpoonup}{\vphantom{a}\smash {#1}}} \)

    \(\newcommand{\avec}{\mathbf a}\) \(\newcommand{\bvec}{\mathbf b}\) \(\newcommand{\cvec}{\mathbf c}\) \(\newcommand{\dvec}{\mathbf d}\) \(\newcommand{\dtil}{\widetilde{\mathbf d}}\) \(\newcommand{\evec}{\mathbf e}\) \(\newcommand{\fvec}{\mathbf f}\) \(\newcommand{\nvec}{\mathbf n}\) \(\newcommand{\pvec}{\mathbf p}\) \(\newcommand{\qvec}{\mathbf q}\) \(\newcommand{\svec}{\mathbf s}\) \(\newcommand{\tvec}{\mathbf t}\) \(\newcommand{\uvec}{\mathbf u}\) \(\newcommand{\vvec}{\mathbf v}\) \(\newcommand{\wvec}{\mathbf w}\) \(\newcommand{\xvec}{\mathbf x}\) \(\newcommand{\yvec}{\mathbf y}\) \(\newcommand{\zvec}{\mathbf z}\) \(\newcommand{\rvec}{\mathbf r}\) \(\newcommand{\mvec}{\mathbf m}\) \(\newcommand{\zerovec}{\mathbf 0}\) \(\newcommand{\onevec}{\mathbf 1}\) \(\newcommand{\real}{\mathbb R}\) \(\newcommand{\twovec}[2]{\left[\begin{array}{r}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\ctwovec}[2]{\left[\begin{array}{c}#1 \\ #2 \end{array}\right]}\) \(\newcommand{\threevec}[3]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\cthreevec}[3]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \end{array}\right]}\) \(\newcommand{\fourvec}[4]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\cfourvec}[4]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \end{array}\right]}\) \(\newcommand{\fivevec}[5]{\left[\begin{array}{r}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\cfivevec}[5]{\left[\begin{array}{c}#1 \\ #2 \\ #3 \\ #4 \\ #5 \\ \end{array}\right]}\) \(\newcommand{\mattwo}[4]{\left[\begin{array}{rr}#1 \amp #2 \\ #3 \amp #4 \\ \end{array}\right]}\) \(\newcommand{\laspan}[1]{\text{Span}\{#1\}}\) \(\newcommand{\bcal}{\cal B}\) \(\newcommand{\ccal}{\cal C}\) \(\newcommand{\scal}{\cal S}\) \(\newcommand{\wcal}{\cal W}\) \(\newcommand{\ecal}{\cal E}\) \(\newcommand{\coords}[2]{\left\{#1\right\}_{#2}}\) \(\newcommand{\gray}[1]{\color{gray}{#1}}\) \(\newcommand{\lgray}[1]{\color{lightgray}{#1}}\) \(\newcommand{\rank}{\operatorname{rank}}\) \(\newcommand{\row}{\text{Row}}\) \(\newcommand{\col}{\text{Col}}\) \(\renewcommand{\row}{\text{Row}}\) \(\newcommand{\nul}{\text{Nul}}\) \(\newcommand{\var}{\text{Var}}\) \(\newcommand{\corr}{\text{corr}}\) \(\newcommand{\len}[1]{\left|#1\right|}\) \(\newcommand{\bbar}{\overline{\bvec}}\) \(\newcommand{\bhat}{\widehat{\bvec}}\) \(\newcommand{\bperp}{\bvec^\perp}\) \(\newcommand{\xhat}{\widehat{\xvec}}\) \(\newcommand{\vhat}{\widehat{\vvec}}\) \(\newcommand{\uhat}{\widehat{\uvec}}\) \(\newcommand{\what}{\widehat{\wvec}}\) \(\newcommand{\Sighat}{\widehat{\Sigma}}\) \(\newcommand{\lt}{<}\) \(\newcommand{\gt}{>}\) \(\newcommand{\amp}{&}\) \(\definecolor{fillinmathshade}{gray}{0.9}\)
    Learning Objectives
    • Approximate the value of a definite integral by using the midpoint and trapezoidal rules.
    • Determine the absolute and relative error in using a numerical integration technique.
    • Estimate the absolute and relative error using an error-bound formula.
    • Recognize when the midpoint and trapezoidal rules over- or underestimate the true value of an integral.
    • Use Simpson’s rule to approximate the value of a definite integral to a given accuracy.

    The antiderivatives of many functions either cannot be expressed or cannot be expressed easily in closed form (that is, in terms of known functions). Consequently, rather than evaluate definite integrals of these functions directly, we resort to various techniques of numerical integration to approximate their values. In this section, we explore several of these techniques. In addition, we examine the process of estimating the error in using these techniques.

    The Midpoint Rule

    Earlier in this text we defined the definite integral of a function over an interval as the limit of Riemann sums. In general, any Riemann sum of a function \( f(x)\) over an interval \([a,b]\) may be viewed as an estimate of \(\displaystyle ∫^b_af(x)\,dx\). Recall that a Riemann sum of a function \( f(x)\) over an interval \( [a,b]\) is obtained by selecting a partition

    \[ P=\{x_0,x_1,x_2,…,x_n\} \nonumber \]

    where \(\quad a=x_0<x_1<x_2<⋯<x_n=b \)

    and a set

    \[ S=\{x^*_1,x^*_2,…,x^*_n\} \nonumber \]

    where \(x_{i−1}≤x^*_i≤x_i \quad \text{for all} \, i.\)

    The Riemann sum corresponding to the partition \(P\) and the set \(S\) is given by \(\displaystyle \sum^n_{i=1}f(x^*_i)Δx_i\), where \( Δx_i=x_i−x_{i−1},\) the length of the \( i^{\text{th}}\) subinterval.

    The midpoint rule for estimating a definite integral uses a Riemann sum with subintervals of equal width and the midpoints, \( m_i\), of each subinterval in place of \( x^*_i\). Formally, we state a theorem regarding the convergence of the midpoint rule as follows.

    The Midpoint Rule

    Assume that \( f(x)\) is continuous on \([a,b]\). Let \( n\) be a positive integer and \( Δx=\dfrac{b−a}{n}\). If \( [a,b]\) is divided into \( n\) subintervals, each of length \( Δx\), and \( m_i\) is the midpoint of the \( i^{\text{th}}\) subinterval, set

    \[M_n=\sum_{i=1}^nf(m_i)Δx. \nonumber \]

    Then \(\displaystyle \lim_{n→∞}M_n=∫^b_af(x)\,dx.\)

    As we can see in Figure \(\PageIndex{1}\), if \( f(x)≥0\) over \( [a,b]\), then \(\displaystyle \sum^n_{i=1}f(m_i)Δx\) corresponds to the sum of the areas of rectangles approximating the area between the graph of \( f(x)\) and the \(x\)-axis over \([a,b]\). The graph shows the rectangles corresponding to \(M_4\) for a nonnegative function over a closed interval \([a,b].\)

    This figure is a graph of a non-negative function in the first quadrant. The function increases and decreases. The quadrant is divided into a grid. Beginning on the x-axis at the point labeled a = x sub 0, there are rectangles shaded whose heights are approximately the height of the curve. The x-axis is scaled by increments of msub1, x sub 1, m sub 2, x sub 2, m sub 3, x sub 3, m sub 4 and b = x sub 4.
    Figure \(\PageIndex{1}\): The midpoint rule approximates the area between the graph of \(f(x)\) and the \(x\)-axis by summing the areas of rectangles with midpoints that are points on \(f(x)\).
    Example \(\PageIndex{1}\): Using the Midpoint Rule with \(M_4\)

    Use the midpoint rule to estimate \(\displaystyle ∫^1_0x^2\,dx\) using four subintervals. Compare the result with the actual value of this integral.

    Solution: Each subinterval has length \( Δx=\dfrac{1−0}{4}=\dfrac{1}{4}.\) Therefore, the subintervals consist of

    \[\left[0,\tfrac{1}{4}\right],\,\left[\tfrac{1}{4},\tfrac{1}{2}\right],\,\left[\tfrac{1}{2},\tfrac{3}{4}\right],\, \text{and}\, \left[\tfrac{3}{4},1\right].\nonumber \]

    The midpoints of these subintervals are \(\left\{\frac{1}{8},\,\frac{3}{8},\,\frac{5}{8},\, \frac{7}{8}\right\}.\) Thus,

    \[\begin{align*} M_4 &=\frac{1}{4}\cdot f\left(\frac{1}{8}\right)+\frac{1}{4}\cdot f\left(\frac{3}{8}\right)+\frac{1}{4}\cdot f\left(\frac{5}{8}\right)+\frac{1}{4}\cdot f\left(\frac{7}{8}\right) \\[4pt] &=\frac{1}{4}⋅\frac{1}{64}+\frac{1}{4}⋅\frac{9}{64}+\frac{1}{4}⋅\frac{25}{64}+\frac{1}{4}⋅\frac{49}{64}\\[4pt] &=\frac{21}{64} = 0.328125. \end{align*}\]

    Since

    \[ ∫^1_0x^2\,dx=\frac{1}{3},\nonumber \]

    the absolute error in this approximation is:

    \[\left\lvert\dfrac{1}{3}−\dfrac{21}{64}\right\rvert=\dfrac{1}{192}≈0.0052, \nonumber \]

    and we see that the midpoint rule produces an estimate that is somewhat close to the actual value of the definite integral.

    Example \(\PageIndex{2}\): Using the Midpoint Rule with \(M_6\)

    Use \(M_6\) to estimate the length of the curve \(y=\frac{1}{2}x^2\) on \([1,4]\).

    Solution: The length of \(y=\frac{1}{2}x^2\) on \([1,4]\) is

    \[s = ∫^4_1\sqrt{1+\left(\frac{dy}{dx}\right)^2}\,dx.\nonumber \]

    Since \(\dfrac{dy}{dx}=x\), this integral becomes \(\displaystyle ∫^4_1\sqrt{1+x^2}\,dx.\)

    If \([1,4]\) is divided into six subintervals, then each subinterval has length \(Δx=\dfrac{4−1}{6}=\dfrac{1}{2}\) and the midpoints of the subintervals are \(\left\{\frac{5}{4},\frac{7}{4},\frac{9}{4},\frac{11}{4},\frac{13}{4},\frac{15}{4}\right\}\). If we set \(f(x)=\sqrt{1+x^2}\),

    \[\begin{align*} M_6 &=\tfrac{1}{2}\cdot f\left(\frac{5}{4}\right)+\tfrac{1}{2}\cdot f\left(\frac{7}{4}\right)+\frac{1}{2}\cdot f\left(\frac{9}{4}\right)+\frac{1}{2}\cdot f\left(\frac{11}{4}\right)+\frac{1}{2}\cdot f\left(\frac{13}{4}\right)+\frac{1}{2}\cdot f\left(\frac{15}{4}\right) \\[4pt] &≈\frac{1}{2}(1.6008+2.0156+2.4622+2.9262+3.4004+3.8810)=8.1431 \, \text{ units}. \end{align*}\]

    Exercise \(\PageIndex{1}\)

    Use the midpoint rule with \( n=2\) to estimate \(\displaystyle ∫^2_1\frac{1}{x}\,dx.\)

    Hint

    \( Δx=\frac{1}{2}, \quad m_1=\frac{5}{4},\quad \text{and} \quad m_2=\frac{7}{4}.\)

    Answer

    \(\dfrac{24}{35}\approx 0.685714\)

    The Trapezoidal Rule

    We can also approximate the value of a definite integral by using trapezoids rather than rectangles. In Figure \(\PageIndex{2}\), the area beneath the curve is approximated by trapezoids rather than by rectangles.

    This figure is a graph of a non-negative function in the first quadrant. The function increases and decreases. The quadrant is divided into a grid. Beginning on the x-axis at the point labeled a = x sub 0, there are trapezoids shaded whose heights are approximately the height of the curve. The x-axis is scaled by increments of a = x sub 0, x sub 1, x sub 2, x sub 3, and b = x sub 4.
    Figure \(\PageIndex{2}\): Trapezoids may be used to approximate the area under a curve, hence approximating the definite integral.

    The trapezoidal rule for estimating definite integrals uses trapezoids rather than rectangles to approximate the area under a curve. To gain insight into the final form of the rule, consider the trapezoids shown in Figure \(\PageIndex{2}\). We assume that the length of each subinterval is given by \(Δx\). First, recall that the area of a trapezoid with a height of \(h\) and bases of length \(b_1\) and \(b_2\) is given by \(\text{Area}=\frac{1}{2}h(b_1+b_2)\). We see that the first trapezoid has a height \(Δx\) and parallel bases of length \( f(x_0)\) and \( f(x_1)\). Thus, the area of the first trapezoid in Figure \(\PageIndex{2}\) is

    \[ \frac{1}{2}Δx\Big(f(x_0)+f(x_1)\Big).\nonumber \]

    The areas of the remaining three trapezoids are

    \(\dfrac{1}{2}Δx\Big(f(x_1)+f(x_2)\Big),\, \dfrac{1}{2}Δx\Big(f(x_2)+f(x_3)\Big),\) and \( \dfrac{1}{2}Δx\Big(f(x_3)+f(x_4)\Big).\)

    Consequently,

    \[∫^b_af(x)\,dx≈\frac{1}{2}Δx\Big(f(x_0)+f(x_1)\Big)+\frac{1}{2}Δx\Big(f(x_1)+f(x_2)\Big)+\frac{1}{2}Δx\Big(f(x_2)+f(x_3)\Big)+\frac{1}{2}Δx\Big(f(x_3)+f(x_4)\Big).\nonumber \]

    After taking out a common factor of \(\frac{1}{2}Δx\) and combining like terms, we have

    \[ ∫^b_af(x)\,dx≈\frac{Δx}{2}\Big[f(x_0)+2\,f(x_1)+2\,f(x_2)+2\,f(x_3)+f(x_4)\Big].\nonumber \]

    Generalizing, we formally state the following rule.

    The Trapezoidal Rule

    Assume that \(f(x)\) is continuous over \([a,b]\). Let \(n\) be a positive integer and \(Δx=\dfrac{b−a}{n}\). Let \( [a,b]\) be divided into \(n\) subintervals, each of length \(Δx\), with endpoints at \( P=\{x_0,x_1,x_2…,x_n\}.\)

    Set

    \[T_n=\frac{Δx}{2}\Big[f(x_0)+2\, f(x_1)+2\, f(x_2)+⋯+2\, f(x_{n−1})+f(x_n)\Big]. \nonumber \]

    Then, \(\displaystyle \lim_{n→+∞}T_n=∫^b_af(x)\,dx.\)

    Before continuing, let’s make a few observations about the trapezoidal rule. First of all, it is useful to note that

    \(T_n=\dfrac{1}{2}(L_n+R_n)\) where \(\displaystyle L_n=\sum_{i=1}^nf(x_{i−1})Δx\) and \(\displaystyle R_n=\sum_{i=1}^nf(x_i)Δx.\)

    That is, \(L_n\) and \(R_n\) approximate the integral using the left-hand and right-hand endpoints of each subinterval, respectively. In addition, a careful examination of Figure \(\PageIndex{3}\) leads us to make the following observations about using the trapezoidal rules and midpoint rules to estimate the definite integral of a nonnegative function. The trapezoidal rule tends to overestimate the value of a definite integral systematically over intervals where the function is concave up and to underestimate the value of a definite integral systematically over intervals where the function is concave down. On the other hand, the midpoint rule tends to average out these errors somewhat by partially overestimating and partially underestimating the value of the definite integral over these same types of intervals. This leads us to hypothesize that, in general, the midpoint rule tends to be more accurate than the trapezoidal rule.

    This figure has two graphs, both of the same non-negative function in the first quadrant. The function increases and decreases. The quadrant is divided into a grid. The first graph, beginning on the x-axis at the point labeled a = x sub 0, there are trapezoids shaded whose heights are approximately the height of the curve. The x-axis is scaled by increments of a = x sub 0, xsub1, x sub 2, x sub 3, and b = x sub 4. The second graph has on the x-axis at the point labeled a = x sub 0. There are rectangles shaded whose heights are approximately the height of the curve. The x-axis is scaled by increments of m sub 1, x sub 1, m sub 2, x sub 2, m sub 3, x sub 3, m sub 4 and b = x sub 4.
    Figure \(\PageIndex{3}\):The trapezoidal rule tends to be less accurate than the midpoint rule.
    Example \(\PageIndex{3}\): Using the Trapezoidal Rule

    Use the trapezoidal rule to estimate \(\displaystyle ∫^1_0x^2\,dx\) using four subintervals.

    Solution

    The endpoints of the subintervals consist of elements of the set \(P=\left\{0,\frac{1}{4},\, \frac{1}{2},\, \frac{3}{4},1\right\}\) and \(Δx=\frac{1−0}{4}=\frac{1}{4}.\) Thus,

    \[\begin{align*} ∫^1_0x^2dx &≈\frac{1}{2}⋅\frac{1}{4}\Big[f(0)+2\, f\left(\tfrac{1}{4}\right)+2\, f\left(\tfrac{1}{2}\right)+2\, f\left(\tfrac{3}{4}\right)+f(1)\Big] \\[4pt]
    &=\tfrac{1}{8}\big(0+2⋅\tfrac{1}{16}+2⋅\tfrac{1}{4}+2⋅\tfrac{9}{16}+1\big) \\[4pt] &=\frac{11}{32} = 0.34375\end{align*}\]

    Exercise \(\PageIndex{2}\)

    Use the trapezoidal rule with \( n=2\) to estimate \(\displaystyle ∫^2_1\frac{1}{x}\,dx.\)

    Hint

    Set \(Δx=\dfrac{1}{2}.\) The endpoints of the subintervals are the elements of the set \(P=\left\{1,\frac{3}{2},2\right\}.\)

    Answer

    \(\dfrac{17}{24} \approx 0.708333\)

    Absolute and Relative Error

    An important aspect of using these numerical approximation rules consists of calculating the error in using them for estimating the value of a definite integral. We first need to define absolute error and relative error.

    Definition: absolute and relative error

    If \(B\) is our estimate of some quantity having an actual value of \(A\), then the absolute error is given by \(|A−B|\).

    The relative error is the error as a percentage of the actual value and is given by \[\left\lvert\frac{A−B}{A}\right\rvert⋅100\%. \nonumber \]

    Example \(\PageIndex{4}\): Calculating Error in the Midpoint Rule

    Calculate the absolute and relative error in the estimate of \(\displaystyle ∫^1_0x^2\,dx\) using the midpoint rule, found in Example \(\PageIndex{1}\).

    Solution: The calculated value is \(\displaystyle ∫^1_0x^2\,dx=\frac{1}{3}\) and our estimate from the example is \(M_4=\frac{21}{64}\). Thus, the absolute error is given by \(\left\lvert\frac{1}{3}−\frac{21}{64}\right\rvert=\frac{1}{192}≈0.0052.\)

    The relative error is \[\frac{1/192}{1/3}=\frac{1}{64}≈0.015625≈1.6\%.\nonumber \]

    Example \(\PageIndex{5}\): Calculating Error in the Trapezoidal Rule

    Calculate the absolute and relative error in the estimate of \(\displaystyle ∫^1_0x^2\,dx\) using the trapezoidal rule, found in Example \(\PageIndex{3}\).

    Solution: The calculated value is \(\displaystyle ∫^1_0x^2\,dx=\frac{1}{3}\) and our estimate from the example is \(T_4=\frac{11}{32}\). Thus, the absolute error is given by \(\left\lvert\frac{1}{3}−\frac{11}{32}\right\rvert=\frac{1}{96}≈0.0104.\)

    The relative error is given by \[\frac{1/96}{1/3}=0.03125≈3.1\%.\nonumber \]

    Exercise \(\PageIndex{3}\)

    In an earlier checkpoint, we estimated \(\displaystyle ∫^2_1\frac{1}{x}\,dx\) to be \(\frac{24}{35}\) using \(M_2\). The actual value of this integral is \(\ln 2\). Using \(\frac{24}{35}≈0.6857\) and \(\ln 2≈0.6931,\) calculate the absolute error and the relative error.

    Hint

    Use the previous examples as a guide.

    Answer

    absolute error \(\approx 0.0074,\) and relative error \(\approx 1.1\%\)

    Error Bounds on the Midpoint and Trapezoidal Rules

    In the two previous examples, we were able to compare our estimate of an integral with the actual value of the integral; however, we do not typically have this luxury. In general, if we are approximating an integral, we are doing so because we cannot compute the exact value of the integral itself easily. Therefore, it is often helpful to be able to determine an upper bound for the error in an approximation of an integral. The following theorem provides error bounds for the midpoint and trapezoidal rules. The theorem is stated without proof.

    Error Bounds for the Midpoint and Trapezoidal Rules

    Let \(f(x)\) be a continuous function over \([a,b]\), having a second derivative \(f''(x)\) over this interval. If \(M\) is the maximum value of \(|f''(x)|\) over \([a,b]\), then the upper bounds for the error in using \(M_n\) and \(T_n\) to estimate \(\displaystyle ∫^b_af(x)\,dx\) are

    \[\text{Error in}\, M_n≤\frac{M(b−a)^3}{24n^2}\label{MidError} \]

    and

    \[\text{Error in}\, T_n≤\frac{M(b−a)^3}{12n^2} \nonumber \].

    We can use these bounds to determine the value of \(n\) necessary to guarantee that the error in an estimate is less than a specified value.

    Example \(\displaystyle \PageIndex{6}\): Determining the Number of Intervals to Use

    What value of \(n\) should be used to guarantee that an estimate of \(\displaystyle ∫^1_0e^{x^2}\,dx\) is accurate to within \(0.01\) if we use the midpoint rule?

    Solution

    We begin by determining the value of \(M\), the maximum value of \( |f''(x)|\) over \( [0,1]\) for \( f(x)=e^{x^2}\). Since \( f′(x)=2xe^{x^2},\) we have

    \[ f''(x)=2e^{x^2}+4x^2e^{x^2}.\nonumber \]

    Thus,

    \[ |f''(x)|=2e^{x^2}(1+2x^2)≤2⋅e⋅3=6e.\nonumber \]

    From the error-bound Equation \(\ref{MidError}\), we have

    \[ \text{Error in}\, M_n≤\frac{M(b−a)^3}{24n^2}≤\frac{6e(1−0)^3}{24n^2}=\frac{6e}{24n^2}.\nonumber \]

    Now we solve the following inequality for \(n\):

    \[\frac{6e}{24n^2}≤0.01.\nonumber \]

    Thus, \(n≥\sqrt{\frac{600e}{24}}≈8.24.\) Since \(n\) must be an integer satisfying this inequality, a choice of \(n=9\) would guarantee that

    \[ \left\lvert ∫^1_0e^{x^2}\,dx−M_n \right\rvert <0.01.\nonumber \]

    Analysis

    We might have been tempted to round \(8.24\) down and choose \(n=8\), but this would be incorrect because we must have an integer greater than or equal to \(8.24\). We need to keep in mind that the error estimates provide an upper bound only for the error. The actual estimate may, in fact, be a much better approximation than is indicated by the error bound.

    Exercise \(\PageIndex{4}\)

    Use Equation \(\ref{MidError}\) to find an upper bound for the error in using \(M_4\) to estimate \(\displaystyle ∫^1_0x^2\,dx.\)

    Hint

    \(f''(x)=2,\) so \(M=2.\)

    Answer

    \(\dfrac{1}{192}\)

    Simpson’s Rule

    With the midpoint rule, we estimated areas of regions under curves by using rectangles. In a sense, we approximated the curve with piecewise constant functions. With the trapezoidal rule, we approximated the curve by using piecewise linear functions. What if we were, instead, to approximate a curve using piecewise quadratic functions? With Simpson’s rule, we do just this. We partition the interval into an even number of subintervals, each of equal width. Over the first pair of subintervals we approximate \(\displaystyle ∫^{x_2}_{x_0}f(x)\,dx\) with \(\displaystyle ∫^{x_2}_{x_0}p(x)\,dx\), where \(p(x)=Ax^2+Bx+C\) is the quadratic function passing through \((x_0,f(x_0)), \,(x_1,f(x_1)),\) and \((x_2,f(x_2))\) (Figure \(\PageIndex{4}\)). Over the next pair of subintervals we approximate \(\displaystyle ∫^{x_4}_{x_2}f(x)\,dx\) with the integral of another quadratic function passing through \( (x_2,f(x_2)), \,(x_3,f(x_3)),\) and \((x_4,f(x_4)).\) This process is continued with each successive pair of subintervals.

    This figure has two graphs, both of the same non-negative function in the first quadrant. The function increases and decreases. The quadrant is divided into a grid. The first graph, beginning on the x-axis at the point labeled x sub 0, there are trapezoids shaded whose heights are represented by the function p(x), which is a curve following an approximate path of the original graph. The x-axis is scaled by increments of x sub 0, x sub 1, x sub 2. The second graph has on the x-axis at the point labeled x sub 0. There are shaded regions under the curve, divided by x sub 0, x sub 1, x sub 2, x sub 3, and x sub 4. The curve is sectioned into two different parts above the shaded areas. These two parts are labeled p sub 1(x) and p sub 2(x).
    Figure \(\PageIndex{4}\): With Simpson’s rule, we approximate a definite integral by integrating a piecewise quadratic function.

    To understand the formula that we obtain for Simpson’s rule, we begin by deriving a formula for this approximation over the first two subintervals. As we go through the derivation, we need to keep in mind the following relationships:

    \[f(x_0)=p(x_0)=Ax_0^2+Bx_0+C \nonumber \]

    \[f(x_1)=p(x_1)=Ax_1^2+Bx_1+C \nonumber \]

    \[f(x_2)=p(x_2)=Ax_2^2+Bx_2+C \nonumber \]

    \(x_2−x_0=2Δx\), where \(Δx\) is the length of a subinterval.

    \(x_2+x_0=2x_1,\) since \(x_1=\dfrac{(x_2+x_0)}{2}\).

    Thus,

    \[\begin{align*} ∫^{x_2}_{x_0}f(x)\,dx &≈∫^{x_2}_{x_0}p(x)\,dx \\[4pt]
    &=∫^{x_2}_{x_0}(Ax^2+Bx+C)\,dx \\[4pt]
    &=\left(\frac{A}{3}x^3+\frac{B}{2}x^2+Cx\right)\bigg|^{x_2}_{x_0} & &\text{Find the antiderivative.} \\[4pt]
    &=\frac{A}{3}(x_2^3−x_0^3)+\frac{B}{2}(x_2^2−x_0^2)+C(x_2−x_0) & &\text{Evaluate the antiderivative.} \\[4pt]
    &=\frac{A}{3}(x_2−x_0)(x_2^2+x_2x_0+x_0^2)+\frac{B}{2}(x_2−x_0)(x_2+x_0)+C(x_2−x_0) \\[4pt]
    &=\frac{x_2−x_0}{6}\bigg(2A(x_2^2+x_2x_0+x_0^2)+3B(x_2+x_0)+6C\bigg) & &\text{Factor out}\, \frac{x_2−x_0}{6}. \\[4pt]
    &=\frac{Δx}{3}\bigg((Ax_2^2+Bx_2+C)+(Ax_0^2+Bx_0+C)+A(x_2^2+2x_2x_0+x_0^2)+2B(x_2+x_0)+4C\bigg) & &\text{Rearrange the terms. Note:} \enspace Δx = \frac{x_2−x_0}{2} \\[4pt]
    &=\frac{Δx}{3}\big(f(x_2)+f(x_0)+A(x_2+x_0)^2+2B(x_2+x_0)+4C\big) & &\text{Factor and substitute:} \\[4pt]
    & & &\quad f(x_2)=Ax_2^2+Bx_2+C \enspace \text{and} \enspace f(x_0)=Ax_0^2+Bx_0+C. \\[4pt]
    &=\frac{Δx}{3}\big(f(x_2)+f(x_0)+A(2x_1)^2+2B(2x_1)+4C\big) & &\text{Substitute}\, x_2+x_0=2x_1. \\[4pt]
    & & &\quad\text{Note:}\, x_1 = \frac{x_2+x_0}{2} \enspace \text{is the midpoint.} \\[4pt]
    &=\frac{Δx}{3}\big(f(x_2)+4f(x_1)+f(x_0)\big). & &\text{Expand and substitute} \, f(x_1)=Ax_1^2+Bx_1+C. \end{align*}\]

    If we approximate \(\displaystyle ∫^{x_4}_{x_2}f(x)\,dx\) using the same method, we see that we have

    \[∫^{x_4}_{x_2}f(x)\,dx≈\frac{Δx}{3}(f(x_4)+4\,f(x_3)+f(x_2)).\nonumber \]

    Combining these two approximations, we get

    \[∫^{x_4}_{x_0}f(x)\,dx≈\frac{Δx}{3}(f(x_0)+4\,f(x_1)+2\,f(x_2)+4\,f(x_3)+f(x_4)).\nonumber \]

    The pattern continues as we add pairs of subintervals to our approximation. The general rule may be stated as follows.

    Simpson’s Rule

    Assume that \(f(x)\) is continuous over \([a,b]\). Let \(n\) be a positive even integer and \(Δx=\dfrac{b−a}{n}\). Let \([a,b]\) be divided into \(n\) subintervals, each of length \(Δx\), with endpoints at \(P=\{x_0,x_1,x_2,…,x_n\}.\) Set

    \[S_n=\frac{Δx}{3}\Big[f(x_0)+4\,f(x_1)+2\,f(x_2)+4\,f(x_3)+2\,f(x_4)+⋯+2\,f(x_{n−2})+4\,f(x_{n−1})+f(x_n)\Big]. \nonumber \]

    Then,

    \[\lim_{n→+∞}S_n=∫^b_af(x)\,dx.\nonumber \]

    Just as the trapezoidal rule is the average of the left-hand and right-hand rules for estimating definite integrals, Simpson’s rule may be obtained from the midpoint and trapezoidal rules by using a weighted average. It can be shown that \(S_{2n}=\left(\frac{2}{3}\right)M_n+\left(\frac{1}{3}\right)T_n\).

    It is also possible to put a bound on the error when using Simpson’s rule to approximate a definite integral. The bound in the error is given by the following rule:

    Rule: Error Bound for Simpson’s Rule

    Let \(f(x)\) be a continuous function over \([a,b]\) having a fourth derivative, \( f^{(4)}(x)\), over this interval. If \(M\) is the maximum value of \(∣f^{(4)}(x)∣\) over \([a,b]\), then the upper bound for the error in using \(S_n\) to estimate \(\displaystyle ∫^b_af(x)\,dx\) is given by

    \[\text{Error in}\, S_n≤\frac{M(b−a)^5}{180n^4}. \nonumber \]

    Example \(\displaystyle \PageIndex{7}\): Applying Simpson’s Rule 1

    Use \(S_2\) to approximate \(\displaystyle ∫^1_0x^3\,dx\). Estimate a bound for the error in \(S_2\).

    Solution

    Since \([0,1]\) is divided into two intervals, each subinterval has length \(Δx=\frac{1−0}{2}=\frac{1}{2}\). The endpoints of these subintervals are \(\left\{0,\frac{1}{2},1\right\}\). If we set \(f(x)=x^3,\) then

    \[S_2=\frac{1}{3}⋅\frac{1}{2}(f(0)+4\,f(\frac{1}{2})+f(1))=\frac{1}{6}(0+4⋅\frac{1}{8}+1)=\frac{1}{4}.\nonumber \]

    Since \( f^{(4)}(x)=0\) and consequently \(M=0,\) we see that

    Error in \(S_2≤\frac{0(1)^5}{180⋅2^4}=0.\)

    This bound indicates that the value obtained through Simpson’s rule is exact. A quick check will verify that, in fact, \(\displaystyle ∫^1_0x^3\,dx=\frac{1}{4}.\)

    Example \(\displaystyle \PageIndex{8}\): Applying Simpson’s Rule 2

    Use \(S_6\) to estimate the length of the curve \(y=\frac{1}{2}x^2\) over \([1,4].\)

    Solution

    The length of \(y=\frac{1}{2}x^2\) over \([1,4]\) is \(\displaystyle ∫^4_1\sqrt{1+x^2}\,dx\). If we divide \([1,4]\) into six subintervals, then each subinterval has length \(Δx=\frac{4−1}{6}=\frac{1}{2}\), and the endpoints of the subintervals are \( \left\{1,\frac{3}{2},2,\frac{5}{2},3,\frac{7}{2},4\right\}.\) Setting \( f(x)=\sqrt{1+x^2}\),

    \[S_6=\frac{1}{3}⋅\frac{1}{2}(f(1)+4f(\frac{3}{2})+2f(2)+4f(\frac{5}{2})+2f(3)+4f(\frac{7}{2})+f(4)).\nonumber \]

    After substituting, we have

    \[S_6=\frac{1}{6}(1.4142+4⋅1.80278+2⋅2.23607+4⋅2.69258+2⋅3.16228+4⋅3.64005+4.12311)≈8.14594\,\text{units}. \nonumber \]

    Exercise \(\PageIndex{5}\)

    Use \(S_2\) to estimate \(\displaystyle ∫^2_1\frac{1}{x}\,dx.\)

    Hint

    \[S_2=\frac{1}{3}Δx\left(f(x_0)+4f(x_1)+f(x_2)\right) \nonumber \]

    Answer

    \(\frac{25}{36} \approx 0.694444\)


    Key Concepts

    • We can use numerical integration to estimate the values of definite integrals when a closed form of the integral is difficult to find or when an approximate value only of the definite integral is needed.
    • The most commonly used techniques for numerical integration are the midpoint rule, trapezoidal rule, and Simpson’s rule.
    • The midpoint rule approximates the definite integral using rectangular regions whereas the trapezoidal rule approximates the definite integral using trapezoidal approximations.
    • Simpson’s rule approximates the definite integral by first approximating the original function using piecewise quadratic functions.

    Key Equations

    • Midpoint rule

    \(\displaystyle M_n=\sum^n_{i=1}f(m_i)Δx\)

    • Trapezoidal rule

    \(T_n=\frac{Δx}{2}\Big[f(x_0)+2\,f(x_1)+2\,f(x_2)+⋯+2\,f(x_{n−1})+f(x_n)\Big]\)

    • Simpson’s rule

    \(S_n=\frac{Δx}{3}\Big[f(x_0)+4\,f(x_1)+2\,f(x_2)+4\,f(x_3)+2\,f(x_4)+4\,f(x_5)+⋯+2\,f(x_{n−2})+4\,f(x_{n−1})+f(x_n)\Big]\)

    • Error bound for midpoint rule

    Error in \(M_n≤\dfrac{M(b−a)^3}{24n^2}\), where \(M\) is the maximum value of \(|f''(x)|\) over \([a,b]\).

    • Error bound for trapezoidal rule

    Error in \(T_n≤\dfrac{M(b−a)^3}{12n^2}\), where \(M\) is the maximum value of \(|f''(x)|\) over \([a,b]\).

    • Error bound for Simpson’s rule

    Error in \(S_n≤\dfrac{M(b−a)^5}{180n^4}\), where \(M\) is the maximum value of \(∣f^{(4)}(x)∣\) over \([a,b]\).


    Glossary

    absolute error
    if \(B\) is an estimate of some quantity having an actual value of \(A\), then the absolute error is given by \( |A−B|\)
    midpoint rule
    a rule that uses a Riemann sum of the form \(\displaystyle M_n=\sum^n_{i=1}f(m_i)Δx\), where \( m_i\) is the midpoint of the \(i^{\text{th}}\) subinterval to approximate \(\displaystyle ∫^b_af(x)\,dx\)
    numerical integration
    the variety of numerical methods used to estimate the value of a definite integral, including the midpoint rule, trapezoidal rule, and Simpson’s rule
    relative error
    error as a percentage of the actual value, given by \[\text{relative error}=\left|\frac{A−B}{A}\right|⋅100\%\nonumber \]
    Simpson’s rule
    a rule that approximates \(\displaystyle ∫^b_af(x)\,dx\) using the area under a piecewise quadratic function.
    The approximation \(S_n\) to \(\displaystyle ∫^b_af(x)\,dx\) is given by \[S_n=\frac{Δx}{3}\big(f(x_0)+4\,f(x_1)+2\,f(x_2)+4\,f(x_3)+2\,f(x_4)+⋯+2\,f(x_{n−2})+4\,f(x_{n−1})+f(x_n)\big).\nonumber \]
    trapezoidal rule
    a rule that approximates \(\displaystyle ∫^b_af(x)\,dx\) using the area of trapezoids.
    The approximation \(T_n\) to \(\displaystyle ∫^b_af(x)\,dx\) is given by \[T_n=\frac{Δx}{2}\big(f(x_0)+2\, f(x_1)+2\, f(x_2)+⋯+2\, f(x_{n−1})+f(x_n)\big).\nonumber \]

    Contributors and Attributions

    • Gilbert Strang (MIT) and Edwin “Jed” Herman (Harvey Mudd) with many contributing authors. This content by OpenStax is licensed with a CC-BY-SA-NC 4.0 license. Download for free at http://cnx.org.

    • Edited by Paul Seeburger (Monroe Community College). Notes added to development of area under a parabola and typos fixed in original text.

    This page titled 3.6: Numerical Integration is shared under a CC BY-NC-SA 4.0 license and was authored, remixed, and/or curated by Gilbert Strang & Edwin “Jed” Herman (OpenStax) via source content that was edited to the style and standards of the LibreTexts platform.