License Public Domain
Lines 44
Keywords
Num_integration (3)
Permissions
Owner: jsimones
Group Owner: SnortSnort
Viewable by Everyone
Editable by All Siafoo Users
Hide
Don't get spied on – We respect your privacy and provide numerous options to protect it. Join Siafoo Now or Learn More

Numerical Integration Atom Feed 0

# 's
 1#include <iostream>  // Not needed for integration routines
2#include <cmath> // sqrt, exp, acos
3using namespace std;
4
5double PI = acos(-1);
6
7double x_squared(double x) {
8 return x * x;
9}
10
11double normal_dist(double x) {
12 double m = 0, s = 1;
13 return exp((x - m) * (m - x) / 2 / s / s) / s / sqrt(2 * PI);
14}
15
16double rectangle_method(double (*func)(double), double x1, double x2, int n) {
17 /*
18 Integrate func between x1 and x2 using the rectangle method with the midpoint
19 approximation for n samples.
20
21 */
22 double step = (x2 - x1) / n, integral = 0;
23 for (x1 += step / 2; x1 < x2; x1 += step) {
24 integral += (*func)(x1) * step;
25 }
26 return integral;
27}
28
29double simpsons_rule(double (*func)(double), double x1, double x2, int n) {
30 /*
31 Integrate func between x1 and x2 using simpson's rule for n samples.
32
33 */
34 double step = (x2 - x1) / n, integral = 0;
35 for (x1; x1 < x2; x1 += step) {
36 integral += step / 6 * ((*func)(x1) + 4 * (*func)(x1 + step / 2) +
37 (*func)(x1 + step));
38 }
39 return integral;
40}
41
42int main() {
43 // Test using x_squared:
44 double integral1 = rectangle_method(x_squared, 1, 10, 1000);
45 double integral2 = simpsons_rule(x_squared, 1, 10, 1000);
46 cout << integral1 << " " << integral2 << endl;
47
48 // Test using normal_dist:
49 integral1 = rectangle_method(normal_dist, -10, 10, 1000);
50 integral2 = simpsons_rule(normal_dist, -10, 10, 1000);
51 cout << integral1 << " " << integral2 << endl;
52
53 return 0;
54}