<Text-field style="_cstyle324" layout="Heading 1">Riemann, Trapezoids, and Simpson</Text-field>

Note: You may notice differences between this Maple worksheet and the equivalent Mathematica notebook. These differences were introduced to preserve the content of these modules and were necessary because of major functional differences between Maple and Mathematica.

<Text-field style="Heading 1" layout="Heading 1">Introduction</Text-field>

OBJECTIVE: To visualize the process of using Riemann sums, the trapezoid rule, and Simpson's rule for approximating definite integrals and to understand the error associated with each method. To evaluate a definite integral, it is often necessary to use a numerical estimate of the integral in place of an exact value. This occurs in three instances: 1) when there is no simple formula for an antiderivative of the integrand; 2) when the antiderivative of the integrand is difficult to determine and/or evaluate; and, 3) when the integrand function is represented by a table of numeric values rather than by a formula. In this module, we investigate several numerical methods for estimating definite integrals.

<Text-field style="_cstyle256" layout="Heading 1">Technology Guidelines</Text-field>

NOTE: If you have just finished a worksheet, restart Maple before executing a new worksheet. TO OPEN SECTIONS, Click on the PLUS sign at the left hand side of the screen or select Expand All Sections from the View drop down menu.

TO STOP AN EXECUTION Click on STOP button from the toolbar.

ORDER OF EXECUTION Execute commands in the order given. Do not skip any Maple Input lines within a given worksheet

Alternatively, you can execute the entire worksheet by selecting the Execute Worksheet command from the Edit drop down menu.

SAVING WORKSHEETS.

You can save anytime to any directory you choose, and it is wise to save often. EXPERIENCING MAJOR PROBLEMS Save if appropriate, and then shut down Maple and start it up again.

<Text-field style="Heading 1" layout="Heading 1">Part I: Riemann Sums and Errors</Text-field>

