The constants s1, s2, s3, etc. are pre-computed values of 1/3!, 1/5! and so
on. The result is returned to LHS and correction in COR. */
-#define TAYLOR_SINCOS(xx, a, da, cor) \
+#define TAYLOR_SIN(xx, a, da, cor) \
({ \
double t = ((POLYNOMIAL (xx) * (a) - 0.5 * (da)) * (xx) + (da)); \
double res = (a) + t; \
if (xx < 0.01588)
{
/* Taylor series. */
- res = TAYLOR_SINCOS (xx, a, da, cor);
+ res = TAYLOR_SIN (xx, a, da, cor);
cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
retval = (res == res + cor) ? res : sloww (a, da, x);
}
if (xx < 0.01588)
{
/* Taylor series. */
- res = TAYLOR_SINCOS (xx, a, da, cor);
+ res = TAYLOR_SIN (xx, a, da, cor);
cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
retval = (res == res + cor) ? res : bsloww (a, da, x, n);
}
xx = a * a;
if (xx < 0.01588)
{
- res = TAYLOR_SINCOS (xx, a, da, cor);
+ res = TAYLOR_SIN (xx, a, da, cor);
cor = (cor > 0) ? 1.02 * cor + 1.0e-31 : 1.02 * cor - 1.0e-31;
retval = (res == res + cor) ? res : csloww (a, da, x);
}
}
if (xx < 0.01588)
{
- res = TAYLOR_SINCOS (xx, a, da, cor);
+ res = TAYLOR_SIN (xx, a, da, cor);
cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
retval = (res == res + cor) ? res : csloww (a, da, x);
}
}
if (xx < 0.01588)
{
- res = TAYLOR_SINCOS (xx, a, da, cor);
+ res = TAYLOR_SIN (xx, a, da, cor);
cor = (cor > 0) ? 1.02 * cor + eps : 1.02 * cor - eps;
retval = (res == res + cor) ? res : bsloww (a, da, x, n);
}
sloww (double x, double dx, double orig)
{
double y, t, res, cor, w[2], a, da, xn;
- union
- {
- int4 i[2];
- double x;
- } v;
+ mynumber v;
int4 n;
res = TAYLOR_SLOW (x, dx, cor);
- cor =
- (cor >
- 0) ? 1.0005 * cor + ABS (orig) * 3.1e-30 : 1.0005 * cor -
- ABS (orig) * 3.1e-30;
+ if (cor > 0)
+ cor = 1.0005 * cor + ABS (orig) * 3.1e-30;
+ else
+ cor = 1.0005 * cor - ABS (orig) * 3.1e-30;
+
if (res == res + cor)
return res;
else
csloww (double x, double dx, double orig)
{
double y, t, res, cor, w[2], a, da, xn;
- union
- {
- int4 i[2];
- double x;
- } v;
+ mynumber v;
int4 n;
/* Taylor series */