I'm putting the Chapter 4 skills to this numerical integration
exercise, but I always get the results different than those in the
book (BTW, it's an exercise for finance subjects). The results in the
book use 8 decimal places, while what I'm getting here is just 6
places. Each of the results for each of the intervals are very much
closer to each other, that's why I guess it might be the reason for
differences between my results and those in the book. Or maybe I did
something wrong in the codes. Do you know how to increase decimal
places? Thank you so much.
// Integration
// f(x) = exp(-x^2)
// a = 0; b = 2; start from 4-interval initial partition
// tolerance level = 0.5*10^-7
#include "std_lib_facilities.h"
#include <cmath>
double midpoint (double, double, double);
double trapez (double, double, double);
double simps (double, double, double);
double square (double);
int main()
{
double a = 0;
double b = 0;
double n = 4; // try 4-interval initial partition
const double tol = 0.00000005; // tolerance level
double midpoint_i_old = 0;
double midpoint_i_new = 0;
double trapez_i_old = 0;
double trapez_i_new = 0;
double simps_i_old = 0;
double simps_i_new = 0;
double midpoint_approx = 0;
double trapez_approx = 0;
double simps_approx = 0;
cout << "Enter lower bound (a): ";
cin >> a;
cout << "Enter upper bound (b): ";
cin >> b;
midpoint_i_old = midpoint(a, b, n);
trapez_i_old = trapez(a, b, n);
simps_i_old = simps(a, b, n);
n = 2*n; // double interval each time around
midpoint_i_new = midpoint(a, b, n);
trapez_i_new = trapez(a, b, n);
simps_i_new = simps(a, b, n);
cout << endl;
// Midpoint
while (abs(midpoint_i_new - midpoint_i_old)>tol) {
midpoint_i_old = midpoint_i_new;
n *= 2;
midpoint_i_new = midpoint(a, b, n);
cout << "n: " << n << "\tMidpoint: " << midpoint_i_new << endl;
}
midpoint_approx = midpoint_i_new;
cout << "Approx. Midpoint = " << midpoint_approx << endl << endl;
// Trapezoidal
n = 8; // reset n
while (abs(trapez_i_new - trapez_i_old)>tol) {
trapez_i_old = trapez_i_new;
n *= 2;
trapez_i_new = trapez(a, b, n);
cout << "n: " << n << "\tTrapezoidal: " << trapez_i_new << endl;
}
trapez_approx = trapez_i_new;
cout << "Approx. Trapezoidal = " << trapez_approx << endl << endl;
// Simpson
n = 8; // reset n
while (abs(simps_i_new - simps_i_old)>tol) {
simps_i_old = simps_i_new;
n *= 2;
simps_i_new = simps(a, b, n);
cout << "n: " << n << "\tSimpson: " << simps_i_new << endl;
}
simps_approx = simps_i_new;
cout << "Approx. Simpson = " << simps_approx << endl << endl;
}
double square (double x)
{
return x * x;
}
double midpoint (double a, double b, double n)
{
double h = (b - a)/n;
double i_midpoint = 0;
for (int i=1; i<=n; ++i)
i_midpoint = i_midpoint + exp(-square(a+h*(i-0.5)));
return h * i_midpoint;
}
double trapez (double a, double b, double n)
{
double h = (b - a)/n;
double i_trapez = 0.5*exp(-square(a)) + 0.5*exp(-square(b));
for (int i=1; i<=(n-1);++i)
i_trapez = i_trapez + exp(-square(a + i*h));
return h * i_trapez;
}
double simps (double a, double b, double n)
{
double h = (b - a)/n;
double i_simps = (exp(-square(a)))/6.0 + (exp(-square(b)))/6.0;
for (int i=1; i<=(n-1); ++i)
i_simps = i_simps + (exp(-square(a+i*h)))/3.0;
for (int j=1; j<=n; ++j)
i_simps = i_simps + (2.0*exp(-square(a+h*(j-0.5))))/3.0;
return h * i_simps;
}
Thank a lot for the reply. I like your book a lot, an excellent book
for self-study (harder than other beginner books, but a lot more
interesting).