This is the mail archive of the newlib@sourceware.org mailing list for the newlib project.


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]
Other format: [Raw text]

erf(x)+erfc(x)=1 precision problem


While testing the new clisp-cvs on cygwin I came to the following
erf/erfc precision problems with the double definitions.

There are some minor deviations from the glibc2 (linux) results, but a
bigger problem is that erf(-1) + erfc(-1) is unequal to (double) 1.0.
The error is 1.1102230246251565e-16
This looks like just one digit, but it is the only significant digit
deviation in the erf(x) + erfc(x)
test  between [-10, 10].
erfc(27) leads to floating point underflow too, but this is probably acceptable.

Sorry, I'm not qualified at all to offer a patch, just the link to sun's
free erf/erfc implementation [1], which lcc also uses.

Results in detail:
CORRECT is glibc, CLISP is cygwin newlib gcc-3.4.4 on a i686 with HW FPU.

Form: (ERFC -1)
CORRECT: 1.842700792949715d0
CLISP  : 1.8427007929497148d0


Form: (ERFC 3)
CORRECT: 2.2090496998585438d-5
CLISP  : 2.209049699858544d-5


Form: (ERFC 5)
CORRECT: 1.5374597944280351d-12
CLISP  : 1.537459794428035d-12


Form: (ERFC 20)
CORRECT: 5.3958656116079005d-176
CLISP  : 5.395865611607901d-176


Form: (ERFC 22)
CORRECT: 1.6219058609334724d-212
CLISP  : 1.6219058609334726d-212


Form: (ERFC 23)
CORRECT: 4.441265948088058d-232
CLISP  : 4.441265948088057d-232


Form: (ERFC 25)
CORRECT: 8.300172571196522d-274
CLISP  : 8.300172571196523d-274


Form: (LOOP :FOR I :FROM -10 :TO 10 :ALWAYS (= 1 (+ (ERF I) (ERFC I))))
CORRECT: T
CLISP  : NIL


Form: (Y0 1.0)
CORRECT: 0.08825696421567698d0
CLISP  : 0.08825696421567696d0


Form: (Y0 10.0)
CORRECT: 0.055671167283599395d0
CLISP  : 0.05567116728359938d0


Form: (Y1 10.0)
CORRECT: 0.2490154242069538d0
CLISP  : 0.24901542420695386d0


Form: (YN 2 10.0)
CORRECT: -0.0058680824422086345d0
CLISP  : -0.005868082442208612d0

[6]> (loop :for i :from -10 :to 10 :always (print (abs (- (+ (os:erf
i) (os:erfc i)) 1))))

0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
1.1102230246251565d-16
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0
0.0d0

That means that for double precision between -10 and +10 only -1 is imprecise.
glibc is correct here.

(loop :for i :from -100 :to 100 :always (print (list i (abs (- (+
(os:erf i) (os:erfc i)) 1)))))
says that also only -1 has this problem, -100 to +26 is acceptable,
and that erfc(27) leads to floating point underflow.


[1] http://www.velocityreviews.com/forums/t317358-erf-function-in-c.html
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
-- 
Reini Urban
http://phpwiki.org/ http://murbreak.at/


Index Nav: [Date Index] [Subject Index] [Author Index] [Thread Index]
Message Nav: [Date Prev] [Date Next] [Thread Prev] [Thread Next]