The most direct method for estimating a definite integral is to replace it with a Riemann sum. That is, we estimate the integral as follows: LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KEYrLUYjNiotSShtc3Vic3VwR0YkNictSSNtb0dGJDYuUSYmaW50O0YnLyUnaXRhbGljR1EmZmFsc2VGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRicvJSZmZW5jZUdGPC8lKnNlcGFyYXRvckdGPC8lKXN0cmV0Y2h5R0Y8LyUqc3ltbWV0cmljR0Y8LyUobGFyZ2VvcEdGPC8lLm1vdmFibGVsaW1pdHNHRjwvJSdhY2NlbnRHRjwvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZQLUYsNiRRImFGJy9GPlEnaXRhbGljRictRiw2JFEiYkYnRlYvJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMkYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGKy1GIzYmLUYsNiRRImZGJ0ZWLUY3Ni1RMCZBcHBseUZ1bmN0aW9uO0YnRj1GQEZCRkRGRkZIRkpGTEZORlEtSShtZmVuY2VkR0YkNiUtRiM2JC1GLDYkUSJ4RidGVkZWRjpGPUZWRistSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuM2VtRicvJSZkZXB0aEdGYHAvJSpsaW5lYnJlYWtHUSVhdXRvRictRjc2LlEwJkRpZmZlcmVudGlhbEQ7RidGOkY9RkBGQkZERkZGSEZKRkxGTkZRRmhvRlYtRjc2LlEiPUYnRjpGPUZARkJGREZGRkhGSkZML0ZPUSwwLjI3Nzc3NzhlbUYnL0ZSRmBxLUYjNigtSSttdW5kZXJvdmVyR0YkNictRjc2LlEmJlN1bTtGJ0Y6Rj1GQEZCL0ZFUSV0cnVlRidGRi9GSUZbci9GS0ZbckZMRk4vRlJRLDAuMTY2NjY2N2VtRictRiM2Ji1GLDYkUSJpRidGVkZccS1JI21uR0YkNiVRIjFGJ0Y6Rj1GVi1GLDYkUSJuRidGVkZMLyUsYWNjZW50dW5kZXJHRjxGKy1GXHA2JkZecC9GYnBRJDUuMEYnRmRwRmZwLUYjNidGKy1GIzYmRl1vRmBvLUZkbzYlLUYjNiQtSSVtc3ViR0YkNiUtRiw2JFEiY0YnRlYtRiM2JEZickZWRmhuRlZGOkY9RlYtRjc2LVExJkludmlzaWJsZVRpbWVzO0YnRj1GQEZCRkRGRkZIRkpGTEZORlEtRiw2JFEiaEYnRlZGVkYrRlZGK0ZWRitGVg==, where the interval [a, b] is divided into n subintervals, each of length LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJoRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEiPUYnL0Y4USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQi8lKXN0cmV0Y2h5R0ZCLyUqc3ltbWV0cmljR0ZCLyUobGFyZ2VvcEdGQi8lLm1vdmFibGVsaW1pdHNHRkIvJSdhY2NlbnRHRkIvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZRLUkmbWZyYWNHRiQ2KC1GIzYmRistRiM2Ji1GLDYlUSJiRidGNEY3LUY7Ni1RKCZtaW51cztGJ0Y+RkBGQ0ZFRkdGSUZLRk0vRlBRLDAuMjIyMjIyMmVtRicvRlNGXG8tRiw2JVEiYUYnRjRGN0Y+RitGPi1GIzYkLUYsNiVRIm5GJ0Y0RjdGPi8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGW3AvJSliZXZlbGxlZEdGQkY+RitGPg== , and where LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1YkdGJDYlLUkjbWlHRiQ2JVEiY0YnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GIzYkLUYvNiVRImlGJ0YyRjUvRjZRJ25vcm1hbEYnLyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGPQ== is any value of x taken from the LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEiaUYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GLzYlUSN0aEYnRjJGNS8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRicvRjZRJ25vcm1hbEYn subinterval. To assess the error in the approximation, we will perform a numerical experiment wherein we use a Riemann sum to estimate an integral for which we know the exact value. (You should keep in mind, however, that numerical integration is primarily used for estimating integrals in situations where we can't determine a decimal or fraction representation for the exact value of the integral.) Our first goal is to determine what factors affect the error when we use left- or right-hand Riemann sums to estimate the exact value of an integral. The integral we choose for our experiment is LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2J0YrLUYjNiotSShtc3Vic3VwR0YkNictSSNtb0dGJDYtUSYmaW50O0YnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGPy8lKXN0cmV0Y2h5R0Y/LyUqc3ltbWV0cmljR0Y/LyUobGFyZ2VvcEdGPy8lLm1vdmFibGVsaW1pdHNHRj8vJSdhY2NlbnRHRj8vJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkjbW5HRiQ2JFEiMEYnRjotSSZtZnJhY0dGJDYoLUYjNiQtRiw2JVElJnBpO0YnLyUnaXRhbGljR0Y/RjpGOi1GIzYkLUZSNiRRIjJGJ0Y6RjovJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRmNvLyUpYmV2ZWxsZWRHRj8vJTFzdXBlcnNjcmlwdHNoaWZ0R0Zdby8lL3N1YnNjcmlwdHNoaWZ0R0ZURistRiM2Ji1GLDYlUSRjb3NGJ0ZnbkY6LUY3Ni1RMCZBcHBseUZ1bmN0aW9uO0YnRjpGPUZARkJGREZGRkhGSkZMRk8tSShtZmVuY2VkR0YkNiQtRiM2JC1GLDYlUSJ4RicvRmhuUSV0cnVlRicvRjtRJ2l0YWxpY0YnRjpGOkY6RistSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuM2VtRicvJSZkZXB0aEdGZXEvJSpsaW5lYnJlYWtHUSVhdXRvRictRjc2LVEwJkRpZmZlcmVudGlhbEQ7RidGOkY9RkBGQkZERkZGSEZKRkxGT0ZpcEY6LUY3Ni1RIj1GJ0Y6Rj1GQEZCRkRGRkZIRkovRk1RLDAuMjc3Nzc3OGVtRicvRlBGZXItRlI2JEZgb0Y6RjpGK0Y6. First, we look at a graphical depiction of the situation. For this, the riemannsum( ) command is created using the leftbox( ) and rightbox( ) commands from the student package. The command is riemannsum(f, xrange, n, box, incdec). The arguments are the function, f, the independent variable, x, the range of the integral xrange, the number of rectangles, n, the "right" or "left" sum indicator, box, and an "increasing/decreasing" indicator, incdec. If the integrand function increases over the interval of the integral, then the argument for incdec is the word "increasing," and if it decreases then the argument is the word "decreasing." The error is taken as the Riemann estimate minus the exact value of the integral. The special command riemannsum( ) is only available in this module and is not a built-in Maple command. Here's how it works.

restart; with(student): with(plots):

riemannsum:=proc(f, xrange, n, box, incdec) local a, p1, p2; a:=array(1..2): if box=right then a[1]:=rightbox(f, xrange, n, color=black, title=`Cos(x)`,tickmarks=[4,4]): p1:=rightbox(f, xrange, n, color=black, title=`error`,tickmarks=[4,4]): p2:=plot(f, xrange, filled=true, tickmarks=[4,4]): if incdec=increasing then a[2]:=display(p2,p1) else a[2]:=display(p1,p2) fi: display(a); else if box=left then a[1]:=leftbox(f, xrange, n, color=black, title=`Cos(x)`,tickmarks=[4,4]): p1:=leftbox(f, xrange, n, color=black, title=`error`,tickmarks=[4,4]): p2:=plot(f, xrange, filled=true,tickmarks=[4,4]): if incdec=increasing then a[2]:=display(p1,p2) else a[2]:=display(p2,p1) fi: display(a); fi: fi: end:

riemannsum(cos(x), x=0..Pi/2, 5, right,decreasing); riemannsum(cos(x), x=0..Pi/2, 5, left,decreasing);

The area of each small triangular-shaped region on the error graph is the local error associated with each rectangle in the Riemann sum. The sum of the areas of these triangular-shaped regions is the global or total error in the estimate of the integral. In general, the areas of the rectangles can be positive, negative, or 0. Because these areas are signed, we will refer to them hereafter as "signed areas."

<Text-field style="Heading 1" layout="Heading 1">You Try It: Visualizing the Errors</Text-field>

Use the riemansumm( ) command in the preceding section (copied below) to respond to the items that follow.

n:=5: f:=x->cos(x): a:=0: b:=Pi/2: incdec:=decreasing: riemannsum(f(x), x=a..b, n, right,incdec); riemannsum(f(x), x=a..b, n, left,incdec);

a) Increase the number of rectangles in the preceding input cell, and describe qualitatively what happens to the error in the estimate of the integral each time you double the number of rectangles. b) What effect does the slope of the function have on the error? c) The function LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYnLUkjbWlHRiQ2I1EhRictRiM2Jy1GLDYmUSRjb3NGJy8lJWJvbGRHUSZmYWxzZUYnLyUnaXRhbGljR0Y2LyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjkvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZQLUkobWZlbmNlZEdGJDYlLUYjNiUtRiw2JlEieEYnRjQvRjhRJXRydWVGJy9GOlEnaXRhbGljRicvRjpRJWJvbGRGJy8lK2ZvbnR3ZWlnaHRHRmpuRjRGOUZpbkZbb0YrRmluRltv decreases on the interval from 0 to LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYkLUkjbWlHRiQ2JVElJnBpO0YnLyUnaXRhbGljR1EmZmFsc2VGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRidGNy1GIzYkLUkjbW5HRiQ2JFEiMkYnRjdGNy8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGRS8lKWJldmVsbGVkR0Y2Rjc=, and its first derivative is negative. In this case, the right-hand Riemann sum underestimates the integral (i.e., the error is negative), and the left-hand sum overestimates the integral (i.e., the error is positive). Describe what happens to the error when we use left- and right-hand Riemann sums to estimate LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1JJm1mcmFjR0YkNigtRiM2JC1GLDYlUSUmcGk7RicvJSdpdGFsaWNHRj1GOEY4LUYjNiQtRlA2JFEiMkYnRjhGOC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGYW8vJSliZXZlbGxlZEdGPS8lMXN1cGVyc2NyaXB0c2hpZnRHRltvLyUvc3Vic2NyaXB0c2hpZnRHRlJGKy1GIzYmLUYsNiVRJHNpbkYnRmVuRjgtRjU2LVEwJkFwcGx5RnVuY3Rpb247RidGOEY7Rj5GQEZCRkRGRkZIRkpGTS1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRInhGJy9GZm5RJXRydWVGJy9GOVEnaXRhbGljRidGOEY4RjhGKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0ZjcS8lKmxpbmVicmVha0dRJWF1dG9GJy1GNTYtUTAmRGlmZmVyZW50aWFsRDtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNRmdwRjhGK0Y4. d) While the error for each rectangle (i.e., the local error) on the left-hand sum appears to be nearly equal in magnitude to the local error for each corresponding rectangle on the right-hand sum, they aren't exactly the same because the graph of the function curves. Where does the difference in the local errors appear to be largest, and how is the difference in the local errors related to the second derivative of the function? (To see this effect more dramatically, try using two and/or three rectangles in the riemannsum( ) command.) e) Based on your observations in part (d), specify two different ways in which you could use Riemann sums to improve the estimate of the integral by reducing the errors. The trick is to look for ways to estimate the integrals so that the local errors nearly cancel. For your improved methods, which feature of the integrand function would you expect to have the most significant effect on the error? f) Try the riemannsum( ) command on a function that we can't integrate exactly. Two examples are LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KS1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIxRidGOC1GUDYkUSIyRidGOC8lMXN1cGVyc2NyaXB0c2hpZnRHRlUvJS9zdWJzY3JpcHRzaGlmdEdRIjBGJy1JJm1mcmFjR0YkNigtRiM2JEZPRjgtRiM2JC1GLDYlUSJ4RicvJSdpdGFsaWNHUSV0cnVlRicvRjlRJ2l0YWxpY0YnRjgvJS5saW5ldGhpY2tuZXNzR0ZSLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRmhvLyUpYmV2ZWxsZWRHRj1GKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0ZicC8lKmxpbmVicmVha0dRJWF1dG9GJy1GNTYtUTAmRGlmZmVyZW50aWFsRDtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNRlxvRjhGK0Y4 and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KS1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1GUDYkUSIzRidGOC8lMXN1cGVyc2NyaXB0c2hpZnRHUSIyRicvJS9zdWJzY3JpcHRzaGlmdEdGUi1JJW1zdXBHRiQ2JS1GNTYtUS8mRXhwb25lbnRpYWxFO0YnRjhGO0Y+RkBGQkZERkZGSEZKL0ZOUSwwLjExMTExMTFlbUYnLUYjNiYtRjU2LVEqJnVtaW51czA7RidGOEY7Rj5GQEZCRkRGRkZIL0ZLUSwwLjIyMjIyMjJlbUYnL0ZORmNvLUYjNiQtRmZuNiUtRiw2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnL0Y5USdpdGFsaWNGJy1GUDYkRlhGOC9GV0ZSRjhGK0Y4RmNwRistSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuM2VtRicvJSZkZXB0aEdGaXAvJSpsaW5lYnJlYWtHUSVhdXRvRictRjU2LVEwJkRpZmZlcmVudGlhbEQ7RidGOEY7Rj5GQEZCRkRGRkZIRkpGTUZpb0Y4RitGOA==.

