Answer to Question #296528 in C++ for satu

Question

Develop user-friendly codes to implement the Newton-Cotes integration formulae as discussed in class. Use appropriate function names and arguments. Remember to check in which cases a certain algorithm will be applicable. 

Answer 1

#include <iostream>
#include <cmath>
#include <cstdlib>
#include <cmath>
using namespace std;


enum QuadRule {trapesoid, simpson, boole};


double integrate(double fun(double), double a, double b, int n, 
                 QuadRule type=simpson)
{
    double h, x, sum;


    if (type ==  simpson && n%2 != 0) {
        cerr << "integrate error: n must be even" << endl;
        exit(1);
    }
    if (type ==  boole && n%4 != 0) {
        cerr << "integrate error: n must be even" << endl;
        exit(1);
    }

    h = (b - a) / n;
    sum = 0.0;
    for (int i=0; i<=n; i++) {
        x = a + i*h;
        if (type == trapesoid) {
            if (i == 0 || i == n ) {
                sum += fun(x)/2;
            }
            else {
                sum += fun(x);
            }
        }
        else if (type == simpson) {
            if (i == 0 || i == n ) {
                sum += fun(x);
            }
            else if (i%2 == 1) {
                sum += 4*fun(x);
            }
            else {
                sum += 2*fun(x);
            }
        }
        else if (type == boole) {
            if (i == 0 || i == n ) {
                sum += 7*fun(x);
            }
            else if (i%2 == 1) {
                sum += 32*fun(x);
            }
            else if (i%4 == 2) {
                sum += 12*fun(x);
            }
            else {
                sum += 14*fun(x);
            }
        }
    }

    double res;
    if (type == trapesoid) {
        res = sum * h;
    }
    else if (type == simpson) {
        res = sum * h/3;
    }
    else if (type == boole) {
        res = sum * 2*h/45;
    }
    return res;
}


int main() {
    const int n = 20;
    const double pi = 3.141592653589793;
    double q1, q2, q4;

    q1 = integrate(sin, 0.0, pi, n, trapesoid);
    q2 = integrate(sin, 0.0, pi, n, simpson);
    q4 = integrate(sin, 0.0, pi, n, boole);

    cout << "Calculetae integral of Sin(x) from 0 to pi  with " << n << " points" << endl;
    cout << "Trapesoid (Newton-Cortes 1) " << q1 << endl;
    cout << "Simpson   (Newton-Cortes 2) " << q2 << endl;
    cout << "Boole     (Newton-Cortes 4) " << q4 << endl;
    cout << "Precise value               " << 2.0 << endl;

    return 0;

}