Doxygen Source Code Documentation
Main Page Alphabetical List Data Structures File List Data Fields Globals Search
cdf_53.c
Go to the documentation of this file.00001 #include "cdflib.h"
00002 double erf1(double *x)
00003
00004
00005
00006
00007
00008 {
00009 static double c = .564189583547756e0;
00010 static double a[5] = {
00011 .771058495001320e-04,-.133733772997339e-02,.323076579225834e-01,
00012 .479137145607681e-01,.128379167095513e+00
00013 };
00014 static double b[3] = {
00015 .301048631703895e-02,.538971687740286e-01,.375795757275549e+00
00016 };
00017 static double p[8] = {
00018 -1.36864857382717e-07,5.64195517478974e-01,7.21175825088309e+00,
00019 4.31622272220567e+01,1.52989285046940e+02,3.39320816734344e+02,
00020 4.51918953711873e+02,3.00459261020162e+02
00021 };
00022 static double q[8] = {
00023 1.00000000000000e+00,1.27827273196294e+01,7.70001529352295e+01,
00024 2.77585444743988e+02,6.38980264465631e+02,9.31354094850610e+02,
00025 7.90950925327898e+02,3.00459260956983e+02
00026 };
00027 static double r[5] = {
00028 2.10144126479064e+00,2.62370141675169e+01,2.13688200555087e+01,
00029 4.65807828718470e+00,2.82094791773523e-01
00030 };
00031 static double s[4] = {
00032 9.41537750555460e+01,1.87114811799590e+02,9.90191814623914e+01,
00033 1.80124575948747e+01
00034 };
00035 static double erf1,ax,bot,t,top,x2;
00036
00037
00038
00039
00040 ax = fabs(*x);
00041 if(ax > 0.5e0) goto S10;
00042 t = *x**x;
00043 top = (((a[0]*t+a[1])*t+a[2])*t+a[3])*t+a[4]+1.0e0;
00044 bot = ((b[0]*t+b[1])*t+b[2])*t+1.0e0;
00045 erf1 = *x*(top/bot);
00046 return erf1;
00047 S10:
00048 if(ax > 4.0e0) goto S20;
00049 top = ((((((p[0]*ax+p[1])*ax+p[2])*ax+p[3])*ax+p[4])*ax+p[5])*ax+p[6])*ax+p[
00050 7];
00051 bot = ((((((q[0]*ax+q[1])*ax+q[2])*ax+q[3])*ax+q[4])*ax+q[5])*ax+q[6])*ax+q[
00052 7];
00053 erf1 = 0.5e0+(0.5e0-exp(-(*x**x))*top/bot);
00054 if(*x < 0.0e0) erf1 = -erf1;
00055 return erf1;
00056 S20:
00057 if(ax >= 5.8e0) goto S30;
00058 x2 = *x**x;
00059 t = 1.0e0/x2;
00060 top = (((r[0]*t+r[1])*t+r[2])*t+r[3])*t+r[4];
00061 bot = (((s[0]*t+s[1])*t+s[2])*t+s[3])*t+1.0e0;
00062 erf1 = (c-top/(x2*bot))/ax;
00063 erf1 = 0.5e0+(0.5e0-exp(-x2)*erf1);
00064 if(*x < 0.0e0) erf1 = -erf1;
00065 return erf1;
00066 S30:
00067 erf1 = fifdsign(1.0e0,*x);
00068 return erf1;
00069 }
