'****************************************************
'*         Program to demonstrate ASYMERF           *
'* ------------------------------------------------ *
'* Reference: BASIC Scientific Subroutines, Vol. II *
'* by F.R. Ruckdeschel, BYTE/McGRAWW-HILL, 1981 [1].*
'* ------------------------------------------------ *
'* SAMPLE RUN:                                      *
'*                                                  *
'* Find the value of ERF(X)=2*Exp(-X*X)/SQRT(PI)    *
'*                                                  *
'* Input X ? 3                                      *
'* ERF(X)= .9999779 with error estimate= -.00000000 *
'* Number of terms evaluated was 10                 *
'*                                                  *
'* Input X ? 4                                      *
'* ERF(X)= 1.0000000 with error estimate= 0.0000000 *
'* Number of terms evaluated was 17                 *
'*                                                  *
'****************************************************
DEFINT I-N
DEFDBL A-H, O-Z
CLS
PRINT
INPUT " Input X: ", x
GOSUB 1000
PRINT
PRINT USING " ERF(X)= #.#######  with error estimate= #.#######."; y; e
PRINT
PRINT " Number of terms evaluated was "; n
PRINT
END

'***********************************************************
'* Asymptotic series expansion of the integral of          *
'* 2 EXP(-X*X)/(X*SQRT(PI)), the normalized error function *
'* (ASYMERF). This program determines the values of the    *
'* above integrand using an asymptotic series which is     *
'* evaluated to the level of maximum accuracy.             *
'* The integral is from 0 to X. The input parameter, X     *
'* must be > 0. The results are returned in Y and Y1,      *
'* with the error measure in E. The number of terms used   *
'* is returned in N. The error is roughly equal to first   *
'* term neglected in the series summation.                 *
'* ------------------------------------------------------- *
'* Reference: A short table of integrals by B.O. Peirce,   *
'* Ginn and Company, 1957.                                 *
'***********************************************************
1000 n = 1: y = 1
c2 = 1 / (2 * x * x)
1100 y = y - c2
n = n + 2: c1 = c2
c2 = -c1 * n / (2 * x * x)
'Test for divergence - The break point is roughly N=X*X
IF ABS(c2) > ABS(c1) THEN GOTO 1200
'Continue summation
GOTO 1100
1200 n = (n + 1) / 2
e = EXP(-x * x) / (x * 1.772453850905516#)
y1 = y * e
y = 1# - y1
e = e * c2
RETURN

'End of file Asymerf.bas