URL https://opencores.org/ocsvn/openrisc/openrisc/trunk

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-dev/] [or1k-gcc/] [libquadmath/] [math/] [sincosq_kernel.c] - Blame information for rev 740

Details | Compare with Previous | View Log


/* Quad-precision floating point sine and cosine on <-pi/4,pi/4>.
   Copyright (C) 1999 Free Software Foundation, Inc.
   This file is part of the GNU C Library.
   Contributed by Jakub Jelinek <jj@ultra.linux.cz>
 
   The GNU C Library is free software; you can redistribute it and/or
   modify it under the terms of the GNU Lesser General Public
   License as published by the Free Software Foundation; either
   version 2.1 of the License, or (at your option) any later version.
 
   The GNU C Library is distributed in the hope that it will be useful,
   but WITHOUT ANY WARRANTY; without even the implied warranty of
   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
   Lesser General Public License for more details.
 
   You should have received a copy of the GNU Lesser General Public
   License along with the GNU C Library; if not, write to the Free
   Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
   02111-1307 USA.  */
 
#include "quadmath-imp.h"
 
static const __float128 c[] = {
#define ONE c[0]
 1.00000000000000000000000000000000000E+00Q, /* 3fff0000000000000000000000000000 */
 
/* cos x ~ ONE + x^2 ( SCOS1 + SCOS2 * x^2 + ... + SCOS4 * x^6 + SCOS5 * x^8 )
   x in <0,1/256>  */
#define SCOS1 c[1]
#define SCOS2 c[2]
#define SCOS3 c[3]
#define SCOS4 c[4]
#define SCOS5 c[5]
-5.00000000000000000000000000000000000E-01Q, /* bffe0000000000000000000000000000 */
 4.16666666666666666666666666556146073E-02Q, /* 3ffa5555555555555555555555395023 */
-1.38888888888888888888309442601939728E-03Q, /* bff56c16c16c16c16c16a566e42c0375 */
 2.48015873015862382987049502531095061E-05Q, /* 3fefa01a01a019ee02dcf7da2d6d5444 */
-2.75573112601362126593516899592158083E-07Q, /* bfe927e4f5dce637cb0b54908754bde0 */
 
/* cos x ~ ONE + x^2 ( COS1 + COS2 * x^2 + ... + COS7 * x^12 + COS8 * x^14 )
   x in <0,0.1484375>  */
#define COS1 c[6]
#define COS2 c[7]
#define COS3 c[8]
#define COS4 c[9]
#define COS5 c[10]
#define COS6 c[11]
#define COS7 c[12]
#define COS8 c[13]
-4.99999999999999999999999999999999759E-01Q, /* bffdfffffffffffffffffffffffffffb */
 4.16666666666666666666666666651287795E-02Q, /* 3ffa5555555555555555555555516f30 */
-1.38888888888888888888888742314300284E-03Q, /* bff56c16c16c16c16c16c16a463dfd0d */
 2.48015873015873015867694002851118210E-05Q, /* 3fefa01a01a01a01a0195cebe6f3d3a5 */
-2.75573192239858811636614709689300351E-07Q, /* bfe927e4fb7789f5aa8142a22044b51f */
 2.08767569877762248667431926878073669E-09Q, /* 3fe21eed8eff881d1e9262d7adff4373 */
-1.14707451049343817400420280514614892E-11Q, /* bfda9397496922a9601ed3d4ca48944b */
 4.77810092804389587579843296923533297E-14Q, /* 3fd2ae5f8197cbcdcaf7c3fb4523414c */
 
/* sin x ~ ONE * x + x^3 ( SSIN1 + SSIN2 * x^2 + ... + SSIN4 * x^6 + SSIN5 * x^8 )
   x in <0,1/256>  */
#define SSIN1 c[14]
#define SSIN2 c[15]
#define SSIN3 c[16]
#define SSIN4 c[17]
#define SSIN5 c[18]
-1.66666666666666666666666666666666659E-01Q, /* bffc5555555555555555555555555555 */
 8.33333333333333333333333333146298442E-03Q, /* 3ff81111111111111111111110fe195d */
-1.98412698412698412697726277416810661E-04Q, /* bff2a01a01a01a01a019e7121e080d88 */
 2.75573192239848624174178393552189149E-06Q, /* 3fec71de3a556c640c6aaa51aa02ab41 */
-2.50521016467996193495359189395805639E-08Q, /* bfe5ae644ee90c47dc71839de75b2787 */
 
/* sin x ~ ONE * x + x^3 ( SIN1 + SIN2 * x^2 + ... + SIN7 * x^12 + SIN8 * x^14 )
   x in <0,0.1484375>  */
#define SIN1 c[19]
#define SIN2 c[20]
#define SIN3 c[21]
#define SIN4 c[22]
#define SIN5 c[23]
#define SIN6 c[24]
#define SIN7 c[25]
#define SIN8 c[26]
-1.66666666666666666666666666666666538e-01Q, /* bffc5555555555555555555555555550 */
 8.33333333333333333333333333307532934e-03Q, /* 3ff811111111111111111111110e7340 */
-1.98412698412698412698412534478712057e-04Q, /* bff2a01a01a01a01a01a019e7a626296 */
 2.75573192239858906520896496653095890e-06Q, /* 3fec71de3a556c7338fa38527474b8f5 */
-2.50521083854417116999224301266655662e-08Q, /* bfe5ae64567f544e16c7de65c2ea551f */
 1.60590438367608957516841576404938118e-10Q, /* 3fde6124613a811480538a9a41957115 */
-7.64716343504264506714019494041582610e-13Q, /* bfd6ae7f3d5aef30c7bc660b060ef365 */
 2.81068754939739570236322404393398135e-15Q, /* 3fce9510115aabf87aceb2022a9a9180 */
};
 
