Next: PASCAL version Up: Floating point precision Previous: C version

FORTRAN version


 Is this a program restart after failure (1)
 or a start from scratch (0) ?
 A  Paranoid  Program  to  Diagnose  Floating-point Arithmetic
          ... Double-Precision Version  ...
 Lest this program stop prematurely, i.e. before displaying
   "End of Test"   
 try to persuade the computer NOT to terminate execution
 whenever an error such as Over/Underflow or Division by
 Zero occurs, but rather to persevere with a surrogate value
 after, perhaps, displaying some warning.  If persuasion
 avails naught, don't despair but run this program anyway
 to see how many milestones it passes, and then run it
 again.  It should pick up just beyond the error and
 continue.  If it does not, it needs further debugging.

 Users are invited to help debug and augment this program
 so that it will cope with unanticipated and newly found
 compilers and arithmetic pathologies.

 To continue diagnosis, press return.
 Diagnosis resumes after milestone  #    0,    ... page     1


 Please send suggestions and interesting results to
         Richard Karpinski
         Computer Center U-76
         University of California
         San Francisco, CA 94143-0704
         USA

 In doing so, please include the following information:
         Precision:   Double;
         Version: 31 July 1986;
         Computer:

         Compiler:

         Optimization level:

         Other relevant compiler options:


 To continue diagnosis, press return.
 Diagnosis resumes after milestone  #    1,    ... page     2


 BASIC version (C) 1983 by Prof. W. M. Kahan.
 Translated to FORTRAN by T. Quarles and G. Taylor.
 Modified to ANSI 66/ANSI 77 compatible subset by
 Daniel Feenberg and David Gay.
 You may redistribute this program freely if you
 acknowledge the source.


 Running this program should reveal these characteristics:

 b = radix ( 1, 2, 4, 8, 10, 16, 100, 256, or ... ) .
 p = precision, the number of significant  b-digits carried.
 u2 = b/b^p = one ulp (unit in the last place) of 1.000xxx..
 u1 = 1/b^p = one ulp of numbers a little less than 1.0.

 To continue diagnosis, press return.
 Diagnosis resumes after milestone  #    2,    ... page     3

 g1, g2, g3 tell whether adequate guard digits are carried;
 1 = yes, 0 = no;  g1 for mult.,  g2 for div., g3 for subt.
 r1,r2,r3,r4  tell whether arithmetic is rounded or chopped;
 0=chopped, 1=correctly rounded, -1=some other rounding;
 r1 for mult., r2 for div., r3 for add/subt., r4 for sqrt.
 s=1 when a sticky bit is used correctly in rounding; else s=0.
 u0 = an underflow threshold.
 e0 and z0 tell whether underflow is abrupt, gradual or fuzzy
 v = an overflow threshold, roughly.
 v0  tells, roughly, whether  infinity  is represented.
 Comparisons are checked for consistency with subtraction
        and for contamination by pseudo-zeros.
 Sqrt is tested. so is  y^x  for (mostly) integers  x .
 Extra-precise subexpressions are revealed but not yet tested.
 Decimal-binary conversion is not yet tested for accuracy.

 To continue diagnosis, press return.
 Diagnosis resumes after milestone  #    3,    ... page     4

 The program attempts to discriminate among:
     >FLAWs, like lack of a sticky bit, 
     >SERIOUS DEFECTs, like lack of a guard digit, and
     >FAILUREs, like  2+2 = 5 .
 Failures may confound subsequent diagnoses.

 The diagnostic capabilities of this program go beyond an
 earlier program called  "Machar", which can be found at the
 end of the book "Software Manual for the Elementary Functions"
 (1980) by W. J. Cody and W. Waite. Although both programs
 try to discover the radix (b), precision (p) and         
 range (over/underflow thresholds) of the arithmetic, this
 program tries to cope with a wider variety of pathologies
 and to say how well the arithmetic is implemented.
 The program is based upon a conventional radix
 representation for floating-point numbers,
 but also allows for logarithmic encoding (b = 1)
 as used by certain early wang machines.


 To continue diagnosis, press return.


peter@physik3.gwdg.de
Tue Aug 8 13:05:00 GMT+0200 1995