GNU Octave  3.8.0
A high-level interpreted language, primarily intended for numerical computations, mostly compatible with Matlab
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
derf.f
Go to the documentation of this file.
1 *DECK DERF
2  DOUBLE PRECISION FUNCTION derf (X)
3 C***BEGIN PROLOGUE DERF
4 C***PURPOSE Compute the error function.
5 C***LIBRARY SLATEC (FNLIB)
6 C***CATEGORY C8A, L5A1E
7 C***TYPE DOUBLE PRECISION (ERF-S, DERF-D)
8 C***KEYWORDS ERF, ERROR FUNCTION, FNLIB, SPECIAL FUNCTIONS
9 C***AUTHOR Fullerton, W., (LANL)
10 C***DESCRIPTION
11 C
12 C DERF(X) calculates the double precision error function for double
13 C precision argument X.
14 C
15 C Series for ERF on the interval 0. to 1.00000E+00
16 C with weighted error 1.28E-32
17 C log weighted error 31.89
18 C significant figures required 31.05
19 C decimal places required 32.55
20 C
21 C***REFERENCES (NONE)
22 C***ROUTINES CALLED D1MACH, DCSEVL, DERFC, INITDS
23 C***REVISION HISTORY (YYMMDD)
24 C 770701 DATE WRITTEN
25 C 890531 Changed all specific intrinsics to generic. (WRB)
26 C 890531 REVISION DATE from Version 3.2
27 C 891214 Prologue converted to Version 4.0 format. (BAB)
28 C 900727 Added EXTERNAL statement. (WRB)
29 C 920618 Removed space from variable name. (RWC, WRB)
30 C***END PROLOGUE DERF
31  DOUBLE PRECISION x, erfcs(21), sqeps, sqrtpi, xbig, y, d1mach,
32  1 dcsevl, derfc
33  LOGICAL first
34  EXTERNAL derfc
35  SAVE erfcs, sqrtpi, nterf, xbig, sqeps, first
36  DATA erfcs( 1) / -.4904612123 4691808039 9845440333 76 d-1 /
37  DATA erfcs( 2) / -.1422612051 0371364237 8247418996 31 d+0 /
38  DATA erfcs( 3) / +.1003558218 7599795575 7546767129 33 d-1 /
39  DATA erfcs( 4) / -.5768764699 7674847650 8270255091 67 d-3 /
40  DATA erfcs( 5) / +.2741993125 2196061034 4221607914 71 d-4 /
41  DATA erfcs( 6) / -.1104317550 7344507604 1353812959 05 d-5 /
42  DATA erfcs( 7) / +.3848875542 0345036949 9613114981 74 d-7 /
43  DATA erfcs( 8) / -.1180858253 3875466969 6317518015 81 d-8 /
44  DATA erfcs( 9) / +.3233421582 6050909646 4029309533 54 d-10 /
45  DATA erfcs( 10) / -.7991015947 0045487581 6073747085 95 d-12 /
46  DATA erfcs( 11) / +.1799072511 3961455611 9672454866 34 d-13 /
47  DATA erfcs( 12) / -.3718635487 8186926382 3168282094 93 d-15 /
48  DATA erfcs( 13) / +.7103599003 7142529711 6899083946 66 d-17 /
49  DATA erfcs( 14) / -.1261245511 9155225832 4954248533 33 d-18 /
50  DATA erfcs( 15) / +.2091640694 1769294369 1705002666 66 d-20 /
51  DATA erfcs( 16) / -.3253973102 9314072982 3641600000 00 d-22 /
52  DATA erfcs( 17) / +.4766867209 7976748332 3733333333 33 d-24 /
53  DATA erfcs( 18) / -.6598012078 2851343155 1999999999 99 d-26 /
54  DATA erfcs( 19) / +.8655011469 9637626197 3333333333 33 d-28 /
55  DATA erfcs( 20) / -.1078892517 7498064213 3333333333 33 d-29 /
56  DATA erfcs( 21) / +.1281188399 3017002666 6666666666 66 d-31 /
57  DATA sqrtpi / 1.772453850 9055160272 9816748334 115d0 /
58  DATA first /.true./
59 C***FIRST EXECUTABLE STATEMENT DERF
60  IF (first) THEN
61  nterf = initds(erfcs, 21, 0.1*REAL(d1mach(3)))
62  xbig = sqrt(-log(sqrtpi*d1mach(3)))
63  sqeps = sqrt(2.0d0*d1mach(3))
64  ENDIF
65  first = .false.
66 C
67  y = abs(x)
68  IF (y.GT.1.d0) go to 20
69 C
70 C ERF(X) = 1.0 - ERFC(X) FOR -1.0 .LE. X .LE. 1.0
71 C
72  IF (y.LE.sqeps) derf = 2.0d0*x/sqrtpi
73  IF (y.GT.sqeps) derf = x*(1.0d0 + dcsevl(2.d0*x*x-1.d0,
74  1 erfcs, nterf))
75  RETURN
76 C
77 C ERF(X) = 1.0 - ERFC(X) FOR ABS(X) .GT. 1.0
78 C
79  20 IF (y.LE.xbig) derf = sign(1.0d0-derfc(y), x)
80  IF (y.GT.xbig) derf = sign(1.0d0, x)
81 C
82  RETURN
83  END