Find $\int_{0}^{1} x^2 e^x$ using trapezoidal rule with error less than $10^{-3}$

integrationnumerical methodsnumerical-calculus

The question is to find $$\int_{0}^{1} x^2 e^x$$ using trapezoidal rule with error less than $10^{-3}$

Which is another way to say :”Find h” in the first place.
So we know that :

$f(x)=x^2 e^x$

$f’(x)=2xe^x+x^2e^x$

$f’’(x)=2e^x+4xe^x+x^2e^x= e^x(2+4x+x^2)$

So :
$| f’’(x)|=| 2e^x+4xe^x+x^2e^x | \le 7e$
So the upper bound (M)would be 7e

Using $\frac{b-a}{12} h^2 M\le 10^3$
I have to find h and then use trapezoidal rule to find T(h)
But the problem is h would be really small and using trapezoidal rule would be too long to do without using mathematical programs

Any help ?

Best Answer

You were on the right way: just use the error formula to get

$$h^2 \leq \sqrt{\frac{12 \text{Tol}}{7e}}$$

which yields for $\text{Tol} = 10^{-3}$:

$h \approx 0.025$ and hence $n \approx 40$. All in all, just set $n>40$ and you'll observe an absolute error less that $\text{Tol}$ in your simulation.

The following C++ snippet shows you this fact

#include <iostream>
#include <cmath>
#include <iomanip>      // std::setprecision


double trapz(double f(double),int n, int a, int b){
    
    double sum{f(a) + f(b)};
    double h = double(b-a)/double(n);
    for (int i=1; i<n; ++i) {
        sum+=2*f(a+i*h);
//        std::cout << a+i*h <<std::endl; bound checking
    }
    return (0.5*h)*sum;
}

double f(double x){
    
    return x*x*exp(x);
}


int main(){
    double exact{exp(1)- 2};
    int n;
    std::cout <<"Set n: " <<std::endl;
    std::cin >> n;
    
    double res{trapz(f,n,0,1)};
    std::cout << "Numerical integration with n= " << n << " yields: "<< std::setprecision(9) << res
    << " with an error of " << std::scientific <<  std::fabs(exact - res) << std::endl;
    

    return 0;
    
}

If you're not familiar with C/C++, just go to your terminal (in MacOS and Linux it's the same) and type the following in order to compile your source code

g++ -o file_name -std=c++14 file_name.cpp

Then run it with

./file_name

to see the output

With $n=40$ I get:

Numerical integration with n= 40 yields: 0.718706544 with an error of 4.247156195e-04