workbench/libs/mathieeesingtrans/ieeespsin.c

   1 /*
   2     Copyright © 1995-2003, The AROS Development Team. All rights reserved.
   3     $Id$
   4 */
   5
   6 #include "mathieeesingtrans_intern.h"
   7
   8 /*
   9     FUNCTION
  10       Calculate the sine of a given IEEE single precision number in radians
  11
  12     RESULT
  13       IEEE single precision floating point number
  14
  15       flags:
  16         zero     : result is zero
  17         negative : result is negative
  18         overflow : 0
  19
  20     NOTES
  21
  22     EXAMPLE
  23
  24     BUGS
  25
  26     SEE ALSO
  27
  28     INTERNALS
  29       Algorithm for Calculation of sin(y):
  30       <code>
  31          z    = floor ( |y| / pi );
  32          y_1  = |y| - z * pi;        => 0 <= y_1 < pi
  33
  34          if (y_1 > pi/2 ) then y_1 = pi - y_1;
  35
  36          => 0 <= y_1 < pi/2
  37
  38          Res = y - y^3/3! + y^5/5! - y^7/7! + y^9/9! - y^11/11! =
  39              = y(1+y^2(-1/3!+y^2(1/5!+y^2(-1/7!+y^2(1/9!-1/11!y^2)))));
  40
  41          if (y < 0)
  42            Res = -Res;
  43
  44          if (z was an odd number)
  45            Res = -Res;
  46       </code>
  47
  48     HISTORY
  49 */
  50
  51 AROS_LH1(float, IEEESPSin,
  52     AROS_LHA(float, y, D0),
  53     struct Library *, MathIeeeSingTransBase, 6, MathIeeeSingTrans
  54 )
  55 {
  56     AROS_LIBFUNC_INIT
  57
  58     LONG z,Res,ysquared,yabs,tmp;
  59     yabs = y & (IEEESPMantisse_Mask + IEEESPExponent_Mask);
  60
  61     if (IEEESP_Pinfty == yabs)
  62     {
  63         SetSR(Overflow_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  64         return IEEESP_NAN;
  65     }
  66
  67     z = IEEESPFloor(IEEESPDiv(yabs, pi));
  68     tmp  = IEEESPMul(z,pi);
  69     tmp |= IEEESPSign_Mask; /* tmp = -tmp; */
  70     yabs = IEEESPAdd(yabs, tmp);
  71     if (yabs > pio2)
  72     {
  73         yabs |= IEEESPSign_Mask;
  74         yabs  =IEEESPAdd(pi, yabs);
  75     }
  76     ysquared = IEEESPMul(yabs,yabs);
  77     Res = IEEESPMul(yabs,
  78           IEEESPAdd(sinf1,
  79           IEEESPMul(ysquared,
  80           IEEESPAdd(sinf2,
  81           IEEESPMul(ysquared,
  82           IEEESPAdd(sinf3,
  83           IEEESPMul(ysquared,
  84           IEEESPAdd(sinf4,
  85           IEEESPMul(ysquared,
  86           IEEESPAdd(sinf5,
  87           IEEESPMul(ysquared, sinf6)))))))))));
  88
  89     if (y < 0 )                        Res ^= IEEESPSign_Mask;
  90     if (TRUE == intern_IEEESPisodd(z)) Res ^= IEEESPSign_Mask;
  91
  92     if (0 == Res)
  93     {
  94         SetSR(Zero_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  95         return 0;
  96     }
  97
  98     if (Res < 0) SetSR(Negative_Bit, Zero_Bit | Negative_Bit | Overflow_Bit);
  99
 100     return Res;
 101
 102     AROS_LIBFUNC_EXIT
 103 }