XProfan

WindowTitle "   Asym-Erf(), die Gauss-Fehlerfunktion in asymptotischer Näherung "
' Quelle: https://jean-pierre.moreau.pagesperso-orange.fr/Basic/asymerf_bas.txt
' (D) Demo 2016/11 transponiert nach XProfan 11.2a by P. Specht, Vienna / Austria
' Keine wie auch immer geartete Gewähr! No warranty whatsoever!
'****************************************************
'*         Program to demonstrate ASYMERF           *
'* ------------------------------------------------ *
'* Reference: Di base 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
WindowStyle 24:CLS:font 2:Declare x!,y!,e!,n&:Set("decimals",17)
Tests:
Print "\n\n  Die Gauss-Fehlerfunktion Erf(x) per il valore x = ";:Input x!
'Gosub "S1000"
S1000
Print
Print "\n  Erf(x) per il valore x= ";x!;"  beträgt: ";format$("0.######",y!)
PRINT "  mit einer Fehleabschätzung zu ";format$("##0.##########",e!)
PRINT
PRINT "          ";n&;" Terme wurde berücksichtigt.";
WaitInput
Goto "Tests"
'***********************************************************
'* 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.                                 *
'***********************************************************

Proc S1000

    Declare c1!,c2!,y1!
    n& = 1 : y! = 1 : c2! = 1 / (2 * Sqr(x!))
    'G1100:

    Repeat

        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
        Case ABS(c2!) > Abs(c1!) : BREAK':= Goto "G1200"
        ' Continue summation

    Until 0':=  Goto "G1100"

    ' G1200:
    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.prf