Previous: 15.2 Complex Numbers
Up: 15 Engineering Programming Examples
Next: 15.4 Summary
Previous Page: 15.2.3 Impedance of Electrical Circuits
Next Page: 15.4 Summary
Another common operation that arises in engineering and scientific computing is integration. While software application packages exist to perform symbolic integration, some functions do not lend themselves to such a ``closed form'' method. A common computing method for approximating the value of an integral in numeric integration. In this section we will develop a small program implementing Simpson's Rule for numeric integration.
The integral of a function between specified limits
gives the area under the function curve as shown in Figure 15.23.
Numeric methods can
approximate the area under the curve by summing approximate sub-areas under
linearized parts of the function at uniformly sampled points
(Figure 15.24). The smaller the sampling interval, ,
the greater the precision of the computed integral.
An algorithm to evaluate such an integral may be written in terms of
the value of the function at sample points between the two limits.
For example, assume the limits
of integration for function, are
and
.
Then, the function
values between the two limits at intervals of h are:
The total number of samples is
There are many
methods to approximate the value of an integral in terms of these sample
values. Simpson's Rule gives a fairly accurate integral of function
between specified limits
and
:
where,
for
,
is the sampling interval,
and
with
adjusted so that
is an even integer.
Except for
the multiplier
,
the above sum is called the Simpson sum.
Observe that in
the Simpson sum, sample values of the function evaluated at odd
,
i.e.
,
are multiplied by 4, and sample values at even values of
, except for
and
,
are multiplied by 2. Finally, sample values
and
are added
without a multiplier.
We will now slightly modify the concept of a generic sum from Chapter
to implement a
function that numerically evaluates an integral of a specified function between
two limits, i.e.
modify the generic function, sum(), into a generic Simpson sum
function. Since integral computation requires real numbers, we use type
double for all our computation.
The parameters to simpsum are the function pointer, fp,
a lower limit, a,
an update function pointer, up,
a sampling interval, step, and an upper
limit, b. The code is shown in Figure 15.25.
The integer variable, m represents the sample number.
If m is zero, the function
sample is added to the cumulative sum, if it is odd,
the sample value times 4 is
added to the cumulative sum; otherwise, sample value times 2 is added.
Finally,
the sample value at b is added and
the resulting Simpson sum is returned. It is
easy now to implement the function integral()
to compute the integral of a
function between limits a and b.
It merely gets the Simpson sum, multiplies by
step/3 and returns it as seen in Figure 15.26
The update function incr() merely returns the value of its first argument
increased by the value of the second argument, step.
This function is included
in sumutil.c together with other useful functions,
self(), square(), and cube() shown in
Figure 15.27.
Finally,
we write a simple driver that computes integrals of several functions using
integral() shown in Figure 15.28.
The program first reads the sampling interval, h; then
repeatedly reads the
integration limits until EOF.
For each set of limits, it calculates the number
of samples, n, for the specified h.
Since the Simpson sum requires an even number
of samples, n is increased by one if it is odd,
and the sampling interval h is
adjusted to correspond to the even value of n.
Then, the program computes the
integral by calling integral() for three different functions:
a straight line
, a square,
, and a cube,
.
The values of integrals are printed.
The program is in three source files, which must be compiled and linked:
integr.c, sumutil.c, and simputil.c.
Here are two sample sessions with different sampling intervals:
Remember, the smaller the sampling interval, the greater the accuracy of the computed integral. The first session specifies a fairly large sampling interval of 0.1 and the results are not very accurate. The exact answers for the integrals are 0.5, 0.3333, and 0.25. The second session specifies a somewhat better sampling interval 0.01, and the results are quite accurate. A smaller sampling interval would be even better, but would require more computation time. A compromise between accuracy and speed is required in most numeric computations.