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