1 /* 2 * CDDL HEADER START 3 * 4 * The contents of this file are subject to the terms of the 5 * Common Development and Distribution License (the "License"). 6 * You may not use this file except in compliance with the License. 7 * 8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE 9 * or http://www.opensolaris.org/os/licensing. 10 * See the License for the specific language governing permissions 11 * and limitations under the License. 12 * 13 * When distributing Covered Code, include this CDDL HEADER in each 14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE. 15 * If applicable, add the following below this CDDL HEADER, with the 16 * fields enclosed by brackets "[]" replaced with your own identifying 17 * information: Portions Copyright [yyyy] [name of copyright owner] 18 * 19 * CDDL HEADER END 20 */ 21 22 /* 23 * Copyright 2011 Nexenta Systems, Inc. All rights reserved. 24 */ 25 /* 26 * Copyright 2006 Sun Microsystems, Inc. All rights reserved. 27 * Use is subject to license terms. 28 */ 29 30 #include <sys/isa_defs.h> 31 #include <sys/ccompile.h> 32 33 #ifdef _LITTLE_ENDIAN 34 #define HI(x) *(1+(int*)x) 35 #define LO(x) *(unsigned*)x 36 #else 37 #define HI(x) *(int*)x 38 #define LO(x) *(1+(unsigned*)x) 39 #endif 40 41 #ifdef __RESTRICT 42 #define restrict _Restrict 43 #else 44 #define restrict 45 #endif 46 47 extern const double __vlibm_TBL_sincos_hi[], __vlibm_TBL_sincos_lo[]; 48 49 static const double 50 half[2] = { 0.5, -0.5 }, 51 one = 1.0, 52 invpio2 = 0.636619772367581343075535, 53 pio2_1 = 1.570796326734125614166, 54 pio2_2 = 6.077100506303965976596e-11, 55 pio2_3 = 2.022266248711166455796e-21, 56 pio2_3t = 8.478427660368899643959e-32, 57 pp1 = -1.666666666605760465276263943134982554676e-0001, 58 pp2 = 8.333261209690963126718376566146180944442e-0003, 59 qq1 = -4.999999999977710986407023955908711557870e-0001, 60 qq2 = 4.166654863857219350645055881018842089580e-0002, 61 poly1[2]= { -1.666666666666629669805215138920301589656e-0001, 62 -4.999999999999931701464060878888294524481e-0001 }, 63 poly2[2]= { 8.333333332390951295683993455280336376663e-0003, 64 4.166666666394861917535640593963708222319e-0002 }, 65 poly3[2]= { -1.984126237997976692791551778230098403960e-0004, 66 -1.388888552656142867832756687736851681462e-0003 }, 67 poly4[2]= { 2.753403624854277237649987622848330351110e-0006, 68 2.478519423681460796618128289454530524759e-0005 }; 69 70 static const unsigned thresh[2] = { 0x3fc90000, 0x3fc40000 }; 71 72 /* Don't __ the following; acomp will handle it */ 73 extern double fabs(double); 74 extern void __vlibm_vsin_big(int, double *, int, double *, int, int); 75 76 void 77 __vsin(int n, double * restrict x, int stridex, double * restrict y, 78 int stridey) 79 { 80 double x0_or_one[4], x1_or_one[4], x2_or_one[4]; 81 double y0_or_zero[4], y1_or_zero[4], y2_or_zero[4]; 82 double x0, x1, x2, *py0 = 0, *py1 = 0, *py2, *xsave, *ysave; 83 unsigned hx0, hx1, hx2, xsb0, xsb1 = 0, xsb2; 84 int i, biguns, nsave, sxsave, sysave; 85 volatile int v __unused; 86 nsave = n; 87 xsave = x; 88 sxsave = stridex; 89 ysave = y; 90 sysave = stridey; 91 biguns = 0; 92 93 x0 = *x; /* error: 'x0' may be used uninitialized */ 94 do 95 { 96 LOOP0: 97 xsb0 = HI(x); 98 hx0 = xsb0 & ~0x80000000; 99 if (hx0 > 0x3fe921fb) 100 { 101 biguns = 1; 102 goto MEDIUM; 103 } 104 if (hx0 < 0x3e400000) 105 { 106 v = *x; 107 *y = *x; 108 x += stridex; 109 y += stridey; 110 i = 0; 111 if (--n <= 0) 112 break; 113 goto LOOP0; 114 } 115 x0 = *x; 116 py0 = y; 117 x += stridex; 118 y += stridey; 119 i = 1; 120 if (--n <= 0) 121 break; 122 123 LOOP1: 124 xsb1 = HI(x); 125 hx1 = xsb1 & ~0x80000000; 126 if (hx1 > 0x3fe921fb) 127 { 128 biguns = 2; 129 goto MEDIUM; 130 } 131 if (hx1 < 0x3e400000) 132 { 133 v = *x; 134 *y = *x; 135 x += stridex; 136 y += stridey; 137 i = 1; 138 if (--n <= 0) 139 break; 140 goto LOOP1; 141 } 142 x1 = *x; 143 py1 = y; 144 x += stridex; 145 y += stridey; 146 i = 2; 147 if (--n <= 0) 148 break; 149 150 LOOP2: 151 xsb2 = HI(x); 152 hx2 = xsb2 & ~0x80000000; 153 if (hx2 > 0x3fe921fb) 154 { 155 biguns = 3; 156 goto MEDIUM; 157 } 158 if (hx2 < 0x3e400000) 159 { 160 v = *x; 161 *y = *x; 162 x += stridex; 163 y += stridey; 164 i = 2; 165 if (--n <= 0) 166 break; 167 goto LOOP2; 168 } 169 x2 = *x; 170 py2 = y; 171 172 i = (hx0 - 0x3fc90000) >> 31; 173 i |= ((hx1 - 0x3fc90000) >> 30) & 2; 174 i |= ((hx2 - 0x3fc90000) >> 29) & 4; 175 switch (i) 176 { 177 double a0, a1, a2, w0, w1, w2; 178 double t0, t1, t2, z0, z1, z2; 179 unsigned j0, j1, j2; 180 181 case 0: 182 j0 = (xsb0 + 0x4000) & 0xffff8000; 183 j1 = (xsb1 + 0x4000) & 0xffff8000; 184 j2 = (xsb2 + 0x4000) & 0xffff8000; 185 HI(&t0) = j0; 186 HI(&t1) = j1; 187 HI(&t2) = j2; 188 LO(&t0) = 0; 189 LO(&t1) = 0; 190 LO(&t2) = 0; 191 x0 -= t0; 192 x1 -= t1; 193 x2 -= t2; 194 z0 = x0 * x0; 195 z1 = x1 * x1; 196 z2 = x2 * x2; 197 t0 = z0 * (qq1 + z0 * qq2); 198 t1 = z1 * (qq1 + z1 * qq2); 199 t2 = z2 * (qq1 + z2 * qq2); 200 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 201 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 202 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 203 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 204 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 205 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 206 xsb0 = (xsb0 >> 30) & 2; 207 xsb1 = (xsb1 >> 30) & 2; 208 xsb2 = (xsb2 >> 30) & 2; 209 a0 = __vlibm_TBL_sincos_hi[j0+xsb0]; 210 a1 = __vlibm_TBL_sincos_hi[j1+xsb1]; 211 a2 = __vlibm_TBL_sincos_hi[j2+xsb2]; 212 t0 = (__vlibm_TBL_sincos_hi[j0+1] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+xsb0]; 213 t1 = (__vlibm_TBL_sincos_hi[j1+1] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+xsb1]; 214 t2 = (__vlibm_TBL_sincos_hi[j2+1] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+xsb2]; 215 *py0 = a0 + t0; 216 *py1 = a1 + t1; 217 *py2 = a2 + t2; 218 break; 219 220 case 1: 221 j1 = (xsb1 + 0x4000) & 0xffff8000; 222 j2 = (xsb2 + 0x4000) & 0xffff8000; 223 HI(&t1) = j1; 224 HI(&t2) = j2; 225 LO(&t1) = 0; 226 LO(&t2) = 0; 227 x1 -= t1; 228 x2 -= t2; 229 z0 = x0 * x0; 230 z1 = x1 * x1; 231 z2 = x2 * x2; 232 t0 = z0 * (poly3[0] + z0 * poly4[0]); 233 t1 = z1 * (qq1 + z1 * qq2); 234 t2 = z2 * (qq1 + z2 * qq2); 235 t0 = z0 * (poly1[0] + z0 * (poly2[0] + t0)); 236 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 237 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 238 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 239 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 240 xsb1 = (xsb1 >> 30) & 2; 241 xsb2 = (xsb2 >> 30) & 2; 242 a1 = __vlibm_TBL_sincos_hi[j1+xsb1]; 243 a2 = __vlibm_TBL_sincos_hi[j2+xsb2]; 244 t0 = x0 + x0 * t0; 245 t1 = (__vlibm_TBL_sincos_hi[j1+1] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+xsb1]; 246 t2 = (__vlibm_TBL_sincos_hi[j2+1] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+xsb2]; 247 *py0 = t0; 248 *py1 = a1 + t1; 249 *py2 = a2 + t2; 250 break; 251 252 case 2: 253 j0 = (xsb0 + 0x4000) & 0xffff8000; 254 j2 = (xsb2 + 0x4000) & 0xffff8000; 255 HI(&t0) = j0; 256 HI(&t2) = j2; 257 LO(&t0) = 0; 258 LO(&t2) = 0; 259 x0 -= t0; 260 x2 -= t2; 261 z0 = x0 * x0; 262 z1 = x1 * x1; 263 z2 = x2 * x2; 264 t0 = z0 * (qq1 + z0 * qq2); 265 t1 = z1 * (poly3[0] + z1 * poly4[0]); 266 t2 = z2 * (qq1 + z2 * qq2); 267 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 268 t1 = z1 * (poly1[0] + z1 * (poly2[0] + t1)); 269 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 270 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 271 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 272 xsb0 = (xsb0 >> 30) & 2; 273 xsb2 = (xsb2 >> 30) & 2; 274 a0 = __vlibm_TBL_sincos_hi[j0+xsb0]; 275 a2 = __vlibm_TBL_sincos_hi[j2+xsb2]; 276 t0 = (__vlibm_TBL_sincos_hi[j0+1] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+xsb0]; 277 t1 = x1 + x1 * t1; 278 t2 = (__vlibm_TBL_sincos_hi[j2+1] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+xsb2]; 279 *py0 = a0 + t0; 280 *py1 = t1; 281 *py2 = a2 + t2; 282 break; 283 284 case 3: 285 j2 = (xsb2 + 0x4000) & 0xffff8000; 286 HI(&t2) = j2; 287 LO(&t2) = 0; 288 x2 -= t2; 289 z0 = x0 * x0; 290 z1 = x1 * x1; 291 z2 = x2 * x2; 292 t0 = z0 * (poly3[0] + z0 * poly4[0]); 293 t1 = z1 * (poly3[0] + z1 * poly4[0]); 294 t2 = z2 * (qq1 + z2 * qq2); 295 t0 = z0 * (poly1[0] + z0 * (poly2[0] + t0)); 296 t1 = z1 * (poly1[0] + z1 * (poly2[0] + t1)); 297 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 298 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 299 xsb2 = (xsb2 >> 30) & 2; 300 a2 = __vlibm_TBL_sincos_hi[j2+xsb2]; 301 t0 = x0 + x0 * t0; 302 t1 = x1 + x1 * t1; 303 t2 = (__vlibm_TBL_sincos_hi[j2+1] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+xsb2]; 304 *py0 = t0; 305 *py1 = t1; 306 *py2 = a2 + t2; 307 break; 308 309 case 4: 310 j0 = (xsb0 + 0x4000) & 0xffff8000; 311 j1 = (xsb1 + 0x4000) & 0xffff8000; 312 HI(&t0) = j0; 313 HI(&t1) = j1; 314 LO(&t0) = 0; 315 LO(&t1) = 0; 316 x0 -= t0; 317 x1 -= t1; 318 z0 = x0 * x0; 319 z1 = x1 * x1; 320 z2 = x2 * x2; 321 t0 = z0 * (qq1 + z0 * qq2); 322 t1 = z1 * (qq1 + z1 * qq2); 323 t2 = z2 * (poly3[0] + z2 * poly4[0]); 324 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 325 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 326 t2 = z2 * (poly1[0] + z2 * (poly2[0] + t2)); 327 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 328 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 329 xsb0 = (xsb0 >> 30) & 2; 330 xsb1 = (xsb1 >> 30) & 2; 331 a0 = __vlibm_TBL_sincos_hi[j0+xsb0]; 332 a1 = __vlibm_TBL_sincos_hi[j1+xsb1]; 333 t0 = (__vlibm_TBL_sincos_hi[j0+1] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+xsb0]; 334 t1 = (__vlibm_TBL_sincos_hi[j1+1] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+xsb1]; 335 t2 = x2 + x2 * t2; 336 *py0 = a0 + t0; 337 *py1 = a1 + t1; 338 *py2 = t2; 339 break; 340 341 case 5: 342 j1 = (xsb1 + 0x4000) & 0xffff8000; 343 HI(&t1) = j1; 344 LO(&t1) = 0; 345 x1 -= t1; 346 z0 = x0 * x0; 347 z1 = x1 * x1; 348 z2 = x2 * x2; 349 t0 = z0 * (poly3[0] + z0 * poly4[0]); 350 t1 = z1 * (qq1 + z1 * qq2); 351 t2 = z2 * (poly3[0] + z2 * poly4[0]); 352 t0 = z0 * (poly1[0] + z0 * (poly2[0] + t0)); 353 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 354 t2 = z2 * (poly1[0] + z2 * (poly2[0] + t2)); 355 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 356 xsb1 = (xsb1 >> 30) & 2; 357 a1 = __vlibm_TBL_sincos_hi[j1+xsb1]; 358 t0 = x0 + x0 * t0; 359 t1 = (__vlibm_TBL_sincos_hi[j1+1] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+xsb1]; 360 t2 = x2 + x2 * t2; 361 *py0 = t0; 362 *py1 = a1 + t1; 363 *py2 = t2; 364 break; 365 366 case 6: 367 j0 = (xsb0 + 0x4000) & 0xffff8000; 368 HI(&t0) = j0; 369 LO(&t0) = 0; 370 x0 -= t0; 371 z0 = x0 * x0; 372 z1 = x1 * x1; 373 z2 = x2 * x2; 374 t0 = z0 * (qq1 + z0 * qq2); 375 t1 = z1 * (poly3[0] + z1 * poly4[0]); 376 t2 = z2 * (poly3[0] + z2 * poly4[0]); 377 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 378 t1 = z1 * (poly1[0] + z1 * (poly2[0] + t1)); 379 t2 = z2 * (poly1[0] + z2 * (poly2[0] + t2)); 380 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 381 xsb0 = (xsb0 >> 30) & 2; 382 a0 = __vlibm_TBL_sincos_hi[j0+xsb0]; 383 t0 = (__vlibm_TBL_sincos_hi[j0+1] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+xsb0]; 384 t1 = x1 + x1 * t1; 385 t2 = x2 + x2 * t2; 386 *py0 = a0 + t0; 387 *py1 = t1; 388 *py2 = t2; 389 break; 390 391 case 7: 392 z0 = x0 * x0; 393 z1 = x1 * x1; 394 z2 = x2 * x2; 395 t0 = z0 * (poly3[0] + z0 * poly4[0]); 396 t1 = z1 * (poly3[0] + z1 * poly4[0]); 397 t2 = z2 * (poly3[0] + z2 * poly4[0]); 398 t0 = z0 * (poly1[0] + z0 * (poly2[0] + t0)); 399 t1 = z1 * (poly1[0] + z1 * (poly2[0] + t1)); 400 t2 = z2 * (poly1[0] + z2 * (poly2[0] + t2)); 401 t0 = x0 + x0 * t0; 402 t1 = x1 + x1 * t1; 403 t2 = x2 + x2 * t2; 404 *py0 = t0; 405 *py1 = t1; 406 *py2 = t2; 407 break; 408 } 409 410 x += stridex; 411 y += stridey; 412 i = 0; 413 } while (--n > 0); 414 415 if (i > 0) 416 { 417 double a0, a1, w0, w1; 418 double t0, t1, z0, z1; 419 unsigned j0, j1; 420 421 if (i > 1) 422 { 423 if (hx1 < 0x3fc90000) 424 { 425 z1 = x1 * x1; 426 t1 = z1 * (poly3[0] + z1 * poly4[0]); 427 t1 = z1 * (poly1[0] + z1 * (poly2[0] + t1)); 428 t1 = x1 + x1 * t1; 429 *py1 = t1; 430 } 431 else 432 { 433 j1 = (xsb1 + 0x4000) & 0xffff8000; 434 HI(&t1) = j1; 435 LO(&t1) = 0; 436 x1 -= t1; 437 z1 = x1 * x1; 438 t1 = z1 * (qq1 + z1 * qq2); 439 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 440 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 441 xsb1 = (xsb1 >> 30) & 2; 442 a1 = __vlibm_TBL_sincos_hi[j1+xsb1]; 443 t1 = (__vlibm_TBL_sincos_hi[j1+1] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+xsb1]; 444 *py1 = a1 + t1; 445 } 446 } 447 if (hx0 < 0x3fc90000) 448 { 449 z0 = x0 * x0; 450 t0 = z0 * (poly3[0] + z0 * poly4[0]); 451 t0 = z0 * (poly1[0] + z0 * (poly2[0] + t0)); 452 t0 = x0 + x0 * t0; 453 *py0 = t0; 454 } 455 else 456 { 457 j0 = (xsb0 + 0x4000) & 0xffff8000; 458 HI(&t0) = j0; 459 LO(&t0) = 0; 460 x0 -= t0; 461 z0 = x0 * x0; 462 t0 = z0 * (qq1 + z0 * qq2); 463 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 464 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 465 xsb0 = (xsb0 >> 30) & 2; 466 a0 = __vlibm_TBL_sincos_hi[j0+xsb0]; 467 t0 = (__vlibm_TBL_sincos_hi[j0+1] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+xsb0]; 468 *py0 = a0 + t0; 469 } 470 } 471 472 return; 473 474 /* 475 * MEDIUM RANGE PROCESSING 476 * Jump here at first sign of medium range argument. We are a bit 477 * confused due to the jump.. fix up several variables and jump into 478 * the nth loop, same as was being processed above. 479 */ 480 481 MEDIUM: 482 483 x0_or_one[1] = 1.0; 484 x1_or_one[1] = 1.0; 485 x2_or_one[1] = 1.0; 486 x0_or_one[3] = -1.0; 487 x1_or_one[3] = -1.0; 488 x2_or_one[3] = -1.0; 489 y0_or_zero[1] = 0.0; 490 y1_or_zero[1] = 0.0; 491 y2_or_zero[1] = 0.0; 492 y0_or_zero[3] = 0.0; 493 y1_or_zero[3] = 0.0; 494 y2_or_zero[3] = 0.0; 495 496 if (biguns == 3) 497 { 498 biguns = 0; 499 xsb0 = xsb0 >> 31; 500 xsb1 = xsb1 >> 31; 501 goto loop2; 502 } 503 else if (biguns == 2) 504 { 505 xsb0 = xsb0 >> 31; 506 biguns = 0; 507 goto loop1; 508 } 509 biguns = 0; 510 511 do 512 { 513 double fn0, fn1, fn2, a0, a1, a2, w0, w1, w2, y0, y1, y2; 514 unsigned hx; 515 int n0, n1, n2; 516 517 loop0: 518 hx = HI(x); 519 xsb0 = hx >> 31; 520 hx &= ~0x80000000; 521 if (hx < 0x3e400000) 522 { 523 v = *x; 524 *y = *x; 525 x += stridex; 526 y += stridey; 527 i = 0; 528 if (--n <= 0) 529 break; 530 goto loop0; 531 } 532 if (hx > 0x413921fb) 533 { 534 if (hx >= 0x7ff00000) 535 { 536 x0 = *x; 537 *y = x0 - x0; 538 } 539 else 540 biguns = 1; 541 x += stridex; 542 y += stridey; 543 i = 0; 544 if (--n <= 0) 545 break; 546 goto loop0; 547 } 548 x0 = *x; 549 py0 = y; 550 x += stridex; 551 y += stridey; 552 i = 1; 553 if (--n <= 0) 554 break; 555 556 loop1: 557 hx = HI(x); 558 xsb1 = hx >> 31; 559 hx &= ~0x80000000; 560 if (hx < 0x3e400000) 561 { 562 v = *x; 563 *y = *x; 564 x += stridex; 565 y += stridey; 566 i = 1; 567 if (--n <= 0) 568 break; 569 goto loop1; 570 } 571 if (hx > 0x413921fb) 572 { 573 if (hx >= 0x7ff00000) 574 { 575 x1 = *x; 576 *y = x1 - x1; 577 } 578 else 579 biguns = 1; 580 x += stridex; 581 y += stridey; 582 i = 1; 583 if (--n <= 0) 584 break; 585 goto loop1; 586 } 587 x1 = *x; 588 py1 = y; 589 x += stridex; 590 y += stridey; 591 i = 2; 592 if (--n <= 0) 593 break; 594 595 loop2: 596 hx = HI(x); 597 xsb2 = hx >> 31; 598 hx &= ~0x80000000; 599 if (hx < 0x3e400000) 600 { 601 v = *x; 602 *y = *x; 603 x += stridex; 604 y += stridey; 605 i = 2; 606 if (--n <= 0) 607 break; 608 goto loop2; 609 } 610 if (hx > 0x413921fb) 611 { 612 if (hx >= 0x7ff00000) 613 { 614 x2 = *x; 615 *y = x2 - x2; 616 } 617 else 618 biguns = 1; 619 x += stridex; 620 y += stridey; 621 i = 2; 622 if (--n <= 0) 623 break; 624 goto loop2; 625 } 626 x2 = *x; 627 py2 = y; 628 629 n0 = (int) (x0 * invpio2 + half[xsb0]); 630 n1 = (int) (x1 * invpio2 + half[xsb1]); 631 n2 = (int) (x2 * invpio2 + half[xsb2]); 632 fn0 = (double) n0; 633 fn1 = (double) n1; 634 fn2 = (double) n2; 635 n0 &= 3; 636 n1 &= 3; 637 n2 &= 3; 638 a0 = x0 - fn0 * pio2_1; 639 a1 = x1 - fn1 * pio2_1; 640 a2 = x2 - fn2 * pio2_1; 641 w0 = fn0 * pio2_2; 642 w1 = fn1 * pio2_2; 643 w2 = fn2 * pio2_2; 644 x0 = a0 - w0; 645 x1 = a1 - w1; 646 x2 = a2 - w2; 647 y0 = (a0 - x0) - w0; 648 y1 = (a1 - x1) - w1; 649 y2 = (a2 - x2) - w2; 650 a0 = x0; 651 a1 = x1; 652 a2 = x2; 653 w0 = fn0 * pio2_3 - y0; 654 w1 = fn1 * pio2_3 - y1; 655 w2 = fn2 * pio2_3 - y2; 656 x0 = a0 - w0; 657 x1 = a1 - w1; 658 x2 = a2 - w2; 659 y0 = (a0 - x0) - w0; 660 y1 = (a1 - x1) - w1; 661 y2 = (a2 - x2) - w2; 662 a0 = x0; 663 a1 = x1; 664 a2 = x2; 665 w0 = fn0 * pio2_3t - y0; 666 w1 = fn1 * pio2_3t - y1; 667 w2 = fn2 * pio2_3t - y2; 668 x0 = a0 - w0; 669 x1 = a1 - w1; 670 x2 = a2 - w2; 671 y0 = (a0 - x0) - w0; 672 y1 = (a1 - x1) - w1; 673 y2 = (a2 - x2) - w2; 674 xsb0 = HI(&x0); 675 i = ((xsb0 & ~0x80000000) - thresh[n0&1]) >> 31; 676 xsb1 = HI(&x1); 677 i |= (((xsb1 & ~0x80000000) - thresh[n1&1]) >> 30) & 2; 678 xsb2 = HI(&x2); 679 i |= (((xsb2 & ~0x80000000) - thresh[n2&1]) >> 29) & 4; 680 switch (i) 681 { 682 double t0, t1, t2, z0, z1, z2; 683 unsigned j0, j1, j2; 684 685 case 0: 686 j0 = (xsb0 + 0x4000) & 0xffff8000; 687 j1 = (xsb1 + 0x4000) & 0xffff8000; 688 j2 = (xsb2 + 0x4000) & 0xffff8000; 689 HI(&t0) = j0; 690 HI(&t1) = j1; 691 HI(&t2) = j2; 692 LO(&t0) = 0; 693 LO(&t1) = 0; 694 LO(&t2) = 0; 695 x0 = (x0 - t0) + y0; 696 x1 = (x1 - t1) + y1; 697 x2 = (x2 - t2) + y2; 698 z0 = x0 * x0; 699 z1 = x1 * x1; 700 z2 = x2 * x2; 701 t0 = z0 * (qq1 + z0 * qq2); 702 t1 = z1 * (qq1 + z1 * qq2); 703 t2 = z2 * (qq1 + z2 * qq2); 704 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 705 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 706 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 707 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 708 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 709 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 710 xsb0 = (xsb0 >> 30) & 2; 711 xsb1 = (xsb1 >> 30) & 2; 712 xsb2 = (xsb2 >> 30) & 2; 713 n0 ^= (xsb0 & ~(n0 << 1)); 714 n1 ^= (xsb1 & ~(n1 << 1)); 715 n2 ^= (xsb2 & ~(n2 << 1)); 716 xsb0 |= 1; 717 xsb1 |= 1; 718 xsb2 |= 1; 719 a0 = __vlibm_TBL_sincos_hi[j0+n0]; 720 a1 = __vlibm_TBL_sincos_hi[j1+n1]; 721 a2 = __vlibm_TBL_sincos_hi[j2+n2]; 722 t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0]; 723 t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1]; 724 t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2]; 725 *py0 = ( a0 + t0 ); 726 *py1 = ( a1 + t1 ); 727 *py2 = ( a2 + t2 ); 728 break; 729 730 case 1: 731 j0 = n0 & 1; 732 j1 = (xsb1 + 0x4000) & 0xffff8000; 733 j2 = (xsb2 + 0x4000) & 0xffff8000; 734 HI(&t1) = j1; 735 HI(&t2) = j2; 736 LO(&t1) = 0; 737 LO(&t2) = 0; 738 x0_or_one[0] = x0; 739 x0_or_one[2] = -x0; 740 y0_or_zero[0] = y0; 741 y0_or_zero[2] = -y0; 742 x1 = (x1 - t1) + y1; 743 x2 = (x2 - t2) + y2; 744 z0 = x0 * x0; 745 z1 = x1 * x1; 746 z2 = x2 * x2; 747 t0 = z0 * (poly3[j0] + z0 * poly4[j0]); 748 t1 = z1 * (qq1 + z1 * qq2); 749 t2 = z2 * (qq1 + z2 * qq2); 750 t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0)); 751 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 752 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 753 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 754 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 755 xsb1 = (xsb1 >> 30) & 2; 756 xsb2 = (xsb2 >> 30) & 2; 757 n1 ^= (xsb1 & ~(n1 << 1)); 758 n2 ^= (xsb2 & ~(n2 << 1)); 759 xsb1 |= 1; 760 xsb2 |= 1; 761 a1 = __vlibm_TBL_sincos_hi[j1+n1]; 762 a2 = __vlibm_TBL_sincos_hi[j2+n2]; 763 t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0); 764 t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1]; 765 t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2]; 766 *py0 = t0; 767 *py1 = ( a1 + t1 ); 768 *py2 = ( a2 + t2 ); 769 break; 770 771 case 2: 772 j0 = (xsb0 + 0x4000) & 0xffff8000; 773 j1 = n1 & 1; 774 j2 = (xsb2 + 0x4000) & 0xffff8000; 775 HI(&t0) = j0; 776 HI(&t2) = j2; 777 LO(&t0) = 0; 778 LO(&t2) = 0; 779 x1_or_one[0] = x1; 780 x1_or_one[2] = -x1; 781 x0 = (x0 - t0) + y0; 782 y1_or_zero[0] = y1; 783 y1_or_zero[2] = -y1; 784 x2 = (x2 - t2) + y2; 785 z0 = x0 * x0; 786 z1 = x1 * x1; 787 z2 = x2 * x2; 788 t0 = z0 * (qq1 + z0 * qq2); 789 t1 = z1 * (poly3[j1] + z1 * poly4[j1]); 790 t2 = z2 * (qq1 + z2 * qq2); 791 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 792 t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1)); 793 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 794 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 795 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 796 xsb0 = (xsb0 >> 30) & 2; 797 xsb2 = (xsb2 >> 30) & 2; 798 n0 ^= (xsb0 & ~(n0 << 1)); 799 n2 ^= (xsb2 & ~(n2 << 1)); 800 xsb0 |= 1; 801 xsb2 |= 1; 802 a0 = __vlibm_TBL_sincos_hi[j0+n0]; 803 a2 = __vlibm_TBL_sincos_hi[j2+n2]; 804 t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0]; 805 t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1); 806 t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2]; 807 *py0 = ( a0 + t0 ); 808 *py1 = t1; 809 *py2 = ( a2 + t2 ); 810 break; 811 812 case 3: 813 j0 = n0 & 1; 814 j1 = n1 & 1; 815 j2 = (xsb2 + 0x4000) & 0xffff8000; 816 HI(&t2) = j2; 817 LO(&t2) = 0; 818 x0_or_one[0] = x0; 819 x0_or_one[2] = -x0; 820 x1_or_one[0] = x1; 821 x1_or_one[2] = -x1; 822 y0_or_zero[0] = y0; 823 y0_or_zero[2] = -y0; 824 y1_or_zero[0] = y1; 825 y1_or_zero[2] = -y1; 826 x2 = (x2 - t2) + y2; 827 z0 = x0 * x0; 828 z1 = x1 * x1; 829 z2 = x2 * x2; 830 t0 = z0 * (poly3[j0] + z0 * poly4[j0]); 831 t1 = z1 * (poly3[j1] + z1 * poly4[j1]); 832 t2 = z2 * (qq1 + z2 * qq2); 833 t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0)); 834 t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1)); 835 w2 = x2 * (one + z2 * (pp1 + z2 * pp2)); 836 j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 837 xsb2 = (xsb2 >> 30) & 2; 838 n2 ^= (xsb2 & ~(n2 << 1)); 839 xsb2 |= 1; 840 a2 = __vlibm_TBL_sincos_hi[j2+n2]; 841 t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0); 842 t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1); 843 t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2]; 844 *py0 = t0; 845 *py1 = t1; 846 *py2 = ( a2 + t2 ); 847 break; 848 849 case 4: 850 j0 = (xsb0 + 0x4000) & 0xffff8000; 851 j1 = (xsb1 + 0x4000) & 0xffff8000; 852 j2 = n2 & 1; 853 HI(&t0) = j0; 854 HI(&t1) = j1; 855 LO(&t0) = 0; 856 LO(&t1) = 0; 857 x2_or_one[0] = x2; 858 x2_or_one[2] = -x2; 859 x0 = (x0 - t0) + y0; 860 x1 = (x1 - t1) + y1; 861 y2_or_zero[0] = y2; 862 y2_or_zero[2] = -y2; 863 z0 = x0 * x0; 864 z1 = x1 * x1; 865 z2 = x2 * x2; 866 t0 = z0 * (qq1 + z0 * qq2); 867 t1 = z1 * (qq1 + z1 * qq2); 868 t2 = z2 * (poly3[j2] + z2 * poly4[j2]); 869 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 870 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 871 t2 = z2 * (poly1[j2] + z2 * (poly2[j2] + t2)); 872 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 873 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 874 xsb0 = (xsb0 >> 30) & 2; 875 xsb1 = (xsb1 >> 30) & 2; 876 n0 ^= (xsb0 & ~(n0 << 1)); 877 n1 ^= (xsb1 & ~(n1 << 1)); 878 xsb0 |= 1; 879 xsb1 |= 1; 880 a0 = __vlibm_TBL_sincos_hi[j0+n0]; 881 a1 = __vlibm_TBL_sincos_hi[j1+n1]; 882 t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0]; 883 t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1]; 884 t2 = x2_or_one[n2] + (y2_or_zero[n2] + x2_or_one[n2] * t2); 885 *py0 = ( a0 + t0 ); 886 *py1 = ( a1 + t1 ); 887 *py2 = t2; 888 break; 889 890 case 5: 891 j0 = n0 & 1; 892 j1 = (xsb1 + 0x4000) & 0xffff8000; 893 j2 = n2 & 1; 894 HI(&t1) = j1; 895 LO(&t1) = 0; 896 x0_or_one[0] = x0; 897 x0_or_one[2] = -x0; 898 x2_or_one[0] = x2; 899 x2_or_one[2] = -x2; 900 y0_or_zero[0] = y0; 901 y0_or_zero[2] = -y0; 902 x1 = (x1 - t1) + y1; 903 y2_or_zero[0] = y2; 904 y2_or_zero[2] = -y2; 905 z0 = x0 * x0; 906 z1 = x1 * x1; 907 z2 = x2 * x2; 908 t0 = z0 * (poly3[j0] + z0 * poly4[j0]); 909 t1 = z1 * (qq1 + z1 * qq2); 910 t2 = z2 * (poly3[j2] + z2 * poly4[j2]); 911 t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0)); 912 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 913 t2 = z2 * (poly1[j2] + z2 * (poly2[j2] + t2)); 914 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 915 xsb1 = (xsb1 >> 30) & 2; 916 n1 ^= (xsb1 & ~(n1 << 1)); 917 xsb1 |= 1; 918 a1 = __vlibm_TBL_sincos_hi[j1+n1]; 919 t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0); 920 t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1]; 921 t2 = x2_or_one[n2] + (y2_or_zero[n2] + x2_or_one[n2] * t2); 922 *py0 = t0; 923 *py1 = ( a1 + t1 ); 924 *py2 = t2; 925 break; 926 927 case 6: 928 j0 = (xsb0 + 0x4000) & 0xffff8000; 929 j1 = n1 & 1; 930 j2 = n2 & 1; 931 HI(&t0) = j0; 932 LO(&t0) = 0; 933 x1_or_one[0] = x1; 934 x1_or_one[2] = -x1; 935 x2_or_one[0] = x2; 936 x2_or_one[2] = -x2; 937 x0 = (x0 - t0) + y0; 938 y1_or_zero[0] = y1; 939 y1_or_zero[2] = -y1; 940 y2_or_zero[0] = y2; 941 y2_or_zero[2] = -y2; 942 z0 = x0 * x0; 943 z1 = x1 * x1; 944 z2 = x2 * x2; 945 t0 = z0 * (qq1 + z0 * qq2); 946 t1 = z1 * (poly3[j1] + z1 * poly4[j1]); 947 t2 = z2 * (poly3[j2] + z2 * poly4[j2]); 948 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 949 t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1)); 950 t2 = z2 * (poly1[j2] + z2 * (poly2[j2] + t2)); 951 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 952 xsb0 = (xsb0 >> 30) & 2; 953 n0 ^= (xsb0 & ~(n0 << 1)); 954 xsb0 |= 1; 955 a0 = __vlibm_TBL_sincos_hi[j0+n0]; 956 t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0]; 957 t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1); 958 t2 = x2_or_one[n2] + (y2_or_zero[n2] + x2_or_one[n2] * t2); 959 *py0 = ( a0 + t0 ); 960 *py1 = t1; 961 *py2 = t2; 962 break; 963 964 case 7: 965 j0 = n0 & 1; 966 j1 = n1 & 1; 967 j2 = n2 & 1; 968 x0_or_one[0] = x0; 969 x0_or_one[2] = -x0; 970 x1_or_one[0] = x1; 971 x1_or_one[2] = -x1; 972 x2_or_one[0] = x2; 973 x2_or_one[2] = -x2; 974 y0_or_zero[0] = y0; 975 y0_or_zero[2] = -y0; 976 y1_or_zero[0] = y1; 977 y1_or_zero[2] = -y1; 978 y2_or_zero[0] = y2; 979 y2_or_zero[2] = -y2; 980 z0 = x0 * x0; 981 z1 = x1 * x1; 982 z2 = x2 * x2; 983 t0 = z0 * (poly3[j0] + z0 * poly4[j0]); 984 t1 = z1 * (poly3[j1] + z1 * poly4[j1]); 985 t2 = z2 * (poly3[j2] + z2 * poly4[j2]); 986 t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0)); 987 t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1)); 988 t2 = z2 * (poly1[j2] + z2 * (poly2[j2] + t2)); 989 t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0); 990 t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1); 991 t2 = x2_or_one[n2] + (y2_or_zero[n2] + x2_or_one[n2] * t2); 992 *py0 = t0; 993 *py1 = t1; 994 *py2 = t2; 995 break; 996 } 997 998 x += stridex; 999 y += stridey; 1000 i = 0; 1001 } while (--n > 0); 1002 1003 if (i > 0) 1004 { 1005 double fn0, fn1, a0, a1, w0, w1, y0, y1; 1006 double t0, t1, z0, z1; 1007 unsigned j0, j1; 1008 int n0, n1; 1009 1010 if (i > 1) 1011 { 1012 n1 = (int) (x1 * invpio2 + half[xsb1]); 1013 fn1 = (double) n1; 1014 n1 &= 3; 1015 a1 = x1 - fn1 * pio2_1; 1016 w1 = fn1 * pio2_2; 1017 x1 = a1 - w1; 1018 y1 = (a1 - x1) - w1; 1019 a1 = x1; 1020 w1 = fn1 * pio2_3 - y1; 1021 x1 = a1 - w1; 1022 y1 = (a1 - x1) - w1; 1023 a1 = x1; 1024 w1 = fn1 * pio2_3t - y1; 1025 x1 = a1 - w1; 1026 y1 = (a1 - x1) - w1; 1027 xsb1 = HI(&x1); 1028 if ((xsb1 & ~0x80000000) < thresh[n1&1]) 1029 { 1030 j1 = n1 & 1; 1031 x1_or_one[0] = x1; 1032 x1_or_one[2] = -x1; 1033 y1_or_zero[0] = y1; 1034 y1_or_zero[2] = -y1; 1035 z1 = x1 * x1; 1036 t1 = z1 * (poly3[j1] + z1 * poly4[j1]); 1037 t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1)); 1038 t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1); 1039 *py1 = t1; 1040 } 1041 else 1042 { 1043 j1 = (xsb1 + 0x4000) & 0xffff8000; 1044 HI(&t1) = j1; 1045 LO(&t1) = 0; 1046 x1 = (x1 - t1) + y1; 1047 z1 = x1 * x1; 1048 t1 = z1 * (qq1 + z1 * qq2); 1049 w1 = x1 * (one + z1 * (pp1 + z1 * pp2)); 1050 j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 1051 xsb1 = (xsb1 >> 30) & 2; 1052 n1 ^= (xsb1 & ~(n1 << 1)); 1053 xsb1 |= 1; 1054 a1 = __vlibm_TBL_sincos_hi[j1+n1]; 1055 t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1]; 1056 *py1 = ( a1 + t1 ); 1057 } 1058 } 1059 n0 = (int) (x0 * invpio2 + half[xsb0]); 1060 fn0 = (double) n0; 1061 n0 &= 3; 1062 a0 = x0 - fn0 * pio2_1; 1063 w0 = fn0 * pio2_2; 1064 x0 = a0 - w0; 1065 y0 = (a0 - x0) - w0; 1066 a0 = x0; 1067 w0 = fn0 * pio2_3 - y0; 1068 x0 = a0 - w0; 1069 y0 = (a0 - x0) - w0; 1070 a0 = x0; 1071 w0 = fn0 * pio2_3t - y0; 1072 x0 = a0 - w0; 1073 y0 = (a0 - x0) - w0; 1074 xsb0 = HI(&x0); 1075 if ((xsb0 & ~0x80000000) < thresh[n0&1]) 1076 { 1077 j0 = n0 & 1; 1078 x0_or_one[0] = x0; 1079 x0_or_one[2] = -x0; 1080 y0_or_zero[0] = y0; 1081 y0_or_zero[2] = -y0; 1082 z0 = x0 * x0; 1083 t0 = z0 * (poly3[j0] + z0 * poly4[j0]); 1084 t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0)); 1085 t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0); 1086 *py0 = t0; 1087 } 1088 else 1089 { 1090 j0 = (xsb0 + 0x4000) & 0xffff8000; 1091 HI(&t0) = j0; 1092 LO(&t0) = 0; 1093 x0 = (x0 - t0) + y0; 1094 z0 = x0 * x0; 1095 t0 = z0 * (qq1 + z0 * qq2); 1096 w0 = x0 * (one + z0 * (pp1 + z0 * pp2)); 1097 j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3; 1098 xsb0 = (xsb0 >> 30) & 2; 1099 n0 ^= (xsb0 & ~(n0 << 1)); 1100 xsb0 |= 1; 1101 a0 = __vlibm_TBL_sincos_hi[j0+n0]; 1102 t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0]; 1103 *py0 = ( a0 + t0 ); 1104 } 1105 } 1106 1107 if (biguns) 1108 __vlibm_vsin_big(nsave, xsave, sxsave, ysave, sysave, 0x413921fb); 1109 } 1110