ANS FORTH CELEFUNT TEST RESULTS FOR ZEXP(Z) zexp.fs vs. 0.9.3 dnw 19Jan05 These results were obtained with pfe 0.33.58 compiled with gcc 3.3.3 on my MacOS X 10.3.7, PowerMac dual G4 system. Leaving aside the pfe C-primitive tests, Gforth 0.6.2 gives the same results as pfe on the same system, except for the way NaN and Inf are treated. Results from Cody's original Fortran code are included for reference. It was compiled with GNU Fortran (GCC) 3.4 20031015 (experimental) under MacOS X. The "normal" version calls complex.fs vs. 0.8.2, Julian Noble's original code for ZEXP. The "kahan" version calls complex-kahan.fs vs. 0.8.7, which uses the Kahan algorithm without treatment of IEEE 754 underflow/overflow exceptions. The "pfe prim" version calls the pfe complex module vs. 0.8.7, which codes the Kahan algorithm in C with treatment of underflow/overflow exceptions. This module is part of the pfe 0.33.xx series distribution. Each test has three sections corresponding to the magnitude of the relative complex difference of the compared functions, the relative error in the real part, and the relative error in the imaginary part. N= = number out of 2000 random trials that are equal MRE = maximum relative error RMS = rms relative error ULP = estimated units in last place (loss of base 2 significant digits) Floating point numbers have 53 significant bits, with at most 7 significant decimal digits printed out below. The complex numbers following "box" are opposite corners of the box in the complex plane from which random values of z are chosen. exp[z] vs. xp[z-del]*exp[del], del = 1/16 + i*1/16 box 0+i6.25E-02, 1+i1 N= MRE ULP RMS ULP vector g77 506 3.86E-16 1.80 1.24E-16 0.16 normal 482 4.31E-16 1.96 1.24E-16 0.16 kahan 482 4.31E-16 1.96 1.24E-16 0.16 pfe prim 486 4.31E-16 1.96 1.24E-16 0.16 real g77 1010 4.40E-16 1.99 1.25E-16 0.17 normal 995 4.43E-16 2.00 1.26E-16 0.18 kahan 995 4.43E-16 2.00 1.26E-16 0.18 pfe prim 1003 4.43E-16 2.00 1.26E-16 0.18 imaginary g77 955 4.22E-16 1.93 1.28E-16 0.20 normal 945 4.22E-16 1.93 1.28E-16 0.21 kahan 945 4.22E-16 1.93 1.28E-16 0.21 pfe prim 946 4.22E-16 1.93 1.28E-16 0.21 exp[z] vs. exp[z-del]*exp[del], del = 1/16 + i*1/16 box 1.625E+00+i1.625E+00, 3+i3 N= MRE ULP RMS ULP vector g77 513 3.92E-16 1.82 1.23E-16 0.15 normal 524 3.92E-16 1.82 1.22E-16 0.13 kahan 524 3.92E-16 1.82 1.22E-16 0.13 pfe prim 524 3.92E-16 1.82 1.22E-16 0.14 real g77 1017 4.30E-16 1.95 1.24E-16 0.16 normal 1036 5.91E-16 2.41 1.24E-16 0.16 kahan 1036 5.91E-16 2.41 1.24E-16 0.16 pfe prim 1024 4.08E-16 1.88 1.24E-16 0.16 imaginary g77 982 4.52E-16 2.03 1.28E-16 0.20 normal 987 4.29E-16 1.95 1.26E-16 0.18 kahan 987 4.29E-16 1.95 1.26E-16 0.18 pfe prim 984 4.29E-16 1.95 1.26E-16 0.19 exp[z] vs. exp[z-del]*exp[del], del = 1/16 + i*1/16 box 16+i16, 17+i17 N= MRE ULP RMS ULP vector g77 465 4.27E-16 1.94 1.28E-16 0.21 normal 455 4.27E-16 1.94 1.29E-16 0.22 kahan 455 4.27E-16 1.94 1.29E-16 0.22 pfe prim 457 4.27E-16 1.94 1.29E-16 0.22 real g77 967 4.84E-16 2.12 1.30E-16 0.23 normal 957 4.83E-16 2.12 1.31E-16 0.24 kahan 957 4.83E-16 2.12 1.31E-16 0.24 pfe prim 957 4.83E-16 2.12 1.30E-16 0.23 imaginary g77 959 4.57E-16 2.04 1.31E-16 0.24 normal 927 4.57E-16 2.04 1.32E-16 0.25 kahan 927 4.57E-16 2.04 1.32E-16 0.25 pfe prim 929 4.57E-16 2.04 1.32E-16 0.25 g77 Z EXP[Z]*EXP[-Z] - 1 7.832594E-01 + i 1.907325E+00 -1.110223E-16 + i 0.000000E+00 4.153622E-01 + i 1.920231E+00 0.000000E+00 + i 0.000000E+00 2.862098E-02 + i 1.577621E+00 -1.110223E-16 - i 8.673617E-19 1.202422E+00 + i 3.025713E-01 0.000000E+00 + i 0.000000E+00 1.821380E+00 + i 1.672224E+00 -1.110223E-16 - i 1.387779E-17 normal Z EXP[Z]*EXP[-Z] - 1 7.832594E-01 + i 1.907325E+00 0.000000E+00 + i 2.775558E-17 4.153622E-01 + i 1.920231E+00 0.000000E+00 + i 0.000000E+00 2.862098E-02 + i 1.577621E+00 -1.110223E-16 + i 0.000000E+00 1.202422E+00 + i 3.025713E-01 2.220446E-16 + i 2.220446E-16 1.821380E+00 + i 1.672224E+00 -2.220446E-16 - i 1.387779E-17 kahan Z EXP[Z]*EXP[-Z] - 1 7.832594E-01 + i 1.907325E+00 0.000000E+00 + i 2.775558E-17 4.153622E-01 + i 1.920231E+00 0.000000E+00 + i 0.000000E+00 2.862098E-02 + i 1.577621E+00 -1.110223E-16 + i 0.000000E+00 1.202422E+00 + i 3.025713E-01 2.220446E-16 + i 2.220446E-16 1.821380E+00 + i 1.672224E+00 -2.220446E-16 - i 1.387779E-17 pfe prim Z EXP[Z]*EXP[-Z] - 1 7.832594E-01 + i 1.907325E+00 -1.110223E-16 + i 1.863435E-17 4.153622E-01 + i 1.920231E+00 0.000000E+00 - i 3.649306E-18 2.862098E-02 + i 1.577621E+00 0.000000E+00 - i 8.068020E-20 1.202422E+00 + i 3.025713E-01 2.220446E-16 + i 2.697428E-16 1.821380E+00 + i 1.672224E+00 -1.110223E-16 - i 1.166000E-17 EXP[0+i0] - 1 = g77 0.000000E+00 + i 0.000000E+00 normal 0.000000E+00 + i 0.000000E+00 kahan 0.000000E+00 + i 0.000000E+00 pfe prim 0.000000E+00 + i 0.000000E+00 [REAL[EXP[0+i10]] - COS[10]]/COS[10] = g77 0.000000E+00 normal -0.000000E+00 kahan -0.000000E+00 pfe prim -0.000000E+00 [IMAG[EXP[0+i10]] - SIN[10]]/SIN[10] = g77 0.000000E+00 normal -0.000000E+00 kahan -0.000000E+00 pfe prim -0.000000E+00 [EXP[EXP[10+i0]] - EXP[10]]/EXP[10] = g77 0.000000E+00 normal 0.000000E+00 kahan 0.000000E+00 pfe prim 0.000000E+00 Tests near extreme arguments: The real part of the first result below may underflow: EXP[ -7.08000E+02 + i 1.00000E+00 ] = g77 1.78708E-308 + i 2.78321E-308 normal 1.78708E-308 + i 2.78321E-308 kahan 1.78708E-308 + i 2.78321E-308 pfe prim 1.78708E-308 + i 2.78321E-308 EXP[ 7.09000E+02 + i 1.00000E+00 ] = g77 4.44042E+307 + i 6.91555E+307 normal 4.44042E+307 + i 6.91555E+307 kahan 4.44042E+307 + i 6.91555E+307 pfe prim 4.44042E+307 + i 6.91555E+307 There is an argument reduction error if the following are not about the same: g77 EXP[ 3.54500E+02 + i 1.00000E+01 ] = -7.60664E+153 - i 4.93185E+153 EXP[ 1.77250E+02 + i 5.00000E+00 ]^2 = -7.60664E+153 - i 4.93185E+153 normal EXP[ 3.54500E+02 + i 1.00000E+01 ] = -7.60664E+153 - i 4.93185E+153 EXP[ 1.77250E+02 + i 5.00000E+00 ]^2 = -7.60664E+153 - i 4.93185E+153 kahan EXP[ 3.54500E+02 + i 1.00000E+01 ] = -7.60664E+153 - i 4.93185E+153 EXP[ 1.77250E+02 + i 5.00000E+00 ]^2 = -7.60664E+153 - i 4.93185E+153 pfe prim EXP[ 3.54500E+02 + i 1.00000E+01 ] = -7.60664E+153 - i 4.93185E+153 EXP[ 1.77250E+02 + i 5.00000E+00 ]^2 = -7.60664E+153 - i 4.93185E+153 Test of error returns: The real part of the argument below is so negative that the real exponential should underflow to zero, giving 0+i0: EXP[ -6.704E+153 + i 1.000E+00 ] = g77 0.000E+00 + i 0.000E+00 normal 0.000E+00 + i 0.000E+00 kahan 0.000E+00 + i 0.000E+00 pfe prim 0.000E+00 + i 0.000E+00 The imaginary part of the argument below is too large for the sine and cosine computations to make sense: EXP[ 1.000E+01 + i 1.341E+154 ] = g77 2.202E+04 + i 5.339E+02 normal 2.202E+04 + i 5.339E+02 kahan 2.202E+04 + i 5.339E+02 pfe prim 2.202E+04 + i 5.339E+02 The real part of the argument below is so large that the real exponential should overflow, giving Inf-i*Inf and/or an error: EXP[ 6.704E+153 - i 1.000E+00 ] = g77 INF - i INF normal INF - i INF kahan INF - i INF pfe prim INF - i INF