CS 2073 -- Program 5

CS 2073, Fall 2005
Program 5
Numerical Integration
    Week 5: Sep 19-23
Due (on time): 2005-10-03  23:59:59
Due (late):        2005-10-05  23:59:59

Program 5 must be emailed to: nrwagner@cs.utsa.edu
following directions for: running and submitting a C program, with deadlines:

2005-10-03 23:59:59 (that's Monday, 3 October 2005, 11:59:59 pm) for full credit.
2005-10-05 23:59:59 (that's Wednesday, 5 October 2005, 11:59:59 pm) for 75% credit.

Introduction: For this assignment, we want a program that will do numerical integration. You don't really need to know any calculus, since for us the integral of a function will just be the area under its graph, or its average value. This assignment will use three numerical integration methods:

the Trapezoid method,
Simpson's method, and
a Monte-Carlo method.

See Random Numbers for information about generating the random number needed for this program.

See Comparisons: while and for loops for examples of loops that add sequences of numbers.

Details of the three methods: We will be finding the value of the integral of a function f(x), for x from a to b. You will also start with an integer n representing the number of intervals to divide the segment from a to b into. The size of each interval is h = (b - a)/n. Given these starting values, the trapezoid method uses the formula:


Trapezoid-Integral = (approximately)
     (h/2)[f(a) + 2 f(a+h) + 2 f(a+2h) + 2 f(a+3h) + 2 f(a+4h) + . . . + 2 f(a+i h) +
                             . . .  + 2 f(a+(n-2)h) + 2 f(a+(n-1)h) + f(b)]

Similarly, Simpson's method uses the following formula, where we want n to be even, so that the alternating weights of 1, 4, 2, 4 will come out as shown: with 2, 4, 1.


Simpson's-Integral = (approximately)
     (h/3)[f(a) + 4 f(a+h) + 2 f(a+2h) + 4 f(a+3h) + 2 f(a+4h) + . . .
                                     + 2 f(a+(n-2)h) + 4 f(a+(n-1)h) + f(b)]

Finally, a Monte-Carlo method might use:


Monte-Carlo-Integral = (approximately)
      h [f(x₁) + f(x₂) + f(x₃) + . . . + f(x_n-1) +  f(x_n)]

Here the numbers x₁, x₂, . . . , x_n are randomly chosen from the interval from a to b.

Integrals to approximate numerically: You should use the specific function f, and find its approximate integral:

f(x) = 1/(1 + x²), for x from 0 to 1, or:

In the case of f we can determine the exact answer analytically: The indefinite integral is just arctan(x), which is 0 when x = 0, so the exact answer to this part is arctan(1), which is pi/4 = 0.785398163397448. For this function, Simpson's rule is very accurate, while the Trapezoid method is intermediate, and we always expect a Monte Carlo method to be the least accurate.

After determining the integral for f, you should use the same code to find numerical integral approximations for the function g below. This function is used in statistics. There is no formula for the exact answer in this case.

g(x) = e^-x², for x from -1 to 1, or:

More details: So write a program that computes the first function f, using each of the three integration methods, and using n = 10, n = 1000 and n = 100000. (There should be 9 answers altogether.)

After that you should change your program to compute the second function g and produce another 9 answers.

The output should be clearly labeled with the function (f or g above), the values of a and b, the value of n, and the integration method.

Each integration method should be implemented as a call to one of three functions, with parameters a, b, and c.

Hints and suggestions: Don't panic when computing the sum above. You have terms for values of i from 0 to n inclusive. With the Trapezoid rule each term is multiplied by 2 except for the 0th and the nth terms.

With Simpson's rule, the 0th and the nth terms are multiplied by 1. Otherwise, odd-numbered terms are multiplied by 4 and even-numbered terms are multiplied by 2. (Simpson's rule requires that n be an even number.). To check that term number i is odd, you can check:


if (i%2 == 1) { multiplier = 4 }

The Monte-Carlo method is even easier. Just add up n terms of the form f(random point) and then multiply by h = (b - a)/n.

In order to get n to take on the three successive values 10, 1000, and 100000, you can use the for loop:


for (n = 10; n <= 100000; n *= 100) { /* body of loop */ }

Extra credit part: For extra credit, try to find an approximate value for the following. After you have found this approximate value, square it and see if you recognize this number. (Hint: the function g(x) gets extremely small as abs(x) increases.)

Super Hint: In class on Wednesday, 28 Sept. 2005, students were having trouble with this program, so I actually worked out the details of the Trapezoid method for them. In fairness, I should make this code available to everyone, and you can mimic the code for the other parts of the assignment. If you have already done the program, please stick with your original version, assuming you got close to the correct answer. There is no single "correct" way to write these programs. Remember that Simpson's method will be more accurate, and the Monte Carlo method will be much less accurate. Trapezoid Method.

What you should email: Refer to the submissions directions and to deadlines at the top of this page. The text file that you submit should first have Your Name, the Course Number, and the Program Number. The rest of the file should have the following in it, in the order below, and clearly labeled, including at the beginning the appropriate item letters: a, b, c, etc.

Contents of email submission for Program 5:
Last Name, First Name; Course Number; Program Number.

C source for your program, integral.c. (Or whatever you wish to name it.)
Results of a run of the program used to determine the integral of f three ways, each for three values of n.
Results of another run doing the same as b above for the function g.
(Extra credit) Results of a run computing the approximate integral of g from minus infinity to plus infinity. Then the square of this answer.

Revision date: 2005-09-23. (Please use ISO 8601, the International Standard.)