workbench/libs/mathieeedoubtrans/ieeedpacos.c

   1 /*
   2     Copyright © 1995-2003, The AROS Development Team. All rights reserved.
   3     $Id$
   4 */
   5 /*
   6  * ====================================================
   7  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
   8  *
   9  * Developed at SunSoft, a Sun Microsystems, Inc. business.
  10  * Permission to use, copy, modify, and distribute this
  11  * software is freely granted, provided that this notice
  12  * is preserved.
  13  * ====================================================
  14  */
  15
  16 #include "mathieeedoubtrans_intern.h"
  17
  18 /*
  19     FUNCTION
  20       Calculate the arcus cosine of the IEEE double precision number
  21
  22     RESULT
  23       IEEE double precision floating point number
  24
  25       flags:
  26         zero     : result is zero
  27         negative : result is negative
  28         overflow : argument is out of range
  29
  30     NOTES
  31
  32     EXAMPLE
  33
  34     BUGS
  35
  36     SEE ALSO
  37
  38     INTERNALS
  39
  40     HISTORY
  41 */
  42
  43 AROS_LHQUAD1(double, IEEEDPAcos,
  44     AROS_LHAQUAD(double, x, D0, D1),
  45     struct MathIeeeDoubTransBase *, MathIeeeDoubTransBase, 20, MathIeeeDoubTrans
  46 )
  47 {
  48     AROS_LIBFUNC_INIT
  49
  50     QUAD z,p,q,r,w,s,c,df, tmp;
  51     int hx,ix;
  52
  53 #   if defined AROS_64BIT_TYPE || defined __GNUC__
  54 #   else
  55 #       undef two_64
  56 #       undef pio2_lo_64
  57 #       undef one_64
  58         const QUAD two_64     = two,
  59                    pio2_lo_64 = pio2_lo,
  60                    one_64     = oneC;
  61 #   endif
  62
  63     hx = Get_High32of64(x);
  64     ix = hx & 0x7fffffff;
  65     if (ix >= 0x3ff00000)
  66     {    /* |x| >= 1 */
  67         if (((ix-0x3ff00000) | Get_Low32of64(x)) == 0)
  68         {  /* |x|==1 */
  69             if(hx>0)
  70             {
  71                 /*
  72                     if this is a 32bit-compiler we need to define some
  73                     variables, otherwise these are available as 64bit
  74                     constants
  75                 */
  76 #               if defined AROS_64BIT_TYPE || defined __GNUC__
  77 #               else
  78 #                   undef zero_64
  79                     QUAD zero_64 = zero;
  80 #               endif
  81
  82                 return zero_64;
  83             }   /* acos(1) = 0  */
  84             else
  85             {
  86                 /*
  87                     if this is a 32bit-compiler we need to define some
  88                     variables, otherwise these are available as 64bit constants
  89                 */
  90 #               if defined AROS_64BIT_TYPE || defined __GNUC__
  91 #               else
  92 #                   undef pi_64
  93                     QUAD pi_64 = pi;
  94 #               endif
  95
  96                 return pi_64;     /* acos(-1)= pi */
  97             }
  98         }
  99         {
 100             /*
 101                 if this is a 32bit -compiler we need to define some variables,
 102                 otherwise these are available as 64bit constants
 103             */
 104 #           if defined AROS_64BIT_TYPE || defined __GNUC__
 105 #           else
 106 #               undef IEEEDPNAN_64
 107                 QUAD IEEEDPNAN_64;
 108                 Set_Value64C(IEEEDPNAN_64, IEEEDPNAN_Hi, IEEEDPNAN_Lo);
 109 #           endif
 110
 111             return IEEEDPNAN_64; /* acos(|x|>1) is NaN */
 112         }
 113     }
 114     {
 115         /*
 116             if this is a 32bit -compiler we need to define some variables,
 117             otherwise these are available as 64bit constants
 118         */
 119 #       if defined AROS_64BIT_TYPE || defined __GNUC__
 120 #       else
 121 #           undef pS0_64
 122 #           undef pS1_64
 123 #           undef pS2_64
 124 #           undef pS3_64
 125 #           undef pS4_64
 126 #           undef pS5_64
 127 #           undef qS1_64
 128 #           undef qS2_64
 129 #           undef qS3_64
 130 #           undef qS4_64
 131 #           undef pio2_hi_64
 132
 133             /* Ok, let's define some constants */
 134             const QUAD pS0_64 = pS0,
 135                        pS1_64 = pS1,
 136                        pS2_64 = pS2,
 137                        pS3_64 = pS3,
 138                        pS4_64 = pS4,
 139                        pS5_64 = pS5;
 140             const QUAD qS1_64 = qS1,
 141                        qS2_64 = qS2,
 142                        qS3_64 = qS3,
 143                        qS4_64 = qS4;
 144             const QUAD pio2_hi_64 = pio2_hi;
 145 #       endif
 146
 147         if(ix < 0x3fe00000)
 148         {         /* |x| < 0.5 */
 149             if(ix <= 0x3c600000) return pio2_hi_64;/*if|x|<2**-57*/
 150
 151             Set_Value64(z, IEEEDPMul(x, x));
 152             Set_Value64
 153             (
 154                 p,
 155                 IEEEDPMul(z, IEEEDPAdd(pS0_64,
 156                 IEEEDPMul(z, IEEEDPAdd(pS1_64,
 157                 IEEEDPMul(z, IEEEDPAdd(pS2_64,
 158                 IEEEDPMul(z, IEEEDPAdd(pS3_64,
 159                 IEEEDPMul(z, IEEEDPAdd(pS4_64,
 160                 IEEEDPMul(z, pS5_64)))))))))))
 161             );
 162
 163             Set_Value64
 164             (
 165                 q,
 166                 IEEEDPAdd(one_64,
 167                 IEEEDPMul(z, IEEEDPAdd(qS1_64,
 168                 IEEEDPMul(z, IEEEDPAdd(qS2_64,
 169                 IEEEDPMul(z, IEEEDPAdd(qS3_64,
 170                 IEEEDPMul(z, qS4_64))))))))
 171             );
 172             Set_Value64(r, IEEEDPDiv(p, q));
 173
 174             Set_Value64
 175             (
 176                 tmp,
 177                 IEEEDPSub(pio2_hi_64,
 178                 IEEEDPSub(x,
 179                 IEEEDPSub(pio2_lo_64, IEEEDPMul(x, r))))
 180             );
 181             return tmp;
 182         }
 183         else
 184         {
 185             /*
 186                 if this is a 32bit -compiler we need to define some variables,
 187                 otherwise these are available as 64bit constants
 188             */
 189 #           if defined AROS_64BIT_TYPE || defined __GNUC__
 190 #           else
 191 #               undef onehalf_64
 192                 const QUAD onehalf_64 = onehalf;
 193 #           endif
 194
 195             if (hx < 0)
 196             {       /* x < -0.5 */
 197                 /*
 198                     if this is a 32bit-compiler we need to define some
 199                     variables, otherwise these are available as 64bit constants
 200                 */
 201 #               if defined AROS_64BIT_TYPE || defined __GNUC__
 202 #               else
 203 #                   undef pi_64
 204                     const QUAD pi_64 = pi;
 205 #               endif
 206
 207                 Set_Value64(z, IEEEDPMul(IEEEDPAdd(one_64,x),onehalf_64));
 208                 Set_Value64(p, IEEEDPMul(z, IEEEDPAdd(pS0_64,
 209                                IEEEDPMul(z, IEEEDPAdd(pS1_64,
 210                                IEEEDPMul(z, IEEEDPAdd(pS2_64,
 211                                IEEEDPMul(z, IEEEDPAdd(pS3_64,
 212                                IEEEDPMul(z, IEEEDPAdd(pS4_64,
 213                                IEEEDPMul(z, pS5_64))))))))))));
 214                 Set_Value64(q, IEEEDPAdd(one_64,
 215                                IEEEDPMul(z, IEEEDPAdd(qS1_64,
 216                                IEEEDPMul(z, IEEEDPAdd(qS2_64,
 217                                IEEEDPMul(z, IEEEDPAdd(qS3_64,
 218                                IEEEDPMul(z,qS4_64)))))))));
 219                 Set_Value64(s, IEEEDPSqrt(z));
 220                 Set_Value64(r, IEEEDPDiv(p,q));
 221                 Set_Value64(w, IEEEDPSub(IEEEDPMul(r,s),pio2_lo_64));
 222                 return IEEEDPSub(pi_64, IEEEDPMul(two_64, IEEEDPAdd(s,w)));
 223             }
 224             else
 225             {  /* x > 0.5 */
 226                 Set_Value64(z, IEEEDPMul(IEEEDPSub(one_64,x),onehalf_64));
 227                 Set_Value64(s, IEEEDPSqrt(z));
 228                 Set_Value64(df,s);
 229                 AND64QC(df,0xffffffff, 0x0);
 230                 Set_Value64(c, IEEEDPDiv(IEEEDPSub(z,IEEEDPMul(df,df)),IEEEDPAdd(s,df)));
 231                 Set_Value64(p, IEEEDPMul(z, IEEEDPAdd(pS0_64,
 232                                IEEEDPMul(z, IEEEDPAdd(pS1_64,
 233                                IEEEDPMul(z, IEEEDPAdd(pS2_64,
 234                                IEEEDPMul(z, IEEEDPAdd(pS3_64,
 235                                IEEEDPMul(z, IEEEDPAdd(pS4_64,
 236                                IEEEDPMul(z,pS5_64))))))))))));
 237                 Set_Value64(q, IEEEDPAdd(one_64,
 238                                IEEEDPMul(z, IEEEDPAdd(qS1_64,
 239                                IEEEDPMul(z, IEEEDPAdd(qS2_64,
 240                                IEEEDPMul(z, IEEEDPAdd(qS3_64,
 241                                IEEEDPMul(z,qS4_64)))))))));
 242                 Set_Value64(r, IEEEDPDiv(p,q));
 243                 Set_Value64(w, IEEEDPAdd(IEEEDPMul(r,s),c));
 244                 return IEEEDPMul(two_64,IEEEDPAdd(df,w));
 245             }
 246         }
 247     }
 248
 249     AROS_LIBFUNC_EXIT
 250 }