<Text-field style="Heading 1" layout="Heading 1">Part II: Analyzing the Errors</Text-field>

Now that we have made some qualitative observations about Riemann estimates of an integral and the associated errors, we are ready to be more quantitative in our analysis. For this we use the following commands from the student package, leftsum(f, range, n) and rightsum(f, range, n), to calculate left-hand and right-hand Riemann sums. The arguments are f, the integrand function, range, the range of the function to be integrated and n, the number of subintervals. We use these commands, together with the fact that the exact value of the integral LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1JJm1mcmFjR0YkNigtRiM2JC1GLDYlUSUmcGk7RicvJSdpdGFsaWNHRj1GOEY4LUYjNiQtRlA2JFEiMkYnRjhGOC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGYW8vJSliZXZlbGxlZEdGPS8lMXN1cGVyc2NyaXB0c2hpZnRHRltvLyUvc3Vic2NyaXB0c2hpZnRHRlJGKy1GIzYmLUYsNiVRJGNvc0YnRmVuRjgtRjU2LVEwJkFwcGx5RnVuY3Rpb247RidGOEY7Rj5GQEZCRkRGRkZIRkpGTS1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRInhGJy9GZm5RJXRydWVGJy9GOVEnaXRhbGljRidGOEY4RjhGKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0ZjcS8lKmxpbmVicmVha0dRJWF1dG9GJy1GNTYtUTAmRGlmZmVyZW50aWFsRDtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNRmdwRjhGK0Y4 is 1, to build a table of values that includes the step size h, the numerical estimates of the integral, RleftI and RrightI, and the errors in the estimates, Eleft = RelftI-1 and Eright = RrightI-1. We start with two rectangles and double the number of rectangles ten times.

f:=cos(x): a:=0: b:=Pi/2.: exactvalue:=int(f, x=a..b): print(`The exact value of the integral`, exactvalue);

matrix([[n, h, left_est, right_est, left_error, right_error], seq( [2^i, Pi/2/2^i, evalf(leftsum(f, x=a..b, 2^i)), evalf(rightsum(f, x=a..b,2^i)), evalf(leftsum(f, x=a..b, 2^i))-exactvalue, evalf(rightsum(f, x=a..b, 2^i))-exactvalue], i=1..11 )]);

The data in the table should confirm your qualitative observations from the "You Try It: Visualizing the Errors" above, but now we can be more specific. We make the following quantitative observations. 1) As the number of rectangles increases, the error decreases, and, in fact, each time we double the number of rectangles or cut h in half, we cut the error roughly in half. A numerical estimate that exhibits this characteristic is said to have an error of order h, and we designate this as O(h). 2) The left-hand sum overestimates the integral, whereas the right-hand some underestimates it, and, for the same number of rectangles, the errors are roughly equal in magnitude. 3) When the number of rectangles is small, the difference between the errors for the left- and right-hand estimates is larger, showing the effect of the second derivative and the concavity of the graph on the difference in errors. The one effect that is not evident from the numeric data is that the local error is larger when the magnitude of the first derivative, that is, the magnitude of the slope of the graph, is larger. We explore this effect further in Part III below.

<Text-field style="Heading 1" layout="Heading 1">You Try It: On Another Function</Text-field>

Perform a numerical experiment like the one in Part II on the integral LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1JJm1mcmFjR0YkNigtRiM2JC1GLDYlUSUmcGk7RicvJSdpdGFsaWNHRj1GOEY4LUYjNiQtRlA2JFEiMkYnRjhGOC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGYW8vJSliZXZlbGxlZEdGPS8lMXN1cGVyc2NyaXB0c2hpZnRHRltvLyUvc3Vic2NyaXB0c2hpZnRHRlJGKy1GIzYkLUklbXN1cEdGJDYlLUYjNiYtRiw2JVEkc2luRidGZW5GOC1GNTYtUTAmQXBwbHlGdW5jdGlvbjtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNLUkobWZlbmNlZEdGJDYkLUYjNiQtRiw2JVEieEYnL0ZmblEldHJ1ZUYnL0Y5USdpdGFsaWNGJ0Y4RjhGOEZpbi9GZ29GUkY4RistSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuM2VtRicvJSZkZXB0aEdGaXEvJSpsaW5lYnJlYWtHUSVhdXRvRictRjU2LVEwJkRpZmZlcmVudGlhbEQ7RidGOEY7Rj5GQEZCRkRGRkZIRkpGTUZccUY4RitGOA==, and then answer the questions that follow. To help, we include the commands to generate the table of values like the one in Part II.

f:=sin(x)^2: a:=0: b:=Pi/2.: exactvalue:=int(f, x=a..b): print(`The exact value of the integral`,exactvalue);

matrix([[n, h, left_est, right_est, left_err, right_err], seq( [2^i, Pi/2/2^i, evalf(leftsum(f, x=a..b, 2^i)), evalf(rightsum(f, x=a..b,2^i)), evalf(leftsum(f, x=a..b, 2^i))-exactvalue, evalf(rightsum(f, x=a..b, 2^i))-exactvalue], i=1..11 )]);

1. What happens to the error each time the number of rectangles is doubled? What is the order of the errors for each Riemann sum estimate of the integral? 2. Does the left-hand Riemann sum overestimate or underestimate the integral? What about the right-hand Riemann sum? 3. What effect does the second derivative of the integrand have on the difference of the errors for the left and right Riemann sums?

<Text-field style="Heading 1" layout="Heading 1">Part III: The Effect of f ' on the Error</Text-field>

