The final fitting error is shown in
fitting/Welcome.htm and SGLITCH/Welcome.htm#_Changes
The function AIGAU.FOR uses on the order of 50 operations to reach double precision accuracy for all values of x. The test ide is taigau.wpj.
The relationship with Errf(x) is examined in Erf.htm .doc
Figure 1 AiGauss on a linear scale
Figure 2 AiGauss on a log scale
This is the first attempt to make aigau. The natural log is accurate to ~15 digits rather than the function itself, and there is a glitch caused by integrating the Gaussian from a set value at approximately z ~ -2.8. This is shown most dramatically in the first fit to the noise level.
The result of this early approximation, AIGAU15.FOR, is compared to the current Aigau.for by taigau.wpj The error for z< -10 is the result of exponentiaing a log ~ -100, 2-3 digits are lost for a real with a 3 digit exponent.
Figure 3 The lower limit errors come from taking the exponentials of logarithms ~500. The Glitch at 3 is from subtraction and approximation errors involved in integrating the Gaussian from specified points. See aigau15\Fitting Brack.htm for more details.
The code is texp.for. It was compiled by
ENTER Z
-708.396
EXP(Z)
2.2260053186166400D-308
1+EXP(Z)
1.0000000000000000
ENTER Z
-708.397
EXP(Z)
0.0000000000000000
1+EXP(Z)
1.0000000000000000
Firstly note that values that are too low, z < 27, are simply taken to be zero. The lower limit values are discussed in
In this range, the integrand after removing exp(-z2) can be expanded and integrated term by term. The cutoff is a last term of 10-16. This is discussed in detail along with examples in
genint\Expanding the integrand.doc
The code bracket in this document is incorporated in AIGAU.FOR as the function brack.
These began as the fit
Careful examination showed that at the left edge of the fit, z=-6, there is a spike in the error to about 10-14. This led to a re-examination of the data to fit as generated in
A bit of manipulation of the end points of the data segments coupled with the use of Gauss Laguerre integration for values of z < -3, enable the fit to be to Brack(z) in the expression
This is shown in
fitting\Welcome.doc htm -- contains final error plot.
The values above 0 are equal to 1 minus the values below zero. This means that there is no point in having a relative accuracy of 15 digits in the value subtracted from 1. The data was generated by an FFT of the time domain version of AiGauss as described in
fitabs\AiGaussPgen\Welcome.htm
The fit is described in
The set of constants are
-0.6931471805599452D0
1.995915713776837D0
-2.582310623341017D0
1.953688092571124D0
-0.9387916412290680D0
0.2892445861593335D0
-0.5400784330829852D-01
0.5078739222313374D-02
-0.1061251766837922D-03
-1.251590673686981D0
0.7695652871670581D0
-0.2680046704303274D0
0.5332127596813570D-01
-0.5079508799796978D-02
0.1060850657795260D-03
ENTER Z
-36
EXP(Z) 2.3195228302435680D-016
1+EXP(Z) 1.0000000000000000
For values of z2 > 36, AiGauss(z) = 1 to the double precision limit.
The log terms compared to were generated by code in gleg.zip. the file with these results is AIGAUL.IN The comparison code used to generate the plot below is in taigau2.wpj
Figure 4 |exp(aigaul.in)-aigau|/aigau
I believe that the bulk of this error comes from the fact that the exp(log) loses accuracy when the log is larger than 1.
Directory of public_html\aigau
<DIR> aigau15 aigau15/Welcome.doc
<DIR> aiglz aiglz/Welcome.doc
<DIR> fitabs fitabs/MPNlfitAigau.doc
<DIR> fitting fitting/Welcome.doc
<DIR> genint genint\Welcome.doc
<DIR> SGLITCH SGLITCH/Welcome.doc
0 File(s) 0 bytes