#define SINCOSQ_COS_HI 0
#define SINCOSQ_COS_LO 1
#define SINCOSQ_SIN_HI 2
#define SINCOSQ_SIN_LO 3
extern const __float128 __sincosq_table[];
 
void
__quadmath_kernel_sincosq(__float128 x, __float128 y, __float128 *sinx,
                          __float128 *cosx, int iy)
{
  __float128 h, l, z, sin_l, cos_l_m1;
  int64_t ix;
  uint32_t tix, hix, index;
  GET_FLT128_MSW64 (ix, x);
  tix = ((uint64_t)ix) >> 32;
  tix &= ~0x80000000;                   /* tix = |x|'s high 32 bits */
  if (tix < 0x3ffc3000)                 /* |x| < 0.1484375 */
    {
      /* Argument is small enough to approximate it by a Chebyshev
         polynomial of degree 16(17).  */
      if (tix < 0x3fc60000)             /* |x| < 2^-57 */
        if (!((int)x))                  /* generate inexact */
          {
            *sinx = x;
            *cosx = ONE;
            return;
          }
      z = x * x;
      *sinx = x + (x * (z*(SIN1+z*(SIN2+z*(SIN3+z*(SIN4+
                        z*(SIN5+z*(SIN6+z*(SIN7+z*SIN8)))))))));
      *cosx = ONE + (z*(COS1+z*(COS2+z*(COS3+z*(COS4+
                     z*(COS5+z*(COS6+z*(COS7+z*COS8))))))));
    }
  else
    {
      /* So that we don't have to use too large polynomial,  we find
         l and h such that x = l + h,  where fabsl(l) <= 1.0/256 with 83
         possible values for h.  We look up cosl(h) and sinl(h) in
         pre-computed tables,  compute cosl(l) and sinl(l) using a
         Chebyshev polynomial of degree 10(11) and compute
         sinl(h+l) = sinl(h)cosl(l) + cosl(h)sinl(l) and
         cosl(h+l) = cosl(h)cosl(l) - sinl(h)sinl(l).  */
      index = 0x3ffe - (tix >> 16);
      hix = (tix + (0x200 << index)) & (0xfffffc00 << index);
      x = fabsq (x);
      switch (index)
        {
        case 0: index = ((45 << 10) + hix - 0x3ffe0000) >> 8; break;
        case 1: index = ((13 << 11) + hix - 0x3ffd0000) >> 9; break;
        default:
        case 2: index = (hix - 0x3ffc3000) >> 10; break;
        }
 
      SET_FLT128_WORDS64(h, ((uint64_t)hix) << 32, 0);
      if (iy)
        l = y - (h - x);
      else
        l = x - h;
      z = l * l;
      sin_l = l*(ONE+z*(SSIN1+z*(SSIN2+z*(SSIN3+z*(SSIN4+z*SSIN5)))));
      cos_l_m1 = z*(SCOS1+z*(SCOS2+z*(SCOS3+z*(SCOS4+z*SCOS5))));
      z = __sincosq_table [index + SINCOSQ_SIN_HI]
          + (__sincosq_table [index + SINCOSQ_SIN_LO]
             + (__sincosq_table [index + SINCOSQ_SIN_HI] * cos_l_m1)
             + (__sincosq_table [index + SINCOSQ_COS_HI] * sin_l));
      *sinx = (ix < 0) ? -z : z;
      *cosx = __sincosq_table [index + SINCOSQ_COS_HI]
              + (__sincosq_table [index + SINCOSQ_COS_LO]
                 - (__sincosq_table [index + SINCOSQ_SIN_HI] * sin_l
                    - __sincosq_table [index + SINCOSQ_COS_HI] * cos_l_m1));
    }
}