To quantify the effect of the slope on the error, we investigate straight-line functions, LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KEYrLUYjNiYtRiw2JFEiZkYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y3USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQS8lKXN0cmV0Y2h5R0ZBLyUqc3ltbWV0cmljR0ZBLyUobGFyZ2VvcEdGQS8lLm1vdmFibGVsaW1pdHNHRkEvJSdhY2NlbnRHRkEvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZQLUkobWZlbmNlZEdGJDYlLUYjNiQtRiw2JFEieEYnRjZGNi8lJ2l0YWxpY0dGQUY9RjYtRjo2LlEiPUYnRmVuRj1GP0ZCRkRGRkZIRkpGTC9GT1EsMC4yNzc3Nzc4ZW1GJy9GUkZbby1GIzYmLUYsNiRRIm1GJ0Y2LUY6Ni1RMSZJbnZpc2libGVUaW1lcztGJ0Y9Rj9GQkZERkZGSEZKRkxGTkZRRlhGNkYrRjZGK0Y2, with varying slopes. First, let's look at a graphical depiction.

riemannsum:=proc(f, xrange, n, box) local a, p1, p2: a:=array(1..2): if box=right then a[1]:=rightbox(f, xrange, n, color=black, title=`m*x`): p1:=rightbox(f, xrange, n, color=black, title=`error above the curve`): p2:=plot(f, xrange, filled=true): a[2]:=display(p2,p1): display(a); else a[1]:=leftbox(f, xrange, n, color=black, title=`m*x`): p1:=leftbox(f, xrange, n, color=black, title=`error below the curve`): p2:=plot(f, xrange, filled=true): a[2]:=display(p1,p2): display(a); fi: end:

n:=5: m:=2: riemannsum(m*x, x=0..5, n, right); riemannsum(m*x, x=0..5, n, left);

Now we set up a table of numeric values, but in this case we keep the number of rectangles fixed at 100 and change m, the slope of the straight-line graph, starting with LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJtRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEiPUYnL0Y4USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQi8lKXN0cmV0Y2h5R0ZCLyUqc3ltbWV0cmljR0ZCLyUobGFyZ2VvcEdGQi8lLm1vdmFibGVsaW1pdHNHRkIvJSdhY2NlbnRHRkIvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZRLUkjbW5HRiQ2JFElMC4wMUYnRj5GPkYrRj4= and doubling ten times. The exact values of the integrals are LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkmbWZyYWNHRiQ2KC1GIzYmLUkjbW5HRiQ2JFEjMjVGJy8lLG1hdGh2YXJpYW50R1Enbm9ybWFsRictSSNtb0dGJDYtUTEmSW52aXNpYmxlVGltZXM7RidGNC8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGPS8lKXN0cmV0Y2h5R0Y9LyUqc3ltbWV0cmljR0Y9LyUobGFyZ2VvcEdGPS8lLm1vdmFibGVsaW1pdHNHRj0vJSdhY2NlbnRHRj0vJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZMLUkjbWlHRiQ2JVEibUYnLyUnaXRhbGljR1EldHJ1ZUYnL0Y1USdpdGFsaWNGJ0Y0LUYjNiQtRjE2JFEiMkYnRjRGNC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGXG8vJSliZXZlbGxlZEdGPUY0.

unassign('m,n,f,x'): m:=0.1*2^i: f:=m*x: a:=0: b:=5:

exactvalue:=int(f,x=a..b): print(`The exact value of the integral is`, exactvalue); matrix([[slope, exact, left_est, right_est, left_err, right_err], seq( [m(i), exactvalue(i), evalf(leftsum(f, x=a..b, 100)), evalf(rightsum(f, x=a..b,100)), evalf(leftsum(f, x=a..b, 100))-exactvalue, evalf(rightsum(f, x=a..b, 100))-exactvalue], i=0..10 )]);

What happens to the error each time you double the slope of the line? Note that for nonlinear functions, the effect of the slope on the error is local in that it varies from point to point on the graph. But if we double all of the slopes, we would also double the error in using a left- or right-hand Riemann sum to estimate the integral. Note also that the differences in the magnitudes of the errors are all 0. This supports the observation that the difference in the errors is dependent upon the second derivative or the concavity of the graph. One method for improving our numerical estimate of the integral is to add the left- and right-hand estimates and then divide the result by two. The idea that motivates this is that the errors for the two estimates will nearly cancel each other out when we add the left- and right-hand sums together. For linear functions, the second derivative is 0, and the graph has no concavity; therefore, this improvement will give the exact value of the integral.

<Text-field style="Heading 1" layout="Heading 1">You Try It: The Effect of <Equation executable="false" style="2D Comment" input-equation="" display="LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=">LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=</Equation> " on the Error Difference</Text-field>

Design a numerical experiment like the one in Part III to demonstrate the effect of the value of the second derivative on the difference between the errors in the estimate of the integral that are obtained using right- and left-hand Riemann sums. To do this, select a simple function that has a constant second derivative, is easy to integrate, and allows you to change the value of the second derivative. Generate a table of numerical values that shows the values of the second derivative and the difference in the errors for the left and right Riemann sums. Take 1 as the starting value of the second derivative, and double its value ten times. Try functions that are concave down as well as ones that are concave up, and summarize the results of your investigation. (We have put the solution for this problem in the next input section with the code hidden. Try to do it on your own, and then look at our solution.)

<Text-field style="Heading 1" layout="Heading 1">Solution</Text-field>

unassign('m,n,f,x'): m:=2^i: f:=1/2*m*x^2: a:=0: b:=5:exactvalue:=int(f,x=a..b): print(`The exact value of the integral is`, exactvalue); matrix([[`f''(x)`, exact, left_est, right_est, left_err, right_err], seq( [m(i), evalf(exactvalue(i)), evalf(leftsum(f, x=a..b, 100)), evalf(rightsum(f, x=a..b,100)), evalf(leftsum(f, x=a..b, 100))-exactvalue, evalf(rightsum(f, x=a..b, 100))-exactvalue], i=0..10 )]);

<Text-field style="Heading 1" layout="Heading 1">Part IV: Trapezoids and Midpoint Riemann Sums</Text-field>

In this part, we investigate two ways to improve the efficiency of numerical estimates of definite integrals.

<Text-field style="Heading 1" layout="Heading 1">Trapezoids</Text-field>

We achieve the first improvement by adding together the left and right Riemann sum estimates of the integral and dividing by two. The motivation for this approach is that the errors in the left and right sums should nearly cancel each other out. In the next section, we form the new function trapezoid( ) that does this.

trapezoid:=proc(f, a, b, n): (evalf(leftsum(f, x=a..b, n)+rightsum(f, x=a..b, n)))/2: end:

