math-157-section-7-7-numerical-integrals

1666 days ago by basyrova

We consider the integral \int_a^b f(x) dx

To obtain a numerical approximation to the integral, we set up the following

f(x) is the function we are integrating

a and b are the limits of integration

n is the number of subdivisions (it should be even for Simpson's rule to work)

f(x)=x^2; a=1; b=2; n=4; 
       

The following calculates the \Delta x that is used in virtually all formulas

delta=(b-a)/n; 
       

This variable will be used in sigma notation

j=var("j"); 
       

We start with left-hand Riemann sum

left_sum=delta*sum(f(a + j*delta),j,0,n-1); print left_sum; print N(left_sum); 
       
63/32
1.96875000000000
63/32
1.96875000000000

Now we find the right-hand Riemann sum

right_sum=delta*sum(f(a + j*delta),j,1,n); print right_sum; print N(right_sum); 
       
87/32
2.71875000000000
87/32
2.71875000000000

This uses midpoint rule

midpoint_sum=delta*sum(f(a+(j+1/2)*delta),j,0,n-1); print midpoint_sum; print N(midpoint_sum); 
       
149/64
2.32812500000000
149/64
2.32812500000000

Now we use trapezoid rule for approximation

trapezoid_sum=(delta/2)*(f(a) + f(b) + 2*sum(f(a+j*delta),j,1,n-1)); print trapezoid_sum; print N(trapezoid_sum); 
       
75/32
2.34375000000000
75/32
2.34375000000000

The Simpson's rule

simpson_sum=(delta/3)*(f(a) + f(b) + 2*sum(f(a+j*delta),j,1,n-1) + 2*sum(f(a+delta+2*j*delta),j,0,n/2-1)); print simpson_sum; print N(simpson_sum); 
       
7/3
2.33333333333333
7/3
2.33333333333333