diff options
author | Eric Andersen <andersen@codepoet.org> | 2001-11-22 14:04:29 +0000 |
---|---|---|
committer | Eric Andersen <andersen@codepoet.org> | 2001-11-22 14:04:29 +0000 |
commit | 7ce331c01ce6eb7b3f5c715a38a24359da9c6ee2 (patch) | |
tree | 3a7e8476e868ae15f4da1b7ce26b2db6f434468c /libm/double/dawsn.c | |
parent | c117dd5fb183afb1a4790a6f6110d88704be6bf8 (diff) |
Totally rework the math library, this time based on the MacOs X
math library (which is itself based on the math lib from FreeBSD).
-Erik
Diffstat (limited to 'libm/double/dawsn.c')
-rw-r--r-- | libm/double/dawsn.c | 392 |
1 files changed, 0 insertions, 392 deletions
diff --git a/libm/double/dawsn.c b/libm/double/dawsn.c deleted file mode 100644 index 4f8d27a0c..000000000 --- a/libm/double/dawsn.c +++ /dev/null @@ -1,392 +0,0 @@ -/* dawsn.c - * - * Dawson's Integral - * - * - * - * SYNOPSIS: - * - * double x, y, dawsn(); - * - * y = dawsn( x ); - * - * - * - * DESCRIPTION: - * - * Approximates the integral - * - * x - * - - * 2 | | 2 - * dawsn(x) = exp( -x ) | exp( t ) dt - * | | - * - - * 0 - * - * Three different rational approximations are employed, for - * the intervals 0 to 3.25; 3.25 to 6.25; and 6.25 up. - * - * - * ACCURACY: - * - * Relative error: - * arithmetic domain # trials peak rms - * IEEE 0,10 10000 6.9e-16 1.0e-16 - * DEC 0,10 6000 7.4e-17 1.4e-17 - * - * - */ - -/* dawsn.c */ - - -/* -Cephes Math Library Release 2.8: June, 2000 -Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier -*/ - -#include <math.h> -/* Dawson's integral, interval 0 to 3.25 */ -#ifdef UNK -static double AN[10] = { - 1.13681498971755972054E-11, - 8.49262267667473811108E-10, - 1.94434204175553054283E-8, - 9.53151741254484363489E-7, - 3.07828309874913200438E-6, - 3.52513368520288738649E-4, --8.50149846724410912031E-4, - 4.22618223005546594270E-2, --9.17480371773452345351E-2, - 9.99999999999999994612E-1, -}; -static double AD[11] = { - 2.40372073066762605484E-11, - 1.48864681368493396752E-9, - 5.21265281010541664570E-8, - 1.27258478273186970203E-6, - 2.32490249820789513991E-5, - 3.25524741826057911661E-4, - 3.48805814657162590916E-3, - 2.79448531198828973716E-2, - 1.58874241960120565368E-1, - 5.74918629489320327824E-1, - 1.00000000000000000539E0, -}; -#endif -#ifdef DEC -static unsigned short AN[40] = { -0027107,0176630,0075752,0107612, -0030551,0070604,0166707,0127727, -0031647,0002210,0117120,0056376, -0033177,0156026,0141275,0140627, -0033516,0112200,0037035,0165515, -0035270,0150613,0016423,0105634, -0135536,0156227,0023515,0044413, -0037055,0015273,0105147,0064025, -0137273,0163145,0014460,0166465, -0040200,0000000,0000000,0000000, -}; -static unsigned short AD[44] = { -0027323,0067372,0115566,0131320, -0030714,0114432,0074206,0006637, -0032137,0160671,0044203,0026344, -0033252,0146656,0020247,0100231, -0034303,0003346,0123260,0022433, -0035252,0125460,0173041,0155415, -0036144,0113747,0125203,0124617, -0036744,0166232,0143671,0133670, -0037442,0127755,0162625,0000100, -0040023,0026736,0003604,0106265, -0040200,0000000,0000000,0000000, -}; -#endif -#ifdef IBMPC -static unsigned short AN[40] = { -0x51f1,0x0f7d,0xffb3,0x3da8, -0xf5fb,0x9db8,0x2e30,0x3e0d, -0x0ba0,0x13ca,0xe091,0x3e54, -0xb833,0xd857,0xfb82,0x3eaf, -0xbd6a,0x07c3,0xd290,0x3ec9, -0x7174,0x63a2,0x1a31,0x3f37, -0xa921,0xe4e9,0xdb92,0xbf4b, -0xed03,0x714c,0xa357,0x3fa5, -0x1da7,0xa326,0x7ccc,0xbfb7, -0x0000,0x0000,0x0000,0x3ff0, -}; -static unsigned short AD[44] = { -0xd65a,0x536e,0x6ddf,0x3dba, -0xc1b4,0x4f10,0x9323,0x3e19, -0x659c,0x2910,0xfc37,0x3e6b, -0xf013,0xc414,0x59b5,0x3eb5, -0x04a3,0xd4d6,0x60dc,0x3ef8, -0x3b62,0x1ec4,0x5566,0x3f35, -0x7532,0xf550,0x92fc,0x3f6c, -0x36f7,0x58f7,0x9d93,0x3f9c, -0xa008,0xbcb2,0x55fd,0x3fc4, -0x9197,0xc0f0,0x65bb,0x3fe2, -0x0000,0x0000,0x0000,0x3ff0, -}; -#endif -#ifdef MIEEE -static unsigned short AN[40] = { -0x3da8,0xffb3,0x0f7d,0x51f1, -0x3e0d,0x2e30,0x9db8,0xf5fb, -0x3e54,0xe091,0x13ca,0x0ba0, -0x3eaf,0xfb82,0xd857,0xb833, -0x3ec9,0xd290,0x07c3,0xbd6a, -0x3f37,0x1a31,0x63a2,0x7174, -0xbf4b,0xdb92,0xe4e9,0xa921, -0x3fa5,0xa357,0x714c,0xed03, -0xbfb7,0x7ccc,0xa326,0x1da7, -0x3ff0,0x0000,0x0000,0x0000, -}; -static unsigned short AD[44] = { -0x3dba,0x6ddf,0x536e,0xd65a, -0x3e19,0x9323,0x4f10,0xc1b4, -0x3e6b,0xfc37,0x2910,0x659c, -0x3eb5,0x59b5,0xc414,0xf013, -0x3ef8,0x60dc,0xd4d6,0x04a3, -0x3f35,0x5566,0x1ec4,0x3b62, -0x3f6c,0x92fc,0xf550,0x7532, -0x3f9c,0x9d93,0x58f7,0x36f7, -0x3fc4,0x55fd,0xbcb2,0xa008, -0x3fe2,0x65bb,0xc0f0,0x9197, -0x3ff0,0x0000,0x0000,0x0000, -}; -#endif - -/* interval 3.25 to 6.25 */ -#ifdef UNK -static double BN[11] = { - 5.08955156417900903354E-1, --2.44754418142697847934E-1, - 9.41512335303534411857E-2, --2.18711255142039025206E-2, - 3.66207612329569181322E-3, --4.23209114460388756528E-4, - 3.59641304793896631888E-5, --2.14640351719968974225E-6, - 9.10010780076391431042E-8, --2.40274520828250956942E-9, - 3.59233385440928410398E-11, -}; -static double BD[10] = { -/* 1.00000000000000000000E0,*/ --6.31839869873368190192E-1, - 2.36706788228248691528E-1, --5.31806367003223277662E-2, - 8.48041718586295374409E-3, --9.47996768486665330168E-4, - 7.81025592944552338085E-5, --4.55875153252442634831E-6, - 1.89100358111421846170E-7, --4.91324691331920606875E-9, - 7.18466403235734541950E-11, -}; -#endif -#ifdef DEC -static unsigned short BN[44] = { -0040002,0045342,0113762,0004360, -0137572,0120346,0172745,0144046, -0037300,0151134,0123440,0117047, -0136663,0025423,0014755,0046026, -0036157,0177561,0027535,0046744, -0135335,0161052,0071243,0146535, -0034426,0154060,0164506,0135625, -0133420,0005356,0100017,0151334, -0032303,0066137,0024013,0046212, -0131045,0016612,0066270,0047574, -0027435,0177025,0060625,0116363, -}; -static unsigned short BD[40] = { -/*0040200,0000000,0000000,0000000,*/ -0140041,0140101,0174552,0037073, -0037562,0061503,0124271,0160756, -0137131,0151760,0073210,0110534, -0036412,0170562,0117017,0155377, -0135570,0101374,0074056,0037276, -0034643,0145376,0001516,0060636, -0133630,0173540,0121344,0155231, -0032513,0005602,0134516,0007144, -0131250,0150540,0075747,0105341, -0027635,0177020,0012465,0125402, -}; -#endif -#ifdef IBMPC -static unsigned short BN[44] = { -0x411e,0x52fe,0x495c,0x3fe0, -0xb905,0xdebc,0x541c,0xbfcf, -0x13c5,0x94e4,0x1a4b,0x3fb8, -0xa983,0x633d,0x6562,0xbf96, -0xa9bd,0x25eb,0xffee,0x3f6d, -0x79ac,0x4e54,0xbc45,0xbf3b, -0xd773,0x1d28,0xdb06,0x3f02, -0xfa5b,0xd001,0x015d,0xbec2, -0x6991,0xe501,0x6d8b,0x3e78, -0x09f0,0x4d97,0xa3b1,0xbe24, -0xb39e,0xac32,0xbfc2,0x3dc3, -}; -static unsigned short BD[40] = { -/*0x0000,0x0000,0x0000,0x3ff0,*/ -0x47c7,0x3f2d,0x3808,0xbfe4, -0x3c3e,0x7517,0x4c68,0x3fce, -0x122b,0x0ed1,0x3a7e,0xbfab, -0xfb60,0x53c1,0x5e2e,0x3f81, -0xc7d8,0x8f05,0x105f,0xbf4f, -0xcc34,0xc069,0x795f,0x3f14, -0x9b53,0x145c,0x1eec,0xbed3, -0xc1cd,0x5729,0x6170,0x3e89, -0xf15c,0x0f7c,0x1a2c,0xbe35, -0xb560,0x02a6,0xbfc2,0x3dd3, -}; -#endif -#ifdef MIEEE -static unsigned short BN[44] = { -0x3fe0,0x495c,0x52fe,0x411e, -0xbfcf,0x541c,0xdebc,0xb905, -0x3fb8,0x1a4b,0x94e4,0x13c5, -0xbf96,0x6562,0x633d,0xa983, -0x3f6d,0xffee,0x25eb,0xa9bd, -0xbf3b,0xbc45,0x4e54,0x79ac, -0x3f02,0xdb06,0x1d28,0xd773, -0xbec2,0x015d,0xd001,0xfa5b, -0x3e78,0x6d8b,0xe501,0x6991, -0xbe24,0xa3b1,0x4d97,0x09f0, -0x3dc3,0xbfc2,0xac32,0xb39e, -}; -static unsigned short BD[40] = { -/*0x3ff0,0x0000,0x0000,0x0000,*/ -0xbfe4,0x3808,0x3f2d,0x47c7, -0x3fce,0x4c68,0x7517,0x3c3e, -0xbfab,0x3a7e,0x0ed1,0x122b, -0x3f81,0x5e2e,0x53c1,0xfb60, -0xbf4f,0x105f,0x8f05,0xc7d8, -0x3f14,0x795f,0xc069,0xcc34, -0xbed3,0x1eec,0x145c,0x9b53, -0x3e89,0x6170,0x5729,0xc1cd, -0xbe35,0x1a2c,0x0f7c,0xf15c, -0x3dd3,0xbfc2,0x02a6,0xb560, -}; -#endif - -/* 6.25 to infinity */ -#ifdef UNK -static double CN[5] = { --5.90592860534773254987E-1, - 6.29235242724368800674E-1, --1.72858975380388136411E-1, - 1.64837047825189632310E-2, --4.86827613020462700845E-4, -}; -static double CD[5] = { -/* 1.00000000000000000000E0,*/ --2.69820057197544900361E0, - 1.73270799045947845857E0, --3.93708582281939493482E-1, - 3.44278924041233391079E-2, --9.73655226040941223894E-4, -}; -#endif -#ifdef DEC -static unsigned short CN[20] = { -0140027,0030427,0176477,0074402, -0040041,0012617,0112375,0162657, -0137461,0000761,0074120,0135160, -0036607,0004325,0117246,0115525, -0135377,0036345,0064750,0047732, -}; -static unsigned short CD[20] = { -/*0040200,0000000,0000000,0000000,*/ -0140454,0127521,0071653,0133415, -0040335,0144540,0016105,0045241, -0137711,0112053,0155034,0062237, -0037015,0002102,0177442,0074546, -0135577,0036345,0064750,0052152, -}; -#endif -#ifdef IBMPC -static unsigned short CN[20] = { -0xef20,0xffa7,0xe622,0xbfe2, -0xbcb6,0xf29f,0x22b1,0x3fe4, -0x174e,0x2f0a,0x203e,0xbfc6, -0xd36b,0xb3d4,0xe11a,0x3f90, -0x09fb,0xad3d,0xe79c,0xbf3f, -}; -static unsigned short CD[20] = { -/*0x0000,0x0000,0x0000,0x3ff0,*/ -0x76e2,0x2e75,0x95ea,0xc005, -0xa954,0x0388,0xb92c,0x3ffb, -0x8c94,0x7b43,0x3285,0xbfd9, -0x4f2d,0x5fe4,0xa088,0x3fa1, -0x0a8d,0xad3d,0xe79c,0xbf4f, -}; -#endif -#ifdef MIEEE -static unsigned short CN[20] = { -0xbfe2,0xe622,0xffa7,0xef20, -0x3fe4,0x22b1,0xf29f,0xbcb6, -0xbfc6,0x203e,0x2f0a,0x174e, -0x3f90,0xe11a,0xb3d4,0xd36b, -0xbf3f,0xe79c,0xad3d,0x09fb, -}; -static unsigned short CD[20] = { -/*0x3ff0,0x0000,0x0000,0x0000,*/ -0xc005,0x95ea,0x2e75,0x76e2, -0x3ffb,0xb92c,0x0388,0xa954, -0xbfd9,0x3285,0x7b43,0x8c94, -0x3fa1,0xa088,0x5fe4,0x4f2d, -0xbf4f,0xe79c,0xad3d,0x0a8d, -}; -#endif - -#ifdef ANSIPROT -extern double chbevl ( double, void *, int ); -extern double sqrt ( double ); -extern double fabs ( double ); -extern double polevl ( double, void *, int ); -extern double p1evl ( double, void *, int ); -#else -double chbevl(), sqrt(), fabs(), polevl(), p1evl(); -#endif -extern double PI, MACHEP; - -double dawsn( xx ) -double xx; -{ -double x, y; -int sign; - - -sign = 1; -if( xx < 0.0 ) - { - sign = -1; - xx = -xx; - } - -if( xx < 3.25 ) -{ -x = xx*xx; -y = xx * polevl( x, AN, 9 )/polevl( x, AD, 10 ); -return( sign * y ); -} - - -x = 1.0/(xx*xx); - -if( xx < 6.25 ) - { - y = 1.0/xx + x * polevl( x, BN, 10) / (p1evl( x, BD, 10) * xx); - return( sign * 0.5 * y ); - } - - -if( xx > 1.0e9 ) - return( (sign * 0.5)/xx ); - -/* 6.25 to infinity */ -y = 1.0/xx + x * polevl( x, CN, 4) / (p1evl( x, CD, 5) * xx); -return( sign * 0.5 * y ); -} |