This file:
http://ftp.aset.psu.edu/pub/ger/fortran/hdk/FloatingPointHDK.txt
Advice about floating point by Herman D.Knoble, 9 March 2004.
One technique to illustrate to students how significant digits can affect
computations is to use examples which break down very quickly
and drastically. (Worse yet are ones which break down more subtly!)
Back in '93 John R. Ehrman (then at IBM Santa Teresa Labs) and
I did a back to back presentation at Share 75 (See Proceedings of
Share 75) on this topic. John presented an excellent "Floating-point
Numbers and Arithmetic Tutorial showing the various models and
implementations and how arithmetic really breaks down. My presentation,
there was "Reality of REAL Arithmetic's Really Wrong Results". The
following example is from that presentation and has been for
us, a fairly good teaching example that you don't need to
limit how many digits participate in the model to make a point.
program Example1
!
!=== Example 1: Numeric breakdown for IEEE Arithmetic
! (or worse yet for IBM 360 Floating-point arithmetic).
! For example try: A=1000000., X=0.0000001
! and see P computed as: 4534650. in Single Precision
! or: 4768372. in Double Precision.
! Note that P is algebraically equal to 1.0 for all A and X /= 0.
!
! Penn State University Center for Academic Computing
! H. D. Knoble, August 1979.
REAL :: P,A,X
DO
WRITE(*,*) 'Please enter A and X (or press Enter to Quit):'
WRITE(*,*) '(Try values 10000., 0.0001, then 0.0001, 10000.)'
READ(*,*,END=99) A,X
IF (A == 0 .OR. X == 0) GOTO 99
P=((A+X)**2 - A**2 - 2.*A*X)/X**2
WRITE(*,100) P,A,X
100 FORMAT(' P=',G14.7,' A=',G14.7,' X=',G14.7)
ENDDO
99 WRITE(*,*) '** Example 1 Ended **'
End Program Example1
Interesting enough, the bigger the "spread" of A and X in general the bigger
the "wrong" answer. However essentially a decent result is computed even
with a very large spread between A and X if the rolls of A and X are
reversed (X large, A, small). The problem here of course is subtracting nearly
equal quantities which are formed by converting decimal numbers that have
infinite binary expansions; when subtracting nearly equal quantities, the
most significant digits cancel during the subtraction, thus shifting the "noise"
toward the most significant digits. Then dividing by X**2 is like multiplying
the "noise" by a very large number (scaling it up). "Really wrong results"
here illustrate with very little arithmetic so-called "significance loss".
(This presentation shows that decimal fractions like 1/10th have an infinite
binary expansion,like: 0.00011001100110011001100 .
A second part to a student exercise studying such an example (like the one
above or one using the so-called "calculator formula" for Standard Deviation
which subtracts two sums that will be nearly equal for a small Standard Deviation)
is to make a set of "recommendations" to help avoid numerical computational
difficulties. Working with some undergrad students, working as a team,
back in 1993 came up with:
1) Use Proven (refereed or commercial) library programs that are as robust
(give the right result or an error message if not) as possible. For example,
IMSL, NAG mathematical statistical libraries. Avoid blindly taking numerical
code from textbooks and magazines unless they are refereed or presented
by a proven expert in numerical computations.
2) Study (know) the possible domain of values that each program input
quantity will take on and test each algorithm thoroughly, especially on
the limits of that domain. (This of course is an art in itself).
3) Use at least Double Precision Arithmetic, including constants (like
0.1D0). In more critical cases use multiple or quadruple precision arithmetic:
http://www.nersc.gov/~dhbailey/ and
http://www.ozemail.com.au/~milleraj/ respectively.
4) For a given mantissa precision (computer) do not input or output numbers
in ASCII (or EBCDIC) format with more decimal digits than can be handled by the
conversion mapping algorithms for that precision (computer) unless you
document that that's what you are doing.
5) Use relative comparisons in place of absolute comparisons of real
quantities. (Also note that you can also use "fuzzy" comparisons. See
http://ftp.aset.psu.edu/pub/ger/fortran/hdk/eps.f90"
6) Avoid computations which are algebraically undefined (or singular) or which
develop results whose magnitudes approach or exceed machine limits.E.g.,
Blindly setting results of underflows to zero can yield "Really Wrong Results".
7) Avoid mixing Fortran INTEGER and REAL quantities without knowing how they
will be compiled. E.g., X=X+I**(-J) does not increment REAL typed X except when
INTEGER typed I=1 or INTEGER typed J=0 and I < 0.
8) Do not program a numerical model beyond your competence to either avoid
significance loss or monitor it. Be especially careful when less than a pre-specified
number of mantissa digits participate in a sum or when all but a few mantissa
digits cancel during a subtraction. You are responsible for disasters caused
by other people using your code!
I'm sure there are more points to add here, especially along the lines of
ethics implied by the last sentence in (8) above.