The trapezoid method gets its name from the signed area of a trapezoid in each subinterval that results from adding the signed areas of the left and right rectangles in each subinterval and dividing by 2. That is, LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1JJm1mcmFjR0YkNigtRiM2JkYrLUYjNihGKy1GIzYnRistRiM2Ji1GLDYlUSJmRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEwJkFwcGx5RnVuY3Rpb247RicvRkNRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0ZNLyUpc3RyZXRjaHlHRk0vJSpzeW1tZXRyaWNHRk0vJShsYXJnZW9wR0ZNLyUubW92YWJsZWxpbWl0c0dGTS8lJ2FjY2VudEdGTS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRmZuLUkobWZlbmNlZEdGJDYkLUYjNiQtSSVtc3ViR0YkNiUtRiw2JVEiY0YnRj9GQi1GIzYkLUYsNiVRJWxlZnRGJ0Y/RkJGSS8lL3N1YnNjcmlwdHNoaWZ0R1EiMEYnRklGSUZJLUZGNi1RMSZJbnZpc2libGVUaW1lcztGJ0ZJRktGTkZQRlJGVEZWRlhGWkZnbi1GLDYlUSJoRidGP0ZCRkktRkY2LVEiK0YnRklGS0ZORlBGUkZURlZGWC9GZW5RLDAuMjIyMjIyMmVtRicvRmhuRmZwLUYjNidGKy1GIzYmRjxGRS1Gam42JC1GIzYkLUZfbzYlRmFvLUYjNiQtRiw2JVEmcmlnaHRGJ0Y/RkJGSUZpb0ZJRklGSUZccEZfcEZJRitGSUYrRkktRiM2JC1JI21uR0YkNiRRIjJGJ0ZJRkkvJS5saW5ldGhpY2tuZXNzR1EiMUYnLyUrZGVub21hbGlnbkdRJ2NlbnRlckYnLyUpbnVtYWxpZ25HRmJyLyUpYmV2ZWxsZWRHRk0tRkY2LVEiPUYnRklGS0ZORlBGUkZURlZGWC9GZW5RLDAuMjc3Nzc3OGVtRicvRmhuRltzLUYyNigtRiM2Ji1Gam42JC1GIzYoRitGOkZicEZqcEYrRklGSUZccEZfcEZJRmdxRl1yRmByRmNyRmVyRklGK0ZJ , the area of a trapezoid, where LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYkUSJmRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEwJkFwcGx5RnVuY3Rpb247RicvRjVRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y/LyUpc3RyZXRjaHlHRj8vJSpzeW1tZXRyaWNHRj8vJShsYXJnZW9wR0Y/LyUubW92YWJsZWxpbWl0c0dGPy8lJ2FjY2VudEdGPy8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRk4tSShtZmVuY2VkR0YkNiUtRiM2JC1JJW1zdWJHRiQ2JS1GLDYkUSJjRidGNC1GIzYkLUYsNiRRJWxlZnRGJ0Y0RjQvJS9zdWJzY3JpcHRzaGlmdEdRIjBGJ0Y0LyUnaXRhbGljR0Y/RjtGNEYrRjQ= is the length of one parallel side, LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYkUSJmRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LVEwJkFwcGx5RnVuY3Rpb247RicvRjVRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y/LyUpc3RyZXRjaHlHRj8vJSpzeW1tZXRyaWNHRj8vJShsYXJnZW9wR0Y/LyUubW92YWJsZWxpbWl0c0dGPy8lJ2FjY2VudEdGPy8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRk4tSShtZmVuY2VkR0YkNiUtRiM2JC1JJW1zdWJHRiQ2JS1GLDYkUSJjRidGNC1GIzYkLUYsNiRRJnJpZ2h0RidGNEY0LyUvc3Vic2NyaXB0c2hpZnRHUSIwRidGNC8lJ2l0YWxpY0dGP0Y7RjRGK0Y0 is the length of the other parallel side, and h is the distance between the two sides. Now we use trapezoid( ) to estimate LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1JJm1mcmFjR0YkNigtRiM2JC1GLDYlUSUmcGk7RicvJSdpdGFsaWNHRj1GOEY4LUYjNiQtRlA2JFEiMkYnRjhGOC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGYW8vJSliZXZlbGxlZEdGPS8lMXN1cGVyc2NyaXB0c2hpZnRHRltvLyUvc3Vic2NyaXB0c2hpZnRHRlJGKy1GIzYmLUYsNiVRJGNvc0YnRmVuRjgtRjU2LVEwJkFwcGx5RnVuY3Rpb247RidGOEY7Rj5GQEZCRkRGRkZIRkpGTS1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRInhGJy9GZm5RJXRydWVGJy9GOVEnaXRhbGljRidGOEY4RjhGKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0ZjcS8lKmxpbmVicmVha0dRJWF1dG9GJy1GNTYtUTAmRGlmZmVyZW50aWFsRDtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNRmdwRjhGK0Y4 with ten subintervals and compare the result with the left- and right-hand Riemann sums. First, we calculate the exact value of the integral.

unassign('f, a, b, n, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: exactvalue:=int(f,x=a..b): print(`The exact value of the integral is `,exactvalue);

Now we calculate the estimates of the integral .

trapezoid(f, a, b, 10); evalf(leftsum(f, x=a..b, 10)); evalf(rightsum(f, x=a..b, 10));

And we compute the errors.

trapezoid(f, a, b, 10)-exactvalue; evalf(leftsum(f, x=a..b, 10))-exactvalue; evalf(rightsum(f, x=a..b, 10))-exactvalue;

For the same number of subintervals, we get a much more precise estimate of the integral from the trapezoid method.

<Text-field style="Heading 1" layout="Heading 1">Midpoint Riemann Sums</Text-field>

The second improvement we consider is to use the midpoint of each subinterval in a Riemann sum. The effect of doing this is to nearly cancel the error on either side of each midpoint in each subinterval. If the function we integrate is increasing over a subinterval, then the midpoint rectangle will overestimate the actual area on the left side of the midpoint and will underestimate, by nearly the same amount, the actual area on the right side of the midpoint. Again we expect that these errors will nearly cancel each other out. Here is the riemSumMid( ) command.

riemSumMid:=proc(f,a,b,n): evalf(sum('subs(x=(a+(i-1/2)*(b-a)/n),f)*(b-a)/n', i=1..n)): end:

Let's try it on our old friend and compare the result with left and right Riemann sums and with the trapezoid estimate. First, we calculate the exact value of the integral.

unassign('f, a, b, n, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: exactvalue:=evalf(int(f,x=a..b)): print(`The exact value of the integral is `,exactvalue);

Now we calculate the estimates of the integral.

riemSumMid(cos(x), 0, Pi/2, 10); trapezoid(f, a, b, 10); evalf(leftsum(f, x=a..b, 10)); evalf(rightsum(f, x=a..b, 10));

Again, we compute the errors.

riemSumMid(cos(x), 0, Pi/2, 10)-exactvalue; trapezoid(f, a, b, 10)-exactvalue; evalf(leftsum(f, x=a..b, 10))-exactvalue; evalf(rightsum(f, x=a..b, 10))-exactvalue;

For the same number of intervals, we get a much more precise estimate of the integral using the midpoint rule; it is even more precise than the estimate given by the trapezoid rule. What relationship does the error for the midpoint estimate seem to have with that for the trapezoid estimate?

<Text-field style="Heading 1" layout="Heading 1">You Try It: On Some Other Functions</Text-field>

Compare the Riemann left, Riemann right, trapezoid, and midpoint estimates of some integrals that you pick. How does the error for the midpoint rule compare with the error for the trapezoid rule with the same number of subdivisions? Also, try increasing the number of subintervals. What happens to the errors when you double the number of subintervals?

Here are some integrals that you might want to try:

LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1JJm1mcmFjR0YkNigtRiM2JC1GLDYlUSUmcGk7RicvJSdpdGFsaWNHRj1GOEY4LUYjNiQtRlA2JFEiMkYnRjhGOC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGYW8vJSliZXZlbGxlZEdGPS8lMXN1cGVyc2NyaXB0c2hpZnRHRltvLyUvc3Vic2NyaXB0c2hpZnRHRlJGKy1GIzYmLUYsNiVRJHNpbkYnRmVuRjgtRjU2LVEwJkFwcGx5RnVuY3Rpb247RidGOEY7Rj5GQEZCRkRGRkZIRkpGTS1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRInhGJy9GZm5RJXRydWVGJy9GOVEnaXRhbGljRidGOEY4RjhGKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0ZjcS8lKmxpbmVicmVha0dRJWF1dG9GJy1GNTYtUTAmRGlmZmVyZW50aWFsRDtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNRmdwRjhGK0Y4 , LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KS1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtRiM2JS1GNTYtUSomdW1pbnVzMDtGJ0Y4RjtGPkZARkJGREZGRkgvRktRLDAuMjIyMjIyMmVtRicvRk5GVS1JI21uR0YkNiRRIjFGJ0Y4RjhGVy8lMXN1cGVyc2NyaXB0c2hpZnRHUSIyRicvJS9zdWJzY3JpcHRzaGlmdEdRIjBGJy1JJW1zdXBHRiQ2JS1GNTYtUS8mRXhwb25lbnRpYWxFO0YnRjhGO0Y+RkBGQkZERkZGSEZKL0ZOUSwwLjExMTExMTFlbUYnLUYsNiVRInhGJy8lJ2l0YWxpY0dRJXRydWVGJy9GOVEnaXRhbGljRicvRmZuRmpuRistSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuM2VtRicvJSZkZXB0aEdGYXAvJSpsaW5lYnJlYWtHUSVhdXRvRictRjU2LVEwJkRpZmZlcmVudGlhbEQ7RidGOEY7Rj5GQEZCRkRGRkZIRkpGTUZjb0Y4RitGOA== , LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2KS1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1GUDYkUSIxRidGOC8lMXN1cGVyc2NyaXB0c2hpZnRHUSIyRicvJS9zdWJzY3JpcHRzaGlmdEdGUi1JJm1zcXJ0R0YkNiMtRiM2JkZTLUY1Ni1RKCZtaW51cztGJ0Y4RjtGPkZARkJGREZGRkgvRktRLDAuMjIyMjIyMmVtRicvRk5GXm8tRiw2JVEieEYnLyUnaXRhbGljR1EldHJ1ZUYnL0Y5USdpdGFsaWNGJ0Y4RistSSdtc3BhY2VHRiQ2Ji8lJ2hlaWdodEdRJjAuMGV4RicvJSZ3aWR0aEdRJjAuM2VtRicvJSZkZXB0aEdGXXAvJSpsaW5lYnJlYWtHUSVhdXRvRictRjU2LVEwJkRpZmZlcmVudGlhbEQ7RidGOEY7Rj5GQEZCRkRGRkZIRkpGTUZgb0Y4RitGOA==

To help you out, we provide some commands in the two cells that follow.

unassign('f, a, b, n, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: exactvalue:=evalf(int(f,x=a..b)): n:=10: print(`The exact value of the integral is `, exactvalue); print(`The midpoint estimate is`, evalf(riemSumMid(f, a, b, n))); print(`The trapezoid estimate is`, trapezoid(f, a, b, n)); print(`The left Riemann estimate is`, evalf(leftsum(f, x=a..b,n))); print(`The right Riemann estimate is`, evalf(rightsum(f, x=a..b,n)));

print(`The midpoint error is`, riemSumMid(f, a, b, n)-exactvalue); print(`The trapezoid error is`, trapezoid(f, a, b, n)-exactvalue); print(`The left Riemann error is`, evalf(leftsum(f, x=a..b,n))-exactvalue); print(`The right Riemann estimate is`, evalf(rightsum(f, x=a..b,n))-exactvalue);

<Text-field style="Heading 1" layout="Heading 1">Part V: Experiment with Trapezoids and Midpoints</Text-field>

In this part, we do some numerical experiments to determine what factors affect the errors when we use trapezoids and midpoint Riemann sums to estimate a definite integral.

<Text-field style="Heading 1" layout="Heading 1">The Effect of the Step Size on the Error</Text-field>

First, we would like to determine the effect of the step size h on the error. We build a table similar to the one in Part II.

unassign('f, a, b, n,i,j, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: exactvalue:=int(f,x=a..b): print(`The exact value of the integral is `, exactvalue); matrix([[n, h, trap_est, mid_est, trap_err, mid_err], seq( [2^j, Pi/2/2^j, trapezoid(f, a, b, 2^j), riemSumMid(f, a, b, 2^j), trapezoid(f, a, b, 2^j)-exactvalue, riemSumMid(f, a, b, 2^j)-exactvalue], j=1..11 )]);

<Text-field style="Heading 1" layout="Heading 1">The Effect of <Equation executable="false" style="2D Comment" input-equation="" display="LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=">LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUkjbWlHRiQ2JVEiZkYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy9GM1Enbm9ybWFsRic=</Equation> " on the Error</Text-field>

Based on the observations made in Parts II, III, and IV, we conjecture that the error for these methods depend on the magnitude of the second derivative of the integrand function. We test our conjecture by comparing the error in estimating the area under the quadratic functions of the form LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2J0YrLUYjNiYtRiw2JFEiZkYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnL0Y3USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQS8lKXN0cmV0Y2h5R0ZBLyUqc3ltbWV0cmljR0ZBLyUobGFyZ2VvcEdGQS8lLm1vdmFibGVsaW1pdHNHRkEvJSdhY2NlbnRHRkEvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZQLUkobWZlbmNlZEdGJDYlLUYjNiQtRiw2JFEieEYnRjZGNi8lJ2l0YWxpY0dGQUY9RjYtRjo2LlEiPUYnRmVuRj1GP0ZCRkRGRkZIRkpGTC9GT1EsMC4yNzc3Nzc4ZW1GJy9GUkZbby1JJm1mcmFjR0YkNigtRiM2Jy1GLDYkUSJhRidGNi1GOjYtUTEmSW52aXNpYmxlVGltZXM7RidGPUY/RkJGREZGRkhGSkZMRk5GUS1GIzYkLUklbXN1cEdGJDYlRlgtSSNtbkdGJDYlUSIyRidGZW5GPS8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGNkYrRjYtRiM2JEZdcEY2LyUubGluZXRoaWNrbmVzc0dRIjFGJy8lK2Rlbm9tYWxpZ25HUSdjZW50ZXJGJy8lKW51bWFsaWduR0ZbcS8lKWJldmVsbGVkR0ZBRjZGK0Y2 between x = 0 and x = 5 for various values of a. The exact values of the integrals are LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2J0YrLUYjNiktSShtc3Vic3VwR0YkNictSSNtb0dGJDYtUSYmaW50O0YnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGPy8lKXN0cmV0Y2h5R0Y/LyUqc3ltbWV0cmljR0Y/LyUobGFyZ2VvcEdGPy8lLm1vdmFibGVsaW1pdHNHRj8vJSdhY2NlbnRHRj8vJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkjbW5HRiQ2JFEiMEYnRjotRlI2JFEiNUYnRjovJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMkYnLyUvc3Vic2NyaXB0c2hpZnRHRlQtSSZtZnJhY0dGJDYoLUYjNictRiw2JVEiYUYnLyUnaXRhbGljR1EldHJ1ZUYnL0Y7USdpdGFsaWNGJy1GNzYtUTEmSW52aXNpYmxlVGltZXM7RidGOkY9RkBGQkZERkZGSEZKRkxGTy1GIzYkLUklbXN1cEdGJDYlLUYsNiVRInhGJ0Zfb0Ziby1GUjYkRlpGOi9GWUZURjpGK0Y6LUYjNiRGX3BGOi8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGaXAvJSliZXZlbGxlZEdGP0YrLUknbXNwYWNlR0YkNiYvJSdoZWlnaHRHUSYwLjBleEYnLyUmd2lkdGhHUSYwLjNlbUYnLyUmZGVwdGhHRmNxLyUqbGluZWJyZWFrR1ElYXV0b0YnLUY3Ni1RMCZEaWZmZXJlbnRpYWxEO0YnRjpGPUZARkJGREZGRkhGSkZMRk9GXHBGOi1GNzYtUSI9RidGOkY9RkBGQkZERkZGSEZKL0ZNUSwwLjI3Nzc3NzhlbUYnL0ZQRmNyLUZobjYoLUYjNiYtRlI2JFEkMTI1RidGOkZkb0Zcb0Y6LUYjNiQtRlI2JFEiNkYnRjpGOkZkcEZncEZqcEZccUY6RitGOg== , and the second derivative of f(x) is equal to a. We set up a table in which a starts at 0.1 and doubles ten times. The number of subintervals is kept constant at 100.

unassign('ac,n,i,j,f,x'): ac:=0.1*2^j: f:=1/2*ac*x^2: a:=0: b:=5: exactvalue:=int(f,x=a..b): Digits:=(7); print(`The exact value of the integral is`, exactvalue); matrix([[f_2(x), exact, trap_est, mid_est, trap_err, mid_err], seq( [ac(j), evalf(exactvalue(j)), trapezoid(f,a,b, 100), riemSumMid(f, a,b,100), trapezoid(f, a,b, 100)-exactvalue, riemSumMid(f,a,b,100)-exactvalue], j=0..10 )]);

In "You Try It" below, we ask you to summarize the results of these experiments.

<Text-field style="Heading 1" layout="Heading 1">You Try It: Summarize the Experiment's Results</Text-field>

Summarize the results of the numerical experiments in Part V. In your summary, address the following items. a) How is the error for each method related to the number of subintervals n and the interval width h? What is the order of the error for each method? (Note: If cutting the interval width in half reduces the error by a factor of LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1cEdGJDYlLUkobWZlbmNlZEdGJDYkLUkmbWZyYWNHRiQ2KC1GIzYkLUkjbW5HRiQ2JFEiMUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJ0Y6LUYjNiQtRjc2JFEiMkYnRjpGOi8lLmxpbmV0aGlja25lc3NHRjkvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGRi8lKWJldmVsbGVkR1EmZmFsc2VGJ0Y6LUkjbWlHRiQ2JVEibkYnLyUnaXRhbGljR1EldHJ1ZUYnL0Y7USdpdGFsaWNGJy8lMXN1cGVyc2NyaXB0c2hpZnRHUSIwRidGOg==, then the error is of order LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYkLUklbXN1cEdGJDYlLUkjbWlHRiQ2JVEiaEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1GLzYlUSJuRidGMkY1LyUxc3VwZXJzY3JpcHRzaGlmdEdRIjBGJy9GNlEnbm9ybWFsRic=, and we designate this by LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ji1GLDYlUSJPRicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy1JI21vR0YkNi1RMCZBcHBseUZ1bmN0aW9uO0YnRjcvJSZmZW5jZUdGNi8lKnNlcGFyYXRvckdGNi8lKXN0cmV0Y2h5R0Y2LyUqc3ltbWV0cmljR0Y2LyUobGFyZ2VvcEdGNi8lLm1vdmFibGVsaW1pdHNHRjYvJSdhY2NlbnRHRjYvJSdsc3BhY2VHUSYwLjBlbUYnLyUncnNwYWNlR0ZOLUkobWZlbmNlZEdGJDYlLUYjNiQtSSVtc3VwR0YkNiUtRiw2JFEiaEYnL0Y4USdpdGFsaWNGJy1GLDYkUSJuRidGZm4vJTFzdXBlcnNjcmlwdHNoaWZ0R1EiMEYnRmZuRjRGN0ZmbkYrRmZu.) b) How is the error for the trapezoid method related to the error for the midpoint Riemann sum? c) How is the error related to the second derivative of the integrand function? (To answer this completely, repeat the second experiment in Part V for negative values of a.) d) Are your conclusions consistent with the trapezoid error formula that is found in the text? e) From the pattern of errors for the trapezoid method and the midpoint Riemann sum method, specify how you might combine these two estimates of the integral so that the errors will nearly cancel each other. Compare your idea for combining the estimates with those of other students.

<Text-field style="Heading 1" layout="Heading 1">Part VI: Simpson's Method</Text-field>

Noting that the magnitude of the error for the trapezoid method is approximately two times the magnitude of the error for a midpoint Riemann sum with the same number of subintervals, we can make yet another improvement in the efficiency of our numerical estimations of the integrals. If we add two midpoint estimates and one trapezoid estimate and divide the result by three, the errors for the two methods should nearly cancel. The result of this combination is called Simpson's method. To illustrate this, consider the following estimates of the integral LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi1RJiZpbnQ7RicvJSxtYXRodmFyaWFudEdRJ25vcm1hbEYnLyUmZmVuY2VHUSZmYWxzZUYnLyUqc2VwYXJhdG9yR0Y9LyUpc3RyZXRjaHlHRj0vJSpzeW1tZXRyaWNHRj0vJShsYXJnZW9wR0Y9LyUubW92YWJsZWxpbWl0c0dGPS8lJ2FjY2VudEdGPS8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRkwtSSNtbkdGJDYkUSIwRidGOC1JJm1mcmFjR0YkNigtRiM2JC1GLDYlUSUmcGk7RicvJSdpdGFsaWNHRj1GOEY4LUYjNiQtRlA2JFEiMkYnRjhGOC8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGYW8vJSliZXZlbGxlZEdGPS8lMXN1cGVyc2NyaXB0c2hpZnRHRltvLyUvc3Vic2NyaXB0c2hpZnRHRlJGKy1GIzYmLUYsNiVRJGNvc0YnRmVuRjgtRjU2LVEwJkFwcGx5RnVuY3Rpb247RidGOEY7Rj5GQEZCRkRGRkZIRkpGTS1JKG1mZW5jZWRHRiQ2JC1GIzYkLUYsNiVRInhGJy9GZm5RJXRydWVGJy9GOVEnaXRhbGljRidGOEY4RjhGKy1JJ21zcGFjZUdGJDYmLyUnaGVpZ2h0R1EmMC4wZXhGJy8lJndpZHRoR1EmMC4zZW1GJy8lJmRlcHRoR0ZjcS8lKmxpbmVicmVha0dRJWF1dG9GJy1GNTYtUTAmRGlmZmVyZW50aWFsRDtGJ0Y4RjtGPkZARkJGREZGRkhGSkZNRmdwRjhGK0Y4. We compare the midpoint and trapezoid estimates and their errors, and then we combine them as indicated above and calculate the error for the new estimate.