Line No.	Rev	Author	Line
1	740	jeremybenn	`/* Quad-precision floating point sine and cosine on <-pi/4,pi/4>.`
2			`Copyright (C) 1999 Free Software Foundation, Inc.`
3			`This file is part of the GNU C Library.`
4			`Contributed by Jakub Jelinek <jj@ultra.linux.cz>`
5
6			`The GNU C Library is free software; you can redistribute it and/or`
7			`modify it under the terms of the GNU Lesser General Public`
8			`License as published by the Free Software Foundation; either`
9			`version 2.1 of the License, or (at your option) any later version.`
10
11			`The GNU C Library is distributed in the hope that it will be useful,`
12			`but WITHOUT ANY WARRANTY; without even the implied warranty of`
13			`MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU`
14			`Lesser General Public License for more details.`
15
16			`You should have received a copy of the GNU Lesser General Public`
17			`License along with the GNU C Library; if not, write to the Free`
18			`Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA`
19			`02111-1307 USA. */`
20
21			`#include "quadmath-imp.h"`
22
23			`static const __float128 c[] = {`
24			`#define ONE c[0]`
25			`1.00000000000000000000000000000000000E+00Q, /* 3fff0000000000000000000000000000 */`
26
27			`/* cos x ~ ONE + x^2 ( SCOS1 + SCOS2 * x^2 + ... + SCOS4 * x^6 + SCOS5 * x^8 )`
28			`x in <0,1/256> */`
29			`#define SCOS1 c[1]`
30			`#define SCOS2 c[2]`
31			`#define SCOS3 c[3]`
32			`#define SCOS4 c[4]`
33			`#define SCOS5 c[5]`
34			`-5.00000000000000000000000000000000000E-01Q, /* bffe0000000000000000000000000000 */`
35			`4.16666666666666666666666666556146073E-02Q, /* 3ffa5555555555555555555555395023 */`
36			`-1.38888888888888888888309442601939728E-03Q, /* bff56c16c16c16c16c16a566e42c0375 */`
37			`2.48015873015862382987049502531095061E-05Q, /* 3fefa01a01a019ee02dcf7da2d6d5444 */`
38			`-2.75573112601362126593516899592158083E-07Q, /* bfe927e4f5dce637cb0b54908754bde0 */`
39
40			`/* cos x ~ ONE + x^2 ( COS1 + COS2 * x^2 + ... + COS7 * x^12 + COS8 * x^14 )`
41			`x in <0,0.1484375> */`
42			`#define COS1 c[6]`
43			`#define COS2 c[7]`
44			`#define COS3 c[8]`
45			`#define COS4 c[9]`
46			`#define COS5 c[10]`
47			`#define COS6 c[11]`
48			`#define COS7 c[12]`
49			`#define COS8 c[13]`
50			`-4.99999999999999999999999999999999759E-01Q, /* bffdfffffffffffffffffffffffffffb */`
51			`4.16666666666666666666666666651287795E-02Q, /* 3ffa5555555555555555555555516f30 */`
52			`-1.38888888888888888888888742314300284E-03Q, /* bff56c16c16c16c16c16c16a463dfd0d */`
53			`2.48015873015873015867694002851118210E-05Q, /* 3fefa01a01a01a01a0195cebe6f3d3a5 */`
54			`-2.75573192239858811636614709689300351E-07Q, /* bfe927e4fb7789f5aa8142a22044b51f */`
55			`2.08767569877762248667431926878073669E-09Q, /* 3fe21eed8eff881d1e9262d7adff4373 */`
56			`-1.14707451049343817400420280514614892E-11Q, /* bfda9397496922a9601ed3d4ca48944b */`
57			`4.77810092804389587579843296923533297E-14Q, /* 3fd2ae5f8197cbcdcaf7c3fb4523414c */`
58
59			`/* sin x ~ ONE * x + x^3 ( SSIN1 + SSIN2 * x^2 + ... + SSIN4 * x^6 + SSIN5 * x^8 )`
60			`x in <0,1/256> */`
61			`#define SSIN1 c[14]`
62			`#define SSIN2 c[15]`
63			`#define SSIN3 c[16]`
64			`#define SSIN4 c[17]`
65			`#define SSIN5 c[18]`
66			`-1.66666666666666666666666666666666659E-01Q, /* bffc5555555555555555555555555555 */`
67			`8.33333333333333333333333333146298442E-03Q, /* 3ff81111111111111111111110fe195d */`
68			`-1.98412698412698412697726277416810661E-04Q, /* bff2a01a01a01a01a019e7121e080d88 */`
69			`2.75573192239848624174178393552189149E-06Q, /* 3fec71de3a556c640c6aaa51aa02ab41 */`
70			`-2.50521016467996193495359189395805639E-08Q, /* bfe5ae644ee90c47dc71839de75b2787 */`
71
72			`/* sin x ~ ONE * x + x^3 ( SIN1 + SIN2 * x^2 + ... + SIN7 * x^12 + SIN8 * x^14 )`
73			`x in <0,0.1484375> */`
74			`#define SIN1 c[19]`
75			`#define SIN2 c[20]`
76			`#define SIN3 c[21]`
77			`#define SIN4 c[22]`
78			`#define SIN5 c[23]`
79			`#define SIN6 c[24]`
80			`#define SIN7 c[25]`
81			`#define SIN8 c[26]`
82			`-1.66666666666666666666666666666666538e-01Q, /* bffc5555555555555555555555555550 */`
83			`8.33333333333333333333333333307532934e-03Q, /* 3ff811111111111111111111110e7340 */`
84			`-1.98412698412698412698412534478712057e-04Q, /* bff2a01a01a01a01a01a019e7a626296 */`
85			`2.75573192239858906520896496653095890e-06Q, /* 3fec71de3a556c7338fa38527474b8f5 */`
86			`-2.50521083854417116999224301266655662e-08Q, /* bfe5ae64567f544e16c7de65c2ea551f */`
87			`1.60590438367608957516841576404938118e-10Q, /* 3fde6124613a811480538a9a41957115 */`
88			`-7.64716343504264506714019494041582610e-13Q, /* bfd6ae7f3d5aef30c7bc660b060ef365 */`
89			`2.81068754939739570236322404393398135e-15Q, /* 3fce9510115aabf87aceb2022a9a9180 */`
90			`};`
91
92			`#define SINCOSQ_COS_HI 0`
93			`#define SINCOSQ_COS_LO 1`
94			`#define SINCOSQ_SIN_HI 2`
95			`#define SINCOSQ_SIN_LO 3`
96			`extern const __float128 __sincosq_table[];`
97
98			`void`
99			`__quadmath_kernel_sincosq(__float128 x, __float128 y, __float128 *sinx,`
100			`__float128 *cosx, int iy)`
101			`{`
102			`__float128 h, l, z, sin_l, cos_l_m1;`
103			`int64_t ix;`
104			`uint32_t tix, hix, index;`
105			`GET_FLT128_MSW64 (ix, x);`
106			`tix = ((uint64_t)ix) >> 32;`
107			`tix &= ~0x80000000; /* tix = \|x\|'s high 32 bits */`
108			`if (tix < 0x3ffc3000) /* \|x\| < 0.1484375 */`
109			`{`
110			`/* Argument is small enough to approximate it by a Chebyshev`
111			`polynomial of degree 16(17). */`
112			`if (tix < 0x3fc60000) /* \|x\| < 2^-57 */`
113			`if (!((int)x)) /* generate inexact */`
114			`{`
115			`*sinx = x;`
116			`*cosx = ONE;`
117			`return;`
118			`}`
119			`z = x * x;`
120			`sinx = x + (x (z(SIN1+z(SIN2+z(SIN3+z(SIN4+`
121			`z(SIN5+z(SIN6+z(SIN7+zSIN8)))))))));`
122			`cosx = ONE + (z(COS1+z(COS2+z(COS3+z*(COS4+`
123			`z(COS5+z(COS6+z(COS7+zCOS8))))))));`
124			`}`
125			`else`
126			`{`
127			`/* So that we don't have to use too large polynomial, we find`
128			`l and h such that x = l + h, where fabsl(l) <= 1.0/256 with 83`
129			`possible values for h. We look up cosl(h) and sinl(h) in`
130			`pre-computed tables, compute cosl(l) and sinl(l) using a`
131			`Chebyshev polynomial of degree 10(11) and compute`
132			`sinl(h+l) = sinl(h)cosl(l) + cosl(h)sinl(l) and`
133			`cosl(h+l) = cosl(h)cosl(l) - sinl(h)sinl(l). */`
134			`index = 0x3ffe - (tix >> 16);`
135			`hix = (tix + (0x200 << index)) & (0xfffffc00 << index);`
136			`x = fabsq (x);`
137			`switch (index)`
138			`{`
139			`case 0: index = ((45 << 10) + hix - 0x3ffe0000) >> 8; break;`
140			`case 1: index = ((13 << 11) + hix - 0x3ffd0000) >> 9; break;`
141			`default:`
142			`case 2: index = (hix - 0x3ffc3000) >> 10; break;`
143			`}`
144
145			`SET_FLT128_WORDS64(h, ((uint64_t)hix) << 32, 0);`
146			`if (iy)`
147			`l = y - (h - x);`
148			`else`
149			`l = x - h;`
150			`z = l * l;`
151			`sin_l = l(ONE+z(SSIN1+z(SSIN2+z(SSIN3+z(SSIN4+zSSIN5)))));`
152			`cos_l_m1 = z(SCOS1+z(SCOS2+z(SCOS3+z(SCOS4+z*SCOS5))));`
153			`z = __sincosq_table [index + SINCOSQ_SIN_HI]`
154			`+ (__sincosq_table [index + SINCOSQ_SIN_LO]`
155			`+ (__sincosq_table [index + SINCOSQ_SIN_HI] * cos_l_m1)`
156			`+ (__sincosq_table [index + SINCOSQ_COS_HI] * sin_l));`
157			`*sinx = (ix < 0) ? -z : z;`
158			`*cosx = __sincosq_table [index + SINCOSQ_COS_HI]`
159			`+ (__sincosq_table [index + SINCOSQ_COS_LO]`
160			`- (__sincosq_table [index + SINCOSQ_SIN_HI] * sin_l`
161			`- __sincosq_table [index + SINCOSQ_COS_HI] * cos_l_m1));`
162			`}`
163			`}`

Browse

Tools

Subversion Repositories openrisc

[/] [openrisc/] [trunk/] [gnu-dev/] [or1k-gcc/] [libquadmath/] [math/] [sincosq_kernel.c] - Blame information for rev 740