Category : Files from Magazines
Archive   : CUJ9204.ZIP
Filename : 1004024A

 
Output of file : 1004024A contained in archive : CUJ9204.ZIP
/* Listing 1 */
typedef struct {
double X1, X2, X3, X4;
} ARG_D_4; /* vector 4 */
#include "float.h"
ARG_D_4 sin_4(xx1, xx2, xx3, xx4)
double xx1, xx2, xx3, xx4;
/* use where cos of same argument not needed
** 16 digits precision, compare to 15 digits in "dtrig.h"
** T C Prince */
{
union dblfmt {
double dbl;
int idbl;
struct dfmt { /* IEEE p754 */
unsigned int sign:1;
unsigned int ex:11;
} fmt;
} xi1;
double xr, x2, x4, x8;
#ifdef __STDC__
long double pi = 3.1415926535897932385, pml;
#include
#else
register double pi = 3.1415926535897932385, pml;
#define LONG_MIN 0x80000000
#endif
union dblfmt pm, round;
ARG_D_4 res;
#define BIAS DBL_MAX_EXP
#include
#ifndef errno
extern int errno;
#endif
#define ODD(i) ((i)&1)
/* use identity sin(x + n pi) = (-1)^n sin(x)
** to reduce range to -pi/2 < x < pi/2
** pml=rint(xi1/pi) */
#if FLT_ROUNDS != 1
#error "rounding mode not nearest; adjust code"
#endif
#if FLT_RADIX !=2 && FLT_RADIX != 10
#error "code not optimum for accuracy in this RADIX"
#endif
#if DBL_DIG > 16
#error "more terms needed for full accuracy"
#endif
/* shortcut test of sign, not portable to VAX */
round.dbl = 1 / LDBL_EPSILON;
xi1.dbl = xx1;
round.idbl |= xi1.idbl & LONG_MIN;
pml = xx1 / pi + round.dbl;
/* sign reversal may reduce register usage */
xr = pi * (pml -= round.dbl) - xx1;
/* shortcut test for fabs(pml) > INT_MAX */
pm.dbl = pml;
if (pm.fmt.ex > BIAS + 31)
errno = ERANGE;
/* don't wait to calculate xr**2 until sign is fixed;
** another sign reversal is due if pm.dbl is odd */
x2 = xr * xr;
/* first sign reversal compensated in coefficient signs;
** conditional sign fixed by testing odd/even
** first two results are obtained by straight Horner
** polynomial evaluation */
res.X1 = (-.9999999999999999 + x2 * (.1666666666666607
+ x2 * (-.833333333328281e-2 + x2 * (.19841269824875e-3
+ x2 * (-.2755731661057e-5 + x2 * (.25051882036e-7
+ x2 * (-.160481709e-9 + x2 * .7374418e-12)))))))
* (ODD((int) pm.dbl) ? -xr : xr);
/* sin(xi2) */
round.dbl = 1 / LDBL_EPSILON;
xi1.dbl = xx2;
round.idbl |= xi1.idbl & LONG_MIN;
pml = xx2 / pi + round.dbl;
xr = pi * (pml -= round.dbl) - xx2;
pm.dbl = pml;
if (pm.fmt.ex > BIAS + 31)
errno = ERANGE;
x2 = xr * xr;
res.X2 = (-.9999999999999999 + x2 * (.1666666666666607
+ x2 * (-.833333333328281e-2 + x2 * (.19841269824875e-3
+ x2 * (-.2755731661057e-5 + x2 * (.25051882036e-7
+ x2 * (-.160481709e-9 + x2 * .7374418e-12)))))))
* (ODD((int) pm.dbl) ? -xr : xr);
/* sin(xi3) */
round.dbl = 1 / LDBL_EPSILON;
xi1.dbl = xx3;
round.idbl |= xi1.idbl & LONG_MIN;
pml = xx3 / pi + round.dbl;
xr = pi * (pml -= round.dbl) - xx3;
pm.dbl = pml;
if (pm.fmt.ex > BIAS + 31)
errno = ERANGE;
x2 = xr * xr;
x4 = x2 * x2;
/* split into 2 Horner polynomials to increase
** parallelism after 1st result finishes */
res.X3 = (-.9999999999999999 + x2 * (.1666666666666607
+ x2 * (-.833333333328281e-2
+ x2 * .19841269824875e-3))
+ (-.2755731661057e-5 + x2 * (.25051882036e-7
+ x2 * (-.160481709e-9
+ x2 * .7374418e-12))) * x4 * x4) *
(ODD((int) pm.dbl) ? -xr : xr);
/* sin(xi4) */
round.dbl = 1 / LDBL_EPSILON;
xi1.dbl = xx4;
round.idbl |= xi1.idbl & LONG_MIN;
pml = xx4 / pi + round.dbl;
xr = pi * (pml -= round.dbl) - xi1.dbl;
/* errno is set to ERANGE if any of the arguments are too
** large for reasonable range reduction */
pm.dbl = pml;
if (pm.fmt.ex > BIAS + 31)
errno = ERANGE;
x2 = xr * xr;
x4 = x2 * x2;
x8 = x4 * x4;
/* multiply by 1 is K&R way to enforce parentheses */
res.X4 = ((-.9999999999999999 + x2 * (.1666666666666607
+ x2 * (-.833333333328281e-2
+ x2 * .19841269824875e-3))) * 1
+ (-.2755731661057e-5 + x2 * (.25051882036e-7
+ x2 * (-.160481709e-9 + x2 * .7374418e-12))) * x8)
* (ODD((int) pm.dbl) ? -xr : xr);
return res;
}


  3 Responses to “Category : Files from Magazines
Archive   : CUJ9204.ZIP
Filename : 1004024A

  1. Very nice! Thank you for this wonderful archive. I wonder why I found it only now. Long live the BBS file archives!

  2. This is so awesome! 😀 I’d be cool if you could download an entire archive of this at once, though.

  3. But one thing that puzzles me is the “mtswslnkmcjklsdlsbdmMICROSOFT” string. There is an article about it here. It is definitely worth a read: http://www.os2museum.com/wp/mtswslnk/