As we saw in Section
4.6, Ada provides a standard library of elementary mathematical functions,
called Ada.Numerics.Elementary_Functions. The full specification
is given in
Appendix
F; here are a few examples of the functions provided:
Function Purpose Arctan( X) Returns the angle y in radians satisfying X = tan( y) where -pi/2 <= y <= pi/2 Cos( X) Returns the cosine of angle X (in radians) Exp( X) Returns eX where e = 2.71828... Log( X) Returns the natural logarithm of X for X > 0.0 Sin( X) Returns the sine of angle X (in radians) Sqrt( X) Returns the square root of X for X >= 0.0
All functions take arguments of type Float (or a subtype
thereof) and return a value of type Float. The arguments for
Sin and Cos must be expressed in radians, not
degrees. The arguments for Log and Sqrt must be
nonnegative; a negative argument will cause the exception
Ada.Numerics.Argument_Error to be raised.
Example 7.10
The predefined exponentiation operator in Ada does not apply to
floating-point exponents. This means that it is not possible to write
xy directly when x and y are type
Float. Some, but not all, math packages provide such an
exponentiation operation. However, from the study of logarithms we know that
ln(xy) = y x ln(x)
and
z = eln(z)
where e is 2.71828.... So we can derive that
xy = e(y x ln(x))
Assuming that Exp and Ln are available in
Math, this formula can be implemented in Ada as
XToPowerY := Ada.Numerics.Elementary_FunctionsExp( y * Ada.Numerics.Elementary_FunctionsLn( X))
In this book we have been faithful to the convention that all calls to package-provided procedures and functions be prefixed with the name of the package, as in
Ada.Text_IO.New_Line;
Ada.Integer_Text_IO.Get (Item => Next_Num);
Y := Ada.Numerics.Elementary_Functions.Sqrt ( X);
Prefixing
the name of the package is called qualification. There are two main
advantages to doing this. First, the reader can tell at a glance exactly which
package has provided a given operation. Even in a class project there may be
several packages "WITH-ed" by a program, including standard Ada
packages, compiler-provided packages like
Ada.Numerics.Elementary_Functions, and packages you write
yourself
or are supplied by your teacher. Qualification makes it easy to see, for
debugging purposes or for enhancing your program at a later date, just which
operations came from which packages. We will discuss the second advantage in a
moment.
Ada provides a method for avoiding the need to qualify all references to
package-provided operations. This is called the USE clause, looks
just like a context clause, and (in this book) is written at the top of a
program unit along with the context clause. For example,
WITH Ada.Numerics.Elementary_Functions;
USE Ada.Numerics.Elementary_Functions;
might
appear at the top of a program. If a USE clause is present,
qualifying the package references is no longer required, although it is
certainly still permitted. Given a USE clause, the two statements
Y := Ada.Numerics.Elementary_Functions.Sqrt( X);
Y := Sqrt( X);
are
both permitted and have the same meaning. An advantage of the USE
clause is that expressions can be somewhat more compactly written, but of
course the information about just which package provided the Sqrt
operation is lost to the reader.
This direct information about which package provides which operation is also
lost to the compiler if a USE clause is present. This means that
the compiler, in translating an unqualified reference to a package operation,
must search its tables for all the packages mentioned by
USE clauses, and this is a tedious and somewhat time-consuming
task for the compiler. This is the second advantage of qualified reference: It
makes the compiler's job a bit easier.
PROGRAM STYLE
The Proper Use of USE Clauses
USE clause has certain advantages, but
qualification of all references also has advantages. Many experienced Ada
professionals believe that the advantages of qualification outweigh those of
USE clauses, and we tend to agree. Generally we will avoid
USE clauses, continuing the style with which we have begun the
book.
There are certain circumstances, such as the math library, in which the
names of the operations are so obvious, and relate so closely to everyday
mathematics, that the more compact expression notation is desirable, and in
such cases we will write a USE clause. We will also write a
USE clause in cases where the name of the package is not standard
and may vary from compiler to compiler.
As in the case of type conversions, moderation is a virtue in the use of
USE clauses, and we advocate careful, case-by-case consideration
of whether the USE or the qualified reference is more advantageous.
Example 7.11
The function Sqrt (square root) can be used to compute the
roots of a quadratic equation in X of the form
aX2 + bX + c = 0
where a, b, and c are type Float. The two
roots are expressed in algebraic form as
__________
-b + \/ b2 - 4ac
Root1 = ------------------
2a
__________
-b - \/ b2 - 4ac
Root2 = ------------------
2a
The Ada implementation is
IF Disc > 0.0 THEN
Root1 := (-b + Sqrt( Disc)) / ( 2.0 * a);
Root2 := (-b - Sqrt( Disc)) / ( 2.0 * a);
END IF;
where
the variable Disc represents the discriminant
(b2 - 4ac) of the equation.
Example 7.12
Program
7.4 draws a sine curve. It uses the Ada function Sin, provided
by the math package Ada.Numerics.Elementary_Functions, which
returns the trigonometric sine of its parameter, an angle expressed in radians.
Program 7.4
Plotting a Sine Curve Using the Math Package
WITH Ada.Text_IO;
WITH Ada.Float_Ada.Text_IO;
WITH Ada.Numerics;
USE Ada.Numerics;
WITH Ada.Numerics.Elementary_Functions;
USE Ada.Numerics.Elementary_Functions;
PROCEDURE Sine_Curve IS
------------------------------------------------------------------------
--| Plots a sine curve.
--| Author: Michael B. Feldman, The George Washington University
--| Last Modified: July 1995
------------------------------------------------------------------------
RadPerDegree : CONSTANT Float := Pi / 180.0;
-- radians per degree
-- Pi in Ada.Numerics
MinAngle : CONSTANT Float := 0.0; -- smallest angle
MaxAngle : CONSTANT Float := 360.0; -- largest angle
PlotWidth : CONSTANT Integer := 40; -- width of plot
PlotHeight : CONSTANT Integer := 20; -- height of plot
StepAngle : CONSTANT Float
:= (MaxAngle-MinAngle) / Float(PlotHeight);
-- change in angle
Star : CONSTANT Character := '*'; -- symbol being plotted
Blank: CONSTANT Character := ' '; -- to "pad" the '*'
SUBTYPE ColumnRange IS Integer RANGE 0..PlotWidth;
Angle : Float; -- angle in degrees
Radian : Float; -- angle in radians
Scale : Float; -- scale factor
Pad : ColumnRange; -- size of blank padding
BEGIN -- Sine_Curve
Ada.Text_IO.Put(Item => " Sine curve plot");
Ada.Text_IO.New_Line(2);
Scale := Float(PlotWidth / 2);
Angle := MinAngle;
WHILE Angle <= MaxAngle LOOP
Radian := Angle * RadPerDegree;
Pad := Natural(Scale * (1.0 + Sin(Radian)));
Ada.Float_Ada.Text_IO.Put
(Item =>Angle, Fore => 4, Aft => 0, Exp => 0);
-- Display blank padding
Ada.Text_IO.Put(Item => Blank);
FOR BlankCount IN 1 .. Pad LOOP
Ada.Text_IO.Put(Item => Blank);
END LOOP;
Ada.Text_IO.Put(Item => Star); -- Plot * in next column
Ada.Float_Ada.Text_IO.Put
(Item =>Sin(Radian), Fore => 6, Aft => 6, Exp => 0);
Ada.Text_IO.New_Line;
Angle := Angle + StepAngle;
END LOOP;
END Sine_Curve;
Sample
Run
Sine curve plot
0.0 * 0.000000
18.0 * 0.309017
36.0 * 0.587785
54.0 * 0.809017
72.0 * 0.951057
90.0 * 1.000000
108.0 * 0.951057
126.0 * 0.809017
144.0 * 0.587785
162.0 * 0.309017
180.0 * -0.000000
198.0 * -0.309017
216.0 * -0.587785
234.0 * -0.809017
252.0 * -0.951056
270.0 * -1.000000
288.0 * -0.951056
306.0 * -0.809017
324.0 * -0.587785
342.0 * -0.309017
360.0 * 0.000000
Because
degrees are a more intuitive way to represent angles, the outer
FOR loop is executed for values of Angle equal to 0,
18, 36, ... , 360 degrees. This requires a conversion to radians in order to
give Sin a sensible parameter value. We need a conversion
constant, RadPerDegree, which is the value pi/180.0.
We use the value of Pi provided by
Ada.Numerics.
Now let us see how to plot the curve. For each Angle, the first
assignment statement below
Radian := Angle * RadPerDegree; Pad := Natural( Scale * (1.0 + Sin( Radian));computes the number of radians corresponding to
Angle. Then the variable
Pad is assigned a value based on Sin(Radian). This
value increases from 0 when Sin(Radian) is -1.0 to twice the value
of Scale when Sin(Radian) is 1.0. Pad,
the limit variable in the inner FOR loop, determines how many
blanks precede each character * displayed on the screen. In this
way, the position of each * displayed represents the sine of the
current angle. The angle is displayed at the left end of each line; the sine
value is also displayed as a floating-point number after each *.
PROGRAM STYLE
Checking Boundary Values
Pad ranges from 0 to
twice Scale as the sine value goes from -1.0 to 1.0. It is always
a good idea to check the accuracy of these assumptions; this usually can be
done by checking the boundaries of the range as shown
below.
Sin( Radian) is -1.0, Pad is Natural( Scale * ( 1.0 + (-1.0)))
Pad is Natural( 20.0 * 0.0)
Pad is Natural( 0.0) = 0
Sin( Radian) is +1.0, Pad is Natural( Scale * ( 1.0 + 1.0))
Pad is Natural( 20.0 * 2.0)
Pad is Natural( 40.0) = 40
It is also a good idea to check the boundary values for all loop control
variables to see that they make sense. For example, the outer loop control
variable, Angle, has an initial value of MinAngle
( 0.0) and a final value of MaxAngle ( 360.0). The inner loop
control variable, BlankCount, has an initial value of 1 and a
final value of Pad.
PROBLEM SPECIFICATION
The math constant e (whose value is the nonterminating decimal
2.71828...) is the base of the natural logarithms. This value is provided by
Ada as Ada.Numerics.E. Suppose we did not have the numerics
package; develop a program that will compute the value of e. The user
will supply the desired number of decimal places of accuracy.
ANALYSIS
There are a number of mathematical quantities that can be represented using a series approximation, where a series is represented by a summation of an infinite number of terms. For example, e can be determined by evaluating the expression
1 + 1 /1! + 1/2! + 1/3! + ... + 1/n! + ...
where n! is the factorial of n as defined below:
0! = 1
n! = n × (n - 1)! (for n >= 1)
Notice that this is just a different, equivalent, way of defining the same
Factorial that we defined in Section
5.8. Instead of calculating
the factorial for each term in the series, we shall use a different method as
outlined below.
We can get an approximation to the value of e by summing the series for a finite value of n. Obviously, the larger the value of n we use, the more accurate will be the computed result. This expression can be represented using summation notation as
n
---
\ 1/i!
/
---
i=0
where the first term is obtained by substituting 0 for i (1/0! is 1/1),
the second term is obtained by substituting 1 for i (1/1!), etc.
To get an approximation to the desired accuracy, we use a successive
approximations method. Suppose the number of decimal places is given by
Places. Start with a single term, then add terms until two
successive approximations differ by no more than
1/10Places. This last quantity is usually called
epsilon, which in mathematics is used to mean a very small interval. For
example, if we desire 6 decimal places, epsilon is 0.000001 = 1/106.
DESIGN
A general loop can be used to implement the formula above easily. The data requirements and algorithm follow.
Data Requirements
Problem Inputs:
the number of decimal places (Places : Positive)
Problem Outputs:
the approximate value of e (e: Float)
Program Variables:
i, to produce the ith term (i: Natural)
the ith term of the series (ithTerm : Float)
the previously computed estimate of e (eOld: Float)
the number of terms in the series (n: Positive)
the desired accuracy (Epsilon: Float)
DESIGN
The algorithm for this case study follows.
Algorithm
1. Prompt user the value of Places
2. Set Epsilon to
1.0/10.0Places
3. Initialize e to 1.0
4. Initialize the ith term to 1.0
5. LOOP
6. Save previous e in eOld
7 Increment i
8. Compute the ith term in the series
9. Add the ith term to e
10. EXIT WHEN e and eOld differ by
no more than Epsilon
END LOOP;
11. Display the approximate value of e
IMPLEMENTATION
Program 7.5 implements this algorithm.
Program
7.5
Estimating e by Successive Approximations
WITH Ada.Text_IO;
WITH Ada.Integer_Ada.Text_IO;
WITH Ada.Float_Ada.Text_IO;
PROCEDURE Estimate_e IS
------------------------------------------------------------------------
--| Computes the value of e by a series approximation.
--| Number of places of accuracy is specified by user input.
--| Author: Michael B. Feldman, The George Washington University
--| Last Modified: July 1995
------------------------------------------------------------------------
Places : Positive; -- Input - number of decimal places of accuracy
e : Float; -- Output - the value being approximated
eOld : Float; -- previous estimate
i : Natural; -- to produce the i-th term
ithTerm : Float; -- ith term in series
n : Positive; -- number of terms in series
Epsilon : Float; -- desired difference between successive tries
BEGIN -- Estimate_e
Ada.Text_IO.Put(Item => "Enter desired number of decimal places > ");
Ada.Integer_Ada.Text_IO.Get(Item => Places);
Epsilon := 1.0 / (10.0 ** Places);
Ada.Text_IO.Put(Item => "Number of Terms Approximate Value of e");
Ada.Text_IO.New_Line;
Ada.Text_IO.Put(Item => "--------------- ----------------------");
Ada.Text_IO.New_Line;
-- Compute each term and add it to the accumulating sum.
e := 1.0; -- initial sum
ithTerm := 1.0; -- first term
i := 0;
LOOP -- and quit when desired accuracy is achieved
eOld := e; -- save previous approximation
i := i + 1;
ithTerm := ithTerm / Float(i);
e := e + ithTerm; -- find new value
Ada.Integer_Ada.Text_IO.Put(Item => i, Width => 9);
Ada.Float_Ada.Text_IO.Put
(Item => e, Fore => 10, Aft => Places+2, Exp => 0);
Ada.Text_IO.New_Line;
EXIT WHEN ABS (e - eOld) <= Epsilon;
END LOOP;
END Estimate_e;
Sample
Run
Enter desired number of decimal places > 10
Number of Terms Approximate Value of e
--------------- ----------------------------
1 2.000000000000
2 2.500000000000
3 2.666666746140
4 2.708333492279
5 2.716666936874
6 2.718055725098
7 2.718254089355
8 2.718278884888
9 2.718281745911
10 2.718281984329
11 2.718281984329
Inside
the loop, the statement
ithTerm := ithTerm / Float( i);computes the value of the ith term in the series by dividing the previous term by the type
Float representation of the variable
i. The formula
(1 / (i - 1)!) / i = 1 / (i x (i - 1)!) = 1 / i!
shows that this division does indeed produce the next term in the series.
Because 0! is 1, ithTerm must be intialized to 1.0. The statement
e = e + ithTerm;adds the new value of
ithTerm to the sum being
accumulated in e. Trace the execution of this loop to satisfy
yourself that ithTerm takes on the values 1/1!, 1/2!, 1/3!, and so
on, during successive loop iterations.
One of the problems in processing floating-point numbers is that there is sometimes an error in representing floating-point data. Just as there are certain numbers that cannot be represented exactly in the decimal number system (e.g., the fraction 1/3 is 0.333333...), so there are numbers that cannot be represented exactly in floating-point form. The representational error will depend on the number of binary digits (bits) used in the mantissa: The more bits there are, the smaller the error.
The number 0.1 is an example of a real number that has a representational error. The effect of a small error is often magnified through repeated computations. Therefore, the result of adding 0.1 ten times is not exactly 1.0, so the loop below may fail to terminate on some computers:
Trial := 0.0; WHILE Trial /= 1.0 LOOP ... Trial := Trial + 0.1; END LOOP;
If the loop repetition test is changed to Trial <
1.0, the loop may execute ten times on one computer and eleven times on
another. For this reason, it is best to use integer values--which are always
exact--whenever possible in loop repetition tests.
Other problems occur when manipulating very large and very small real
numbers. In adding a large number and a small number, the larger number may
"cancel out" the smaller number (a cancellation error). If
X is much larger than Y, X + Y and
X may have the same value (e.g., 1000.0 + 0.0001234 is equal to
1000.0 on some computers).
For this reason, you can sometimes obtain more accurate results by carefuly
selecting the order in which computations are performed. For example, in
computing the value of e in the preceding case study, the terms of
the series
1 + 1/1! + 1/2! + ... + 1/n!
were generated in left-to-right order and added to a sum being accumulated
in e. When n is large, the value of 1/n! is very
small, so the effect of adding a very small term to a sum that is larger than
2.0 may be lost. If the terms were generated and summed in right-to left-order
instead, the computation result might be more accurate.
If two very small numbers are multiplied, the result may be too small to be
represented accurately and will become zero. This is called arithmetic
underflow. Similarly, if two very large numbers are multiplied, the result
may be too large to be represented. This is called arithmetic overflow
and, in Ada, causes Constraint_Error to be raised. Arithmetic
underflow and overflow can also occur when processing very large and small
integer values.
a. _______ \/ U + V × W2 b. logn( XY) c. _______ \/ X - Y2 d. | XY - W/Z|
Sqrt(ABS(Integer(-15.8)))
Ada.Numerics.E. Look up this value and check the results of
Program 7.5 against it for several different numbers of decimal places of
accuracy.
X to the nearest two decimal places. (Hint: You have to
multiply by 100.0 before rounding.)
1 + x + x2/2! + x3/3! + ... + xn/n! + ...
Write a program to compute and print the value of this series for any
x and any n. Compare the result to Exp( x) (available
in the Ada math library) and print a message O.K. or Not
O.K., depending on whether the difference between these results exceeds
0.001. How many terms--that is, what value of n--seem to provide good
results without making the computation take too many steps?
Copyright © 1996 by Addison-Wesley Publishing Company, Inc.