unassign('f, a, b, n, i,j, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: n:=10: exactvalue:=evalf(int(f,x=a..b)): trap:=trapezoid(f,a,b,n): mid:=riemSumMid(f,a,b,n): simps:=(trap+2*mid)/3: print(`The exact value of the integral is `, exactvalue); print(`The trapezoid estimate is`, trap); print(`The trapezoid error is`, trap-exactvalue); print(`The midpoint estimate is`, mid); print(`The midpoint error is`, mid-exactvalue); print(`The Simpson estimate is`, simps); print(`The error for Simpson's estimate is`, simps-exactvalue);

Since Simpson's method requires three function evaluations in each subinterval, one at each end point for the trapezoid plus one at the mid point, we count the number of subintervals differently. If we take the dividing points for the subintervals to be all of the places where the integrand function is evaluated, then there are actually twice as many subintervals for Simpson's method than there would be for the trapezoid method or the midpoint Riemann sum. One consequence of this is that we need an even number of subintervals for Simpson's method to work. We use Maple's built-in simpson( ) function to estimate the definite integral.

Let's use it to estimate LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYmLUkjbWlHRiQ2I1EhRictRiM2Ki1JKG1zdWJzdXBHRiQ2Jy1JI21vR0YkNi5RJiZpbnQ7RicvJSdpdGFsaWNHUSZmYWxzZUYnLyUsbWF0aHZhcmlhbnRHUSdub3JtYWxGJy8lJmZlbmNlR0Y6LyUqc2VwYXJhdG9yR0Y6LyUpc3RyZXRjaHlHRjovJSpzeW1tZXRyaWNHRjovJShsYXJnZW9wR0Y6LyUubW92YWJsZWxpbWl0c0dGOi8lJ2FjY2VudEdGOi8lJ2xzcGFjZUdRJjAuMGVtRicvJSdyc3BhY2VHRk4tSSNtbkdGJDYlUSIwRidGOEY7LUkmbWZyYWNHRiQ2KC1GIzYkLUYsNiVRJSZwaTtGJ0Y4RjsvRjxRJ2l0YWxpY0YnLUYjNiQtRlI2JVEiMkYnRjhGO0Znbi8lLmxpbmV0aGlja25lc3NHUSIxRicvJStkZW5vbWFsaWduR1EnY2VudGVyRicvJSludW1hbGlnbkdGY28vJSliZXZlbGxlZEdGOi8lMXN1cGVyc2NyaXB0c2hpZnRHRl1vLyUvc3Vic2NyaXB0c2hpZnRHRlRGKy1GIzYmLUYsNiVRJGNvc0YnRjhGOy1GNTYtUTAmQXBwbHlGdW5jdGlvbjtGJ0Y7Rj5GQEZCRkRGRkZIRkpGTEZPLUkobWZlbmNlZEdGJDYlLUYjNiQtRiw2JFEieEYnRmduRmduRjhGO0ZnbkYrLUknbXNwYWNlR0YkNiYvJSdoZWlnaHRHUSYwLjBleEYnLyUmd2lkdGhHUSYwLjNlbUYnLyUmZGVwdGhHRmFxLyUqbGluZWJyZWFrR1ElYXV0b0YnLUY1Ni5RMCZEaWZmZXJlbnRpYWxEO0YnRjhGO0Y+RkBGQkZERkZGSEZKRkxGT0ZpcEZnbkYrRmdu and then calculate the errors, comparing the results with those obtained using trapezoids and midpoint Riemann sums.

unassign('f, a, b, n, i,j, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: n:=6: exactvalue:=evalf(int(f,x=a..b)): print(`The exact value of the integral is `, exactvalue); print(`The Simpson estimate of the integral is`, evalf(simpson(f,x=a..b,n))); print(`The trapezoid estimate of the integral is`, trap); print(`The midpoint estimate is`, riemSumMid(f,a,b,n));

Now we compare the errors.

print(`The error in the Simpson estimate is`, evalf(simpson(f,x=a..b,n))-exactvalue); print(`The error in the trapezoid estimate is`, trapezoid(f,a,b,n)-exactvalue); print(`The error in the midpoint is`, riemSumMid(f,a,b,n)-exactvalue);

<Text-field style="Heading 1" layout="Heading 1"> You Try It: Experiment with Simpson</Text-field>

Perform some numerical experiments like those in the preceding parts to show that the order of the error for Simpson's method is fourth. Also show that the local error is proportional to the fourth derivative of the integrand function. (We have put the solutions for this problem in the next two sections. Try to do it on your own, and then look at our solutions.)

<Text-field style="Heading 1" layout="Heading 1">Solution 1</Text-field>

unassign('f, a, b, n, i,j, exactvalue'): f:=cos(x): a:=0: b:=Pi/2: exactvalue:=evalf(int(f,x=a..b)): Digits:=(20): print(`The exact value of the integral is `, exactvalue); matrix([[n, h, Simpson_estimate, Simpson_error], seq( [2^j, Pi/2/2^j, evalf(simpson(f,x=a..b, 2^j)), evalf(simpson(f,x=a..b,2^j))-exactvalue], j=1..11 )]);

<Text-field style="Heading 1" layout="Heading 1">Solution 2</Text-field>

unassign('f, a, b, n, i,j, exactvalue'): m:=0.1*2^j: f:=m/24*x^4: a:=0: b:=5: exactvalue:=int(f,x=a..b): Digits:=(12): print(`The exact value of the integral is `, exactvalue); matrix([[f_4, h, exact, Simpson_estimate, Simpson_error], seq( [m, Pi/2/2^j, exactvalue, evalf(simpson(f,x=a..b, 100)), evalf(simpson(f,x=a..b,100)-exactvalue)], j=0..10 )]);