125c28e8Piotr Jasiukajtis/*
225c28e8Piotr Jasiukajtis * CDDL HEADER START
325c28e8Piotr Jasiukajtis *
425c28e8Piotr Jasiukajtis * The contents of this file are subject to the terms of the
525c28e8Piotr Jasiukajtis * Common Development and Distribution License (the "License").
625c28e8Piotr Jasiukajtis * You may not use this file except in compliance with the License.
725c28e8Piotr Jasiukajtis *
825c28e8Piotr Jasiukajtis * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
925c28e8Piotr Jasiukajtis * or http://www.opensolaris.org/os/licensing.
1025c28e8Piotr Jasiukajtis * See the License for the specific language governing permissions
1125c28e8Piotr Jasiukajtis * and limitations under the License.
1225c28e8Piotr Jasiukajtis *
1325c28e8Piotr Jasiukajtis * When distributing Covered Code, include this CDDL HEADER in each
1425c28e8Piotr Jasiukajtis * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
1525c28e8Piotr Jasiukajtis * If applicable, add the following below this CDDL HEADER, with the
1625c28e8Piotr Jasiukajtis * fields enclosed by brackets "[]" replaced with your own identifying
1725c28e8Piotr Jasiukajtis * information: Portions Copyright [yyyy] [name of copyright owner]
1825c28e8Piotr Jasiukajtis *
1925c28e8Piotr Jasiukajtis * CDDL HEADER END
2025c28e8Piotr Jasiukajtis */
2125c28e8Piotr Jasiukajtis
2225c28e8Piotr Jasiukajtis/*
2325c28e8Piotr Jasiukajtis * Copyright 2011 Nexenta Systems, Inc.  All rights reserved.
2425c28e8Piotr Jasiukajtis */
2525c28e8Piotr Jasiukajtis/*
2625c28e8Piotr Jasiukajtis * Copyright 2006 Sun Microsystems, Inc.  All rights reserved.
2725c28e8Piotr Jasiukajtis * Use is subject to license terms.
2825c28e8Piotr Jasiukajtis */
2925c28e8Piotr Jasiukajtis
3025c28e8Piotr Jasiukajtis#include <sys/isa_defs.h>
3125c28e8Piotr Jasiukajtis#include <sys/ccompile.h>
3225c28e8Piotr Jasiukajtis
3325c28e8Piotr Jasiukajtis#ifdef _LITTLE_ENDIAN
3425c28e8Piotr Jasiukajtis#define HI(x)	*(1+(int*)x)
3525c28e8Piotr Jasiukajtis#define LO(x)	*(unsigned*)x
3625c28e8Piotr Jasiukajtis#else
3725c28e8Piotr Jasiukajtis#define HI(x)	*(int*)x
3825c28e8Piotr Jasiukajtis#define LO(x)	*(1+(unsigned*)x)
3925c28e8Piotr Jasiukajtis#endif
4025c28e8Piotr Jasiukajtis
4125c28e8Piotr Jasiukajtis#ifdef __RESTRICT
4225c28e8Piotr Jasiukajtis#define restrict _Restrict
4325c28e8Piotr Jasiukajtis#else
4425c28e8Piotr Jasiukajtis#define restrict
4525c28e8Piotr Jasiukajtis#endif
4625c28e8Piotr Jasiukajtis
4725c28e8Piotr Jasiukajtis/*
4825c28e8Piotr Jasiukajtis * vcos.1.c
4925c28e8Piotr Jasiukajtis *
5025c28e8Piotr Jasiukajtis * Vector cosine function.  Just slight modifications to vsin.8.c, mainly
5125c28e8Piotr Jasiukajtis * in the primary range part.
5225c28e8Piotr Jasiukajtis *
5325c28e8Piotr Jasiukajtis * Modification to primary range processing.  If an argument that does not
5425c28e8Piotr Jasiukajtis * fall in the primary range is encountered, then processing is continued
5525c28e8Piotr Jasiukajtis * in the medium range.
5625c28e8Piotr Jasiukajtis *
5725c28e8Piotr Jasiukajtis */
5825c28e8Piotr Jasiukajtis
5925c28e8Piotr Jasiukajtisextern const double __vlibm_TBL_sincos_hi[], __vlibm_TBL_sincos_lo[];
6025c28e8Piotr Jasiukajtis
6125c28e8Piotr Jasiukajtisstatic const double
6225c28e8Piotr Jasiukajtis	half[2]	= { 0.5, -0.5 },
6325c28e8Piotr Jasiukajtis	one		= 1.0,
6425c28e8Piotr Jasiukajtis	invpio2 = 0.636619772367581343075535,  /* 53 bits of pi/2 */
6525c28e8Piotr Jasiukajtis	pio2_1	= 1.570796326734125614166,	/* first 33 bits of pi/2 */
6625c28e8Piotr Jasiukajtis	pio2_2	= 6.077100506303965976596e-11, /* second 33 bits of pi/2 */
6725c28e8Piotr Jasiukajtis	pio2_3	= 2.022266248711166455796e-21, /* third 33 bits of pi/2 */
6825c28e8Piotr Jasiukajtis	pio2_3t	= 8.478427660368899643959e-32, /* pi/2 - pio2_3 */
6925c28e8Piotr Jasiukajtis	pp1		= -1.666666666605760465276263943134982554676e-0001,
7025c28e8Piotr Jasiukajtis	pp2		=  8.333261209690963126718376566146180944442e-0003,
7125c28e8Piotr Jasiukajtis	qq1		= -4.999999999977710986407023955908711557870e-0001,
7225c28e8Piotr Jasiukajtis	qq2		=  4.166654863857219350645055881018842089580e-0002,
7325c28e8Piotr Jasiukajtis	poly1[2]= { -1.666666666666629669805215138920301589656e-0001,
7425c28e8Piotr Jasiukajtis				-4.999999999999931701464060878888294524481e-0001 },
7525c28e8Piotr Jasiukajtis	poly2[2]= {  8.333333332390951295683993455280336376663e-0003,
7625c28e8Piotr Jasiukajtis				 4.166666666394861917535640593963708222319e-0002 },
7725c28e8Piotr Jasiukajtis	poly3[2]= { -1.984126237997976692791551778230098403960e-0004,
7825c28e8Piotr Jasiukajtis				-1.388888552656142867832756687736851681462e-0003 },
7925c28e8Piotr Jasiukajtis	poly4[2]= {  2.753403624854277237649987622848330351110e-0006,
8025c28e8Piotr Jasiukajtis				 2.478519423681460796618128289454530524759e-0005 };
8125c28e8Piotr Jasiukajtis
8225c28e8Piotr Jasiukajtisstatic const unsigned thresh[2] = { 0x3fc90000, 0x3fc40000 };
8325c28e8Piotr Jasiukajtis
8425c28e8Piotr Jasiukajtis/* Don't __ the following; acomp will handle it */
8525c28e8Piotr Jasiukajtisextern double fabs(double);
8625c28e8Piotr Jasiukajtisextern void __vlibm_vcos_big(int, double *, int, double *, int, int);
8725c28e8Piotr Jasiukajtis
8825c28e8Piotr Jasiukajtis/*
8925c28e8Piotr Jasiukajtis * y[i*stridey] := cos( x[i*stridex] ), for i = 0..n.
9025c28e8Piotr Jasiukajtis *
9125c28e8Piotr Jasiukajtis * Calls __vlibm_vcos_big to handle all elts which have abs >~ 1.647e+06.
9225c28e8Piotr Jasiukajtis * Argument reduction is done here for elts pi/4 < arg < 1.647e+06.
9325c28e8Piotr Jasiukajtis *
9425c28e8Piotr Jasiukajtis * elts < 2^-27 use the approximation 1.0 ~ cos(x).
9525c28e8Piotr Jasiukajtis */
9625c28e8Piotr Jasiukajtisvoid
9725c28e8Piotr Jasiukajtis__vcos(int n, double * restrict x, int stridex, double * restrict y,
9825c28e8Piotr Jasiukajtis	int stridey)
9925c28e8Piotr Jasiukajtis{
10025c28e8Piotr Jasiukajtis	double		x0_or_one[4], x1_or_one[4], x2_or_one[4];
10125c28e8Piotr Jasiukajtis	double		y0_or_zero[4], y1_or_zero[4], y2_or_zero[4];
10225c28e8Piotr Jasiukajtis	double		x0, x1, x2, *py0 = 0, *py1 = 0, *py2, *xsave, *ysave;
10325c28e8Piotr Jasiukajtis	unsigned	hx0, hx1, hx2, xsb0, xsb1 = 0, xsb2;
10425c28e8Piotr Jasiukajtis	int		i, biguns, nsave, sxsave, sysave;
105c81a25eToomas Soome	volatile int	v __unused;
10625c28e8Piotr Jasiukajtis	nsave = n;
10725c28e8Piotr Jasiukajtis	xsave = x;
10825c28e8Piotr Jasiukajtis	sxsave = stridex;
10925c28e8Piotr Jasiukajtis	ysave = y;
11025c28e8Piotr Jasiukajtis	sysave = stridey;
11125c28e8Piotr Jasiukajtis	biguns = 0;
11225c28e8Piotr Jasiukajtis
1138b8bfb2Toomas Soome	x0 = *x;	/* 'x0' may be used uninitialized */
11425c28e8Piotr Jasiukajtis	do /* MAIN LOOP */
11525c28e8Piotr Jasiukajtis	{
11625c28e8Piotr Jasiukajtis		/* Gotos here so _break_ exits MAIN LOOP. */
11725c28e8Piotr JasiukajtisLOOP0:  /* Find first arg in right range. */
11825c28e8Piotr Jasiukajtis		xsb0 = HI(x); /* get most significant word */
11925c28e8Piotr Jasiukajtis		hx0 = xsb0 & ~0x80000000; /* mask off sign bit */
12025c28e8Piotr Jasiukajtis		if (hx0 > 0x3fe921fb) {
12125c28e8Piotr Jasiukajtis			/* Too big: arg reduction needed, so leave for second part */
12225c28e8Piotr Jasiukajtis			biguns = 1;
12325c28e8Piotr Jasiukajtis			goto MEDIUM;
12425c28e8Piotr Jasiukajtis		}
12525c28e8Piotr Jasiukajtis		if (hx0 < 0x3e400000) {
12625c28e8Piotr Jasiukajtis			/* Too small.  cos x ~ 1. */
12725c28e8Piotr Jasiukajtis			v = *x;
12825c28e8Piotr Jasiukajtis			*y = 1.0;
12925c28e8Piotr Jasiukajtis			x += stridex;
13025c28e8Piotr Jasiukajtis			y += stridey;
13125c28e8Piotr Jasiukajtis			i = 0;
13225c28e8Piotr Jasiukajtis			if (--n <= 0)
13325c28e8Piotr Jasiukajtis				break;
13425c28e8Piotr Jasiukajtis			goto LOOP0;
13525c28e8Piotr Jasiukajtis		}
13625c28e8Piotr Jasiukajtis		x0 = *x;
13725c28e8Piotr Jasiukajtis		py0 = y;
13825c28e8Piotr Jasiukajtis		x += stridex;
13925c28e8Piotr Jasiukajtis		y += stridey;
14025c28e8Piotr Jasiukajtis		i = 1;
14125c28e8Piotr Jasiukajtis		if (--n <= 0)
14225c28e8Piotr Jasiukajtis			break;
14325c28e8Piotr Jasiukajtis
14425c28e8Piotr JasiukajtisLOOP1: /* Get second arg, same as above. */
14525c28e8Piotr Jasiukajtis		xsb1 = HI(x);
14625c28e8Piotr Jasiukajtis		hx1 = xsb1 & ~0x80000000;
14725c28e8Piotr Jasiukajtis		if (hx1 > 0x3fe921fb)
14825c28e8Piotr Jasiukajtis		{
14925c28e8Piotr Jasiukajtis			biguns = 2;
15025c28e8Piotr Jasiukajtis			goto MEDIUM;
15125c28e8Piotr Jasiukajtis		}
15225c28e8Piotr Jasiukajtis		if (hx1 < 0x3e400000)
15325c28e8Piotr Jasiukajtis		{
15425c28e8Piotr Jasiukajtis			v = *x;
15525c28e8Piotr Jasiukajtis			*y = 1.0;
15625c28e8Piotr Jasiukajtis			x += stridex;
15725c28e8Piotr Jasiukajtis			y += stridey;
15825c28e8Piotr Jasiukajtis			i = 1;
15925c28e8Piotr Jasiukajtis			if (--n <= 0)
16025c28e8Piotr Jasiukajtis				break;
16125c28e8Piotr Jasiukajtis			goto LOOP1;
16225c28e8Piotr Jasiukajtis		}
16325c28e8Piotr Jasiukajtis		x1 = *x;
16425c28e8Piotr Jasiukajtis		py1 = y;
16525c28e8Piotr Jasiukajtis		x += stridex;
16625c28e8Piotr Jasiukajtis		y += stridey;
16725c28e8Piotr Jasiukajtis		i = 2;
16825c28e8Piotr Jasiukajtis		if (--n <= 0)
16925c28e8Piotr Jasiukajtis			break;
17025c28e8Piotr Jasiukajtis
17125c28e8Piotr JasiukajtisLOOP2: /* Get third arg, same as above. */
17225c28e8Piotr Jasiukajtis		xsb2 = HI(x);
17325c28e8Piotr Jasiukajtis		hx2 = xsb2 & ~0x80000000;
17425c28e8Piotr Jasiukajtis		if (hx2 > 0x3fe921fb)
17525c28e8Piotr Jasiukajtis		{
17625c28e8Piotr Jasiukajtis			biguns = 3;
17725c28e8Piotr Jasiukajtis			goto MEDIUM;
17825c28e8Piotr Jasiukajtis		}
17925c28e8Piotr Jasiukajtis		if (hx2 < 0x3e400000)
18025c28e8Piotr Jasiukajtis		{
18125c28e8Piotr Jasiukajtis			v = *x;
18225c28e8Piotr Jasiukajtis			*y = 1.0;
18325c28e8Piotr Jasiukajtis			x += stridex;
18425c28e8Piotr Jasiukajtis			y += stridey;
18525c28e8Piotr Jasiukajtis			i = 2;
18625c28e8Piotr Jasiukajtis			if (--n <= 0)
18725c28e8Piotr Jasiukajtis				break;
18825c28e8Piotr Jasiukajtis			goto LOOP2;
18925c28e8Piotr Jasiukajtis		}
19025c28e8Piotr Jasiukajtis		x2 = *x;
19125c28e8Piotr Jasiukajtis		py2 = y;
19225c28e8Piotr Jasiukajtis
19325c28e8Piotr Jasiukajtis		/*
19425c28e8Piotr Jasiukajtis		 * 0x3fc40000 = 5/32 ~ 0.15625
19525c28e8Piotr Jasiukajtis		 * Get msb after subtraction.  Will be 1 only if
19625c28e8Piotr Jasiukajtis		 * hx0 - 5/32 is negative.
19725c28e8Piotr Jasiukajtis		 */
19825c28e8Piotr Jasiukajtis		i = (hx0 - 0x3fc40000) >> 31;
19925c28e8Piotr Jasiukajtis		i |= ((hx1 - 0x3fc40000) >> 30) & 2;
20025c28e8Piotr Jasiukajtis		i |= ((hx2 - 0x3fc40000) >> 29) & 4;
20125c28e8Piotr Jasiukajtis		switch (i)
20225c28e8Piotr Jasiukajtis		{
20325c28e8Piotr Jasiukajtis			double		a0, a1, a2, w0, w1, w2;
20425c28e8Piotr Jasiukajtis			double		t0, t1, t2, z0, z1, z2;
20525c28e8Piotr Jasiukajtis			unsigned	j0, j1, j2;
20625c28e8Piotr Jasiukajtis
20725c28e8Piotr Jasiukajtis		case 0: /* All are > 5/32 */
20825c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
20925c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
21025c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
21125c28e8Piotr Jasiukajtis			HI(&t0) = j0;
21225c28e8Piotr Jasiukajtis			HI(&t1) = j1;
21325c28e8Piotr Jasiukajtis			HI(&t2) = j2;
21425c28e8Piotr Jasiukajtis			LO(&t0) = 0;
21525c28e8Piotr Jasiukajtis			LO(&t1) = 0;
21625c28e8Piotr Jasiukajtis			LO(&t2) = 0;
21725c28e8Piotr Jasiukajtis			x0 -= t0;
21825c28e8Piotr Jasiukajtis			x1 -= t1;
21925c28e8Piotr Jasiukajtis			x2 -= t2;
22025c28e8Piotr Jasiukajtis			z0 = x0 * x0;
22125c28e8Piotr Jasiukajtis			z1 = x1 * x1;
22225c28e8Piotr Jasiukajtis			z2 = x2 * x2;
22325c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
22425c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
22525c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
22625c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
22725c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
22825c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
22925c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
23025c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
23125c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
23225c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
23325c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
23425c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
23525c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+1]; /* cos_hi(t) */
23625c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+1];
23725c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+1];
23825c28e8Piotr Jasiukajtis			   /*   cos_lo(t)			 sin_hi(t) */
23925c28e8Piotr Jasiukajtis			t0 = __vlibm_TBL_sincos_lo[j0+1] - (__vlibm_TBL_sincos_hi[j0+xsb0]*w0 - a0*t0);
24025c28e8Piotr Jasiukajtis			t1 = __vlibm_TBL_sincos_lo[j1+1] - (__vlibm_TBL_sincos_hi[j1+xsb1]*w1 - a1*t1);
24125c28e8Piotr Jasiukajtis			t2 = __vlibm_TBL_sincos_lo[j2+1] - (__vlibm_TBL_sincos_hi[j2+xsb2]*w2 - a2*t2);
24225c28e8Piotr Jasiukajtis
24325c28e8Piotr Jasiukajtis			*py0 = a0 + t0;
24425c28e8Piotr Jasiukajtis			*py1 = a1 + t1;
24525c28e8Piotr Jasiukajtis			*py2 = a2 + t2;
24625c28e8Piotr Jasiukajtis			break;
24725c28e8Piotr Jasiukajtis
24825c28e8Piotr Jasiukajtis		case 1:
24925c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
25025c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
25125c28e8Piotr Jasiukajtis			HI(&t1) = j1;
25225c28e8Piotr Jasiukajtis			HI(&t2) = j2;
25325c28e8Piotr Jasiukajtis			LO(&t1) = 0;
25425c28e8Piotr Jasiukajtis			LO(&t2) = 0;
25525c28e8Piotr Jasiukajtis			x1 -= t1;
25625c28e8Piotr Jasiukajtis			x2 -= t2;
25725c28e8Piotr Jasiukajtis			z0 = x0 * x0;
25825c28e8Piotr Jasiukajtis			z1 = x1 * x1;
25925c28e8Piotr Jasiukajtis			z2 = x2 * x2;
26025c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[1] + z0 * poly4[1]);
26125c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
26225c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
26325c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[1] + z0 * (poly2[1] + t0));
26425c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
26525c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
26625c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
26725c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
26825c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
26925c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
27025c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+1];
27125c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+1];
27225c28e8Piotr Jasiukajtis			t1 = __vlibm_TBL_sincos_lo[j1+1] - (__vlibm_TBL_sincos_hi[j1+xsb1]*w1 - a1*t1);
27325c28e8Piotr Jasiukajtis			t2 = __vlibm_TBL_sincos_lo[j2+1] - (__vlibm_TBL_sincos_hi[j2+xsb2]*w2 - a2*t2);
27425c28e8Piotr Jasiukajtis			*py0 = one + t0;
27525c28e8Piotr Jasiukajtis			*py1 = a1 + t1;
27625c28e8Piotr Jasiukajtis			*py2 = a2 + t2;
27725c28e8Piotr Jasiukajtis			break;
27825c28e8Piotr Jasiukajtis
27925c28e8Piotr Jasiukajtis		case 2:
28025c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
28125c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
28225c28e8Piotr Jasiukajtis			HI(&t0) = j0;
28325c28e8Piotr Jasiukajtis			HI(&t2) = j2;
28425c28e8Piotr Jasiukajtis			LO(&t0) = 0;
28525c28e8Piotr Jasiukajtis			LO(&t2) = 0;
28625c28e8Piotr Jasiukajtis			x0 -= t0;
28725c28e8Piotr Jasiukajtis			x2 -= t2;
28825c28e8Piotr Jasiukajtis			z0 = x0 * x0;
28925c28e8Piotr Jasiukajtis			z1 = x1 * x1;
29025c28e8Piotr Jasiukajtis			z2 = x2 * x2;
29125c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
29225c28e8Piotr Jasiukajtis			t1 = z1 * (poly3[1] + z1 * poly4[1]);
29325c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
29425c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
29525c28e8Piotr Jasiukajtis			t1 = z1 * (poly1[1] + z1 * (poly2[1] + t1));
29625c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
29725c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
29825c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
29925c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
30025c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
30125c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+1];
30225c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+1];
30325c28e8Piotr Jasiukajtis			t0 = __vlibm_TBL_sincos_lo[j0+1] - (__vlibm_TBL_sincos_hi[j0+xsb0]*w0 - a0*t0);
30425c28e8Piotr Jasiukajtis			t2 = __vlibm_TBL_sincos_lo[j2+1] - (__vlibm_TBL_sincos_hi[j2+xsb2]*w2 - a2*t2);
30525c28e8Piotr Jasiukajtis			*py0 = a0 + t0;
30625c28e8Piotr Jasiukajtis			*py1 = one + t1;
30725c28e8Piotr Jasiukajtis			*py2 = a2 + t2;
30825c28e8Piotr Jasiukajtis			break;
30925c28e8Piotr Jasiukajtis
31025c28e8Piotr Jasiukajtis		case 3:
31125c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
31225c28e8Piotr Jasiukajtis			HI(&t2) = j2;
31325c28e8Piotr Jasiukajtis			LO(&t2) = 0;
31425c28e8Piotr Jasiukajtis			x2 -= t2;
31525c28e8Piotr Jasiukajtis			z0 = x0 * x0;
31625c28e8Piotr Jasiukajtis			z1 = x1 * x1;
31725c28e8Piotr Jasiukajtis			z2 = x2 * x2;
31825c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[1] + z0 * poly4[1]);
31925c28e8Piotr Jasiukajtis			t1 = z1 * (poly3[1] + z1 * poly4[1]);
32025c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
32125c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[1] + z0 * (poly2[1] + t0));
32225c28e8Piotr Jasiukajtis			t1 = z1 * (poly1[1] + z1 * (poly2[1] + t1));
32325c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
32425c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
32525c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
32625c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+1];
32725c28e8Piotr Jasiukajtis			t2 = __vlibm_TBL_sincos_lo[j2+1] - (__vlibm_TBL_sincos_hi[j2+xsb2]*w2 - a2*t2);
32825c28e8Piotr Jasiukajtis			*py0 = one + t0;
32925c28e8Piotr Jasiukajtis			*py1 = one + t1;
33025c28e8Piotr Jasiukajtis			*py2 = a2 + t2;
33125c28e8Piotr Jasiukajtis			break;
33225c28e8Piotr Jasiukajtis
33325c28e8Piotr Jasiukajtis		case 4:
33425c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
33525c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
33625c28e8Piotr Jasiukajtis			HI(&t0) = j0;
33725c28e8Piotr Jasiukajtis			HI(&t1) = j1;
33825c28e8Piotr Jasiukajtis			LO(&t0) = 0;
33925c28e8Piotr Jasiukajtis			LO(&t1) = 0;
34025c28e8Piotr Jasiukajtis			x0 -= t0;
34125c28e8Piotr Jasiukajtis			x1 -= t1;
34225c28e8Piotr Jasiukajtis			z0 = x0 * x0;
34325c28e8Piotr Jasiukajtis			z1 = x1 * x1;
34425c28e8Piotr Jasiukajtis			z2 = x2 * x2;
34525c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
34625c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
34725c28e8Piotr Jasiukajtis			t2 = z2 * (poly3[1] + z2 * poly4[1]);
34825c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
34925c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
35025c28e8Piotr Jasiukajtis			t2 = z2 * (poly1[1] + z2 * (poly2[1] + t2));
35125c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
35225c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
35325c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
35425c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
35525c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+1];
35625c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+1];
35725c28e8Piotr Jasiukajtis			t0 = __vlibm_TBL_sincos_lo[j0+1] - (__vlibm_TBL_sincos_hi[j0+xsb0]*w0 - a0*t0);
35825c28e8Piotr Jasiukajtis			t1 = __vlibm_TBL_sincos_lo[j1+1] - (__vlibm_TBL_sincos_hi[j1+xsb1]*w1 - a1*t1);
35925c28e8Piotr Jasiukajtis			*py0 = a0 + t0;
36025c28e8Piotr Jasiukajtis			*py1 = a1 + t1;
36125c28e8Piotr Jasiukajtis			*py2 = one + t2;
36225c28e8Piotr Jasiukajtis			break;
36325c28e8Piotr Jasiukajtis
36425c28e8Piotr Jasiukajtis		case 5:
36525c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
36625c28e8Piotr Jasiukajtis			HI(&t1) = j1;
36725c28e8Piotr Jasiukajtis			LO(&t1) = 0;
36825c28e8Piotr Jasiukajtis			x1 -= t1;
36925c28e8Piotr Jasiukajtis			z0 = x0 * x0;
37025c28e8Piotr Jasiukajtis			z1 = x1 * x1;
37125c28e8Piotr Jasiukajtis			z2 = x2 * x2;
37225c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[1] + z0 * poly4[1]);
37325c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
37425c28e8Piotr Jasiukajtis			t2 = z2 * (poly3[1] + z2 * poly4[1]);
37525c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[1] + z0 * (poly2[1] + t0));
37625c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
37725c28e8Piotr Jasiukajtis			t2 = z2 * (poly1[1] + z2 * (poly2[1] + t2));
37825c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
37925c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
38025c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+1];
38125c28e8Piotr Jasiukajtis			t1 = __vlibm_TBL_sincos_lo[j1+1] - (__vlibm_TBL_sincos_hi[j1+xsb1]*w1 - a1*t1);
38225c28e8Piotr Jasiukajtis			*py0 = one + t0;
38325c28e8Piotr Jasiukajtis			*py1 = a1 + t1;
38425c28e8Piotr Jasiukajtis			*py2 = one + t2;
38525c28e8Piotr Jasiukajtis			break;
38625c28e8Piotr Jasiukajtis
38725c28e8Piotr Jasiukajtis		case 6:
38825c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
38925c28e8Piotr Jasiukajtis			HI(&t0) = j0;
39025c28e8Piotr Jasiukajtis			LO(&t0) = 0;
39125c28e8Piotr Jasiukajtis			x0 -= t0;
39225c28e8Piotr Jasiukajtis			z0 = x0 * x0;
39325c28e8Piotr Jasiukajtis			z1 = x1 * x1;
39425c28e8Piotr Jasiukajtis			z2 = x2 * x2;
39525c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
39625c28e8Piotr Jasiukajtis			t1 = z1 * (poly3[1] + z1 * poly4[1]);
39725c28e8Piotr Jasiukajtis			t2 = z2 * (poly3[1] + z2 * poly4[1]);
39825c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
39925c28e8Piotr Jasiukajtis			t1 = z1 * (poly1[1] + z1 * (poly2[1] + t1));
40025c28e8Piotr Jasiukajtis			t2 = z2 * (poly1[1] + z2 * (poly2[1] + t2));
40125c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
40225c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
40325c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+1];
40425c28e8Piotr Jasiukajtis			t0 = __vlibm_TBL_sincos_lo[j0+1] - (__vlibm_TBL_sincos_hi[j0+xsb0]*w0 - a0*t0);
40525c28e8Piotr Jasiukajtis			*py0 = a0 + t0;
40625c28e8Piotr Jasiukajtis			*py1 = one + t1;
40725c28e8Piotr Jasiukajtis			*py2 = one + t2;
40825c28e8Piotr Jasiukajtis			break;
40925c28e8Piotr Jasiukajtis
41025c28e8Piotr Jasiukajtis		case 7: /* All are < 5/32 */
41125c28e8Piotr Jasiukajtis			z0 = x0 * x0;
41225c28e8Piotr Jasiukajtis			z1 = x1 * x1;
41325c28e8Piotr Jasiukajtis			z2 = x2 * x2;
41425c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[1] + z0 * poly4[1]);
41525c28e8Piotr Jasiukajtis			t1 = z1 * (poly3[1] + z1 * poly4[1]);
41625c28e8Piotr Jasiukajtis			t2 = z2 * (poly3[1] + z2 * poly4[1]);
41725c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[1] + z0 * (poly2[1] + t0));
41825c28e8Piotr Jasiukajtis			t1 = z1 * (poly1[1] + z1 * (poly2[1] + t1));
41925c28e8Piotr Jasiukajtis			t2 = z2 * (poly1[1] + z2 * (poly2[1] + t2));
42025c28e8Piotr Jasiukajtis			*py0 = one + t0;
42125c28e8Piotr Jasiukajtis			*py1 = one + t1;
42225c28e8Piotr Jasiukajtis			*py2 = one + t2;
42325c28e8Piotr Jasiukajtis			break;
42425c28e8Piotr Jasiukajtis		}
42525c28e8Piotr Jasiukajtis
42625c28e8Piotr Jasiukajtis		x += stridex;
42725c28e8Piotr Jasiukajtis		y += stridey;
42825c28e8Piotr Jasiukajtis		i = 0;
42925c28e8Piotr Jasiukajtis	} while (--n > 0); /* END MAIN LOOP */
43025c28e8Piotr Jasiukajtis
43125c28e8Piotr Jasiukajtis	/*
43225c28e8Piotr Jasiukajtis	 * CLEAN UP last 0, 1, or 2 elts.
43325c28e8Piotr Jasiukajtis	 */
43425c28e8Piotr Jasiukajtis	if (i > 0) /* Clean up elts at tail.  i < 3. */
43525c28e8Piotr Jasiukajtis	{
43625c28e8Piotr Jasiukajtis		double		a0, a1, w0, w1;
43725c28e8Piotr Jasiukajtis		double		t0, t1, z0, z1;
43825c28e8Piotr Jasiukajtis		unsigned	j0, j1;
43925c28e8Piotr Jasiukajtis
44025c28e8Piotr Jasiukajtis		if (i > 1)
44125c28e8Piotr Jasiukajtis		{
44225c28e8Piotr Jasiukajtis			if (hx1 < 0x3fc40000)
44325c28e8Piotr Jasiukajtis			{
44425c28e8Piotr Jasiukajtis				z1 = x1 * x1;
44525c28e8Piotr Jasiukajtis				t1 = z1 * (poly3[1] + z1 * poly4[1]);
44625c28e8Piotr Jasiukajtis				t1 = z1 * (poly1[1] + z1 * (poly2[1] + t1));
44725c28e8Piotr Jasiukajtis				t1 = one + t1;
44825c28e8Piotr Jasiukajtis				*py1 = t1;
44925c28e8Piotr Jasiukajtis			}
45025c28e8Piotr Jasiukajtis			else
45125c28e8Piotr Jasiukajtis			{
45225c28e8Piotr Jasiukajtis				j1 = (xsb1 + 0x4000) & 0xffff8000;
45325c28e8Piotr Jasiukajtis				HI(&t1) = j1;
45425c28e8Piotr Jasiukajtis				LO(&t1) = 0;
45525c28e8Piotr Jasiukajtis				x1 -= t1;
45625c28e8Piotr Jasiukajtis				z1 = x1 * x1;
45725c28e8Piotr Jasiukajtis				t1 = z1 * (qq1 + z1 * qq2);
45825c28e8Piotr Jasiukajtis				w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
45925c28e8Piotr Jasiukajtis				j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
46025c28e8Piotr Jasiukajtis				xsb1 = (xsb1 >> 30) & 2;
46125c28e8Piotr Jasiukajtis				a1 = __vlibm_TBL_sincos_hi[j1+1];
46225c28e8Piotr Jasiukajtis				t1 = __vlibm_TBL_sincos_lo[j1+1]
46325c28e8Piotr Jasiukajtis					- (__vlibm_TBL_sincos_hi[j1+xsb1]*w1 - a1*t1);
46425c28e8Piotr Jasiukajtis				*py1 = a1 + t1;
46525c28e8Piotr Jasiukajtis			}
46625c28e8Piotr Jasiukajtis		}
46725c28e8Piotr Jasiukajtis		if (hx0 < 0x3fc40000)
46825c28e8Piotr Jasiukajtis		{
46925c28e8Piotr Jasiukajtis			z0 = x0 * x0;
47025c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[1] + z0 * poly4[1]);
47125c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[1] + z0 * (poly2[1] + t0));
47225c28e8Piotr Jasiukajtis			t0 = one + t0;
47325c28e8Piotr Jasiukajtis			*py0 = t0;
47425c28e8Piotr Jasiukajtis		}
47525c28e8Piotr Jasiukajtis		else
47625c28e8Piotr Jasiukajtis		{
47725c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
47825c28e8Piotr Jasiukajtis			HI(&t0) = j0;
47925c28e8Piotr Jasiukajtis			LO(&t0) = 0;
48025c28e8Piotr Jasiukajtis			x0 -= t0;
48125c28e8Piotr Jasiukajtis			z0 = x0 * x0;
48225c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
48325c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
48425c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
48525c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
48625c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+1];
48725c28e8Piotr Jasiukajtis			t0 = __vlibm_TBL_sincos_lo[j0+1] - (__vlibm_TBL_sincos_hi[j0+xsb0]*w0 - a0*t0);
48825c28e8Piotr Jasiukajtis			*py0 = a0 + t0;
48925c28e8Piotr Jasiukajtis		}
49025c28e8Piotr Jasiukajtis	} /* END CLEAN UP */
49125c28e8Piotr Jasiukajtis
49225c28e8Piotr Jasiukajtis	return;
49325c28e8Piotr Jasiukajtis
49425c28e8Piotr Jasiukajtis	/*
49525c28e8Piotr Jasiukajtis	 * Take care of BIGUNS.
49625c28e8Piotr Jasiukajtis	 *
49725c28e8Piotr Jasiukajtis	 * We have jumped here in the middle of processing after having
49825c28e8Piotr Jasiukajtis	 * encountered a medium range argument.  Therefore things are in a
49925c28e8Piotr Jasiukajtis	 * bit of a tizzy.
50025c28e8Piotr Jasiukajtis	 */
50125c28e8Piotr Jasiukajtis
50225c28e8Piotr JasiukajtisMEDIUM:
50325c28e8Piotr Jasiukajtis
50425c28e8Piotr Jasiukajtis	x0_or_one[1] = 1.0;
50525c28e8Piotr Jasiukajtis	x1_or_one[1] = 1.0;
50625c28e8Piotr Jasiukajtis	x2_or_one[1] = 1.0;
50725c28e8Piotr Jasiukajtis	x0_or_one[3] = -1.0;
50825c28e8Piotr Jasiukajtis	x1_or_one[3] = -1.0;
50925c28e8Piotr Jasiukajtis	x2_or_one[3] = -1.0;
51025c28e8Piotr Jasiukajtis	y0_or_zero[1] = 0.0;
51125c28e8Piotr Jasiukajtis	y1_or_zero[1] = 0.0;
51225c28e8Piotr Jasiukajtis	y2_or_zero[1] = 0.0;
51325c28e8Piotr Jasiukajtis	y0_or_zero[3] = 0.0;
51425c28e8Piotr Jasiukajtis	y1_or_zero[3] = 0.0;
51525c28e8Piotr Jasiukajtis	y2_or_zero[3] = 0.0;
51625c28e8Piotr Jasiukajtis
51725c28e8Piotr Jasiukajtis	if (biguns == 3)
51825c28e8Piotr Jasiukajtis	{
51925c28e8Piotr Jasiukajtis		biguns = 0;
52025c28e8Piotr Jasiukajtis		xsb0 = xsb0 >> 31;
52125c28e8Piotr Jasiukajtis		xsb1 = xsb1 >> 31;
52225c28e8Piotr Jasiukajtis		goto loop2;
52325c28e8Piotr Jasiukajtis	}
52425c28e8Piotr Jasiukajtis	else if (biguns == 2)
52525c28e8Piotr Jasiukajtis	{
52625c28e8Piotr Jasiukajtis		xsb0 = xsb0 >> 31;
52725c28e8Piotr Jasiukajtis		biguns = 0;
52825c28e8Piotr Jasiukajtis		goto loop1;
52925c28e8Piotr Jasiukajtis	}
53025c28e8Piotr Jasiukajtis	biguns = 0;
53125c28e8Piotr Jasiukajtis
53225c28e8Piotr Jasiukajtis	do
53325c28e8Piotr Jasiukajtis	{
53425c28e8Piotr Jasiukajtis		double		fn0, fn1, fn2, a0, a1, a2, w0, w1, w2, y0, y1, y2;
53525c28e8Piotr Jasiukajtis		unsigned	hx;
53625c28e8Piotr Jasiukajtis		int			n0, n1, n2;
53725c28e8Piotr Jasiukajtis
53825c28e8Piotr Jasiukajtis		/*
53925c28e8Piotr Jasiukajtis		 * Find 3 more to work on: Not already done, not too big.
54025c28e8Piotr Jasiukajtis		 */
54125c28e8Piotr Jasiukajtis
54225c28e8Piotr Jasiukajtisloop0:
54325c28e8Piotr Jasiukajtis		hx = HI(x);
54425c28e8Piotr Jasiukajtis		xsb0 = hx >> 31;
54525c28e8Piotr Jasiukajtis		hx &= ~0x80000000;
54625c28e8Piotr Jasiukajtis		if (hx > 0x413921fb) /* (1.6471e+06) Too big: leave it. */
54725c28e8Piotr Jasiukajtis		{
54825c28e8Piotr Jasiukajtis			if (hx >= 0x7ff00000) /* Inf or NaN */
54925c28e8Piotr Jasiukajtis			{
55025c28e8Piotr Jasiukajtis				x0 = *x;
55125c28e8Piotr Jasiukajtis				*y = x0 - x0;
55225c28e8Piotr Jasiukajtis			}
55325c28e8Piotr Jasiukajtis			else
55425c28e8Piotr Jasiukajtis				biguns = 1;
55525c28e8Piotr Jasiukajtis			x += stridex;
55625c28e8Piotr Jasiukajtis			y += stridey;
55725c28e8Piotr Jasiukajtis			i = 0;
55825c28e8Piotr Jasiukajtis			if (--n <= 0)
55925c28e8Piotr Jasiukajtis				break;
56025c28e8Piotr Jasiukajtis			goto loop0;
56125c28e8Piotr Jasiukajtis		}
56225c28e8Piotr Jasiukajtis		x0 = *x;
56325c28e8Piotr Jasiukajtis		py0 = y;
56425c28e8Piotr Jasiukajtis		x += stridex;
56525c28e8Piotr Jasiukajtis		y += stridey;
56625c28e8Piotr Jasiukajtis		i = 1;
56725c28e8Piotr Jasiukajtis		if (--n <= 0)
56825c28e8Piotr Jasiukajtis			break;
56925c28e8Piotr Jasiukajtis
57025c28e8Piotr Jasiukajtisloop1:
57125c28e8Piotr Jasiukajtis		hx = HI(x);
57225c28e8Piotr Jasiukajtis		xsb1 = hx >> 31;
57325c28e8Piotr Jasiukajtis		hx &= ~0x80000000;
57425c28e8Piotr Jasiukajtis		if (hx > 0x413921fb)
57525c28e8Piotr Jasiukajtis		{
57625c28e8Piotr Jasiukajtis			if (hx >= 0x7ff00000)
57725c28e8Piotr Jasiukajtis			{
57825c28e8Piotr Jasiukajtis				x1 = *x;
57925c28e8Piotr Jasiukajtis				*y = x1 - x1;
58025c28e8Piotr Jasiukajtis			}
58125c28e8Piotr Jasiukajtis			else
58225c28e8Piotr Jasiukajtis				biguns = 1;
58325c28e8Piotr Jasiukajtis			x += stridex;
58425c28e8Piotr Jasiukajtis			y += stridey;
58525c28e8Piotr Jasiukajtis			i = 1;
58625c28e8Piotr Jasiukajtis			if (--n <= 0)
58725c28e8Piotr Jasiukajtis				break;
58825c28e8Piotr Jasiukajtis			goto loop1;
58925c28e8Piotr Jasiukajtis		}
59025c28e8Piotr Jasiukajtis		x1 = *x;
59125c28e8Piotr Jasiukajtis		py1 = y;
59225c28e8Piotr Jasiukajtis		x += stridex;
59325c28e8Piotr Jasiukajtis		y += stridey;
59425c28e8Piotr Jasiukajtis		i = 2;
59525c28e8Piotr Jasiukajtis		if (--n <= 0)
59625c28e8Piotr Jasiukajtis			break;
59725c28e8Piotr Jasiukajtis
59825c28e8Piotr Jasiukajtisloop2:
59925c28e8Piotr Jasiukajtis		hx = HI(x);
60025c28e8Piotr Jasiukajtis		xsb2 = hx >> 31;
60125c28e8Piotr Jasiukajtis		hx &= ~0x80000000;
60225c28e8Piotr Jasiukajtis		if (hx > 0x413921fb)
60325c28e8Piotr Jasiukajtis		{
60425c28e8Piotr Jasiukajtis			if (hx >= 0x7ff00000)
60525c28e8Piotr Jasiukajtis			{
60625c28e8Piotr Jasiukajtis				x2 = *x;
60725c28e8Piotr Jasiukajtis				*y = x2 - x2;
60825c28e8Piotr Jasiukajtis			}
60925c28e8Piotr Jasiukajtis			else
61025c28e8Piotr Jasiukajtis				biguns = 1;
61125c28e8Piotr Jasiukajtis			x += stridex;
61225c28e8Piotr Jasiukajtis			y += stridey;
61325c28e8Piotr Jasiukajtis			i = 2;
61425c28e8Piotr Jasiukajtis			if (--n <= 0)
61525c28e8Piotr Jasiukajtis				break;
61625c28e8Piotr Jasiukajtis			goto loop2;
61725c28e8Piotr Jasiukajtis		}
61825c28e8Piotr Jasiukajtis		x2 = *x;
61925c28e8Piotr Jasiukajtis		py2 = y;
62025c28e8Piotr Jasiukajtis
62125c28e8Piotr Jasiukajtis		n0 = (int) (x0 * invpio2 + half[xsb0]);
62225c28e8Piotr Jasiukajtis		n1 = (int) (x1 * invpio2 + half[xsb1]);
62325c28e8Piotr Jasiukajtis		n2 = (int) (x2 * invpio2 + half[xsb2]);
62425c28e8Piotr Jasiukajtis		fn0 = (double) n0;
62525c28e8Piotr Jasiukajtis		fn1 = (double) n1;
62625c28e8Piotr Jasiukajtis		fn2 = (double) n2;
62725c28e8Piotr Jasiukajtis		n0 = (n0 + 1) & 3; /* Add 1 (before the mod) to make sin into cos */
62825c28e8Piotr Jasiukajtis		n1 = (n1 + 1) & 3;
62925c28e8Piotr Jasiukajtis		n2 = (n2 + 1) & 3;
63025c28e8Piotr Jasiukajtis		a0 = x0 - fn0 * pio2_1;
63125c28e8Piotr Jasiukajtis		a1 = x1 - fn1 * pio2_1;
63225c28e8Piotr Jasiukajtis		a2 = x2 - fn2 * pio2_1;
63325c28e8Piotr Jasiukajtis		w0 = fn0 * pio2_2;
63425c28e8Piotr Jasiukajtis		w1 = fn1 * pio2_2;
63525c28e8Piotr Jasiukajtis		w2 = fn2 * pio2_2;
63625c28e8Piotr Jasiukajtis		x0 = a0 - w0;
63725c28e8Piotr Jasiukajtis		x1 = a1 - w1;
63825c28e8Piotr Jasiukajtis		x2 = a2 - w2;
63925c28e8Piotr Jasiukajtis		y0 = (a0 - x0) - w0;
64025c28e8Piotr Jasiukajtis		y1 = (a1 - x1) - w1;
64125c28e8Piotr Jasiukajtis		y2 = (a2 - x2) - w2;
64225c28e8Piotr Jasiukajtis		a0 = x0;
64325c28e8Piotr Jasiukajtis		a1 = x1;
64425c28e8Piotr Jasiukajtis		a2 = x2;
64525c28e8Piotr Jasiukajtis		w0 = fn0 * pio2_3 - y0;
64625c28e8Piotr Jasiukajtis		w1 = fn1 * pio2_3 - y1;
64725c28e8Piotr Jasiukajtis		w2 = fn2 * pio2_3 - y2;
64825c28e8Piotr Jasiukajtis		x0 = a0 - w0;
64925c28e8Piotr Jasiukajtis		x1 = a1 - w1;
65025c28e8Piotr Jasiukajtis		x2 = a2 - w2;
65125c28e8Piotr Jasiukajtis		y0 = (a0 - x0) - w0;
65225c28e8Piotr Jasiukajtis		y1 = (a1 - x1) - w1;
65325c28e8Piotr Jasiukajtis		y2 = (a2 - x2) - w2;
65425c28e8Piotr Jasiukajtis		a0 = x0;
65525c28e8Piotr Jasiukajtis		a1 = x1;
65625c28e8Piotr Jasiukajtis		a2 = x2;
65725c28e8Piotr Jasiukajtis		w0 = fn0 * pio2_3t - y0;
65825c28e8Piotr Jasiukajtis		w1 = fn1 * pio2_3t - y1;
65925c28e8Piotr Jasiukajtis		w2 = fn2 * pio2_3t - y2;
66025c28e8Piotr Jasiukajtis		x0 = a0 - w0;
66125c28e8Piotr Jasiukajtis		x1 = a1 - w1;
66225c28e8Piotr Jasiukajtis		x2 = a2 - w2;
66325c28e8Piotr Jasiukajtis		y0 = (a0 - x0) - w0;
66425c28e8Piotr Jasiukajtis		y1 = (a1 - x1) - w1;
66525c28e8Piotr Jasiukajtis		y2 = (a2 - x2) - w2;
66625c28e8Piotr Jasiukajtis		xsb0 = HI(&x0);
66725c28e8Piotr Jasiukajtis		i = ((xsb0 & ~0x80000000) - thresh[n0&1]) >> 31;
66825c28e8Piotr Jasiukajtis		xsb1 = HI(&x1);
66925c28e8Piotr Jasiukajtis		i |= (((xsb1 & ~0x80000000) - thresh[n1&1]) >> 30) & 2;
67025c28e8Piotr Jasiukajtis		xsb2 = HI(&x2);
67125c28e8Piotr Jasiukajtis		i |= (((xsb2 & ~0x80000000) - thresh[n2&1]) >> 29) & 4;
67225c28e8Piotr Jasiukajtis		switch (i)
67325c28e8Piotr Jasiukajtis		{
67425c28e8Piotr Jasiukajtis			double		t0, t1, t2, z0, z1, z2;
67525c28e8Piotr Jasiukajtis			unsigned	j0, j1, j2;
67625c28e8Piotr Jasiukajtis
67725c28e8Piotr Jasiukajtis		case 0:
67825c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
67925c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
68025c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
68125c28e8Piotr Jasiukajtis			HI(&t0) = j0;
68225c28e8Piotr Jasiukajtis			HI(&t1) = j1;
68325c28e8Piotr Jasiukajtis			HI(&t2) = j2;
68425c28e8Piotr Jasiukajtis			LO(&t0) = 0;
68525c28e8Piotr Jasiukajtis			LO(&t1) = 0;
68625c28e8Piotr Jasiukajtis			LO(&t2) = 0;
68725c28e8Piotr Jasiukajtis			x0 = (x0 - t0) + y0;
68825c28e8Piotr Jasiukajtis			x1 = (x1 - t1) + y1;
68925c28e8Piotr Jasiukajtis			x2 = (x2 - t2) + y2;
69025c28e8Piotr Jasiukajtis			z0 = x0 * x0;
69125c28e8Piotr Jasiukajtis			z1 = x1 * x1;
69225c28e8Piotr Jasiukajtis			z2 = x2 * x2;
69325c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
69425c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
69525c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
69625c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
69725c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
69825c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
69925c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
70025c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
70125c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
70225c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
70325c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
70425c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
70525c28e8Piotr Jasiukajtis			n0 ^= (xsb0 & ~(n0 << 1));
70625c28e8Piotr Jasiukajtis			n1 ^= (xsb1 & ~(n1 << 1));
70725c28e8Piotr Jasiukajtis			n2 ^= (xsb2 & ~(n2 << 1));
70825c28e8Piotr Jasiukajtis			xsb0 |= 1;
70925c28e8Piotr Jasiukajtis			xsb1 |= 1;
71025c28e8Piotr Jasiukajtis			xsb2 |= 1;
71125c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+n0];
71225c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+n1];
71325c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+n2];
71425c28e8Piotr Jasiukajtis			t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0];
71525c28e8Piotr Jasiukajtis			t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1];
71625c28e8Piotr Jasiukajtis			t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2];
71725c28e8Piotr Jasiukajtis			*py0 = ( a0 + t0 );
71825c28e8Piotr Jasiukajtis			*py1 = ( a1 + t1 );
71925c28e8Piotr Jasiukajtis			*py2 = ( a2 + t2 );
72025c28e8Piotr Jasiukajtis			break;
72125c28e8Piotr Jasiukajtis
72225c28e8Piotr Jasiukajtis		case 1:
72325c28e8Piotr Jasiukajtis			j0 = n0 & 1;
72425c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
72525c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
72625c28e8Piotr Jasiukajtis			HI(&t1) = j1;
72725c28e8Piotr Jasiukajtis			HI(&t2) = j2;
72825c28e8Piotr Jasiukajtis			LO(&t1) = 0;
72925c28e8Piotr Jasiukajtis			LO(&t2) = 0;
73025c28e8Piotr Jasiukajtis			x0_or_one[0] = x0;
73125c28e8Piotr Jasiukajtis			x0_or_one[2] = -x0;
73225c28e8Piotr Jasiukajtis			y0_or_zero[0] = y0;
73325c28e8Piotr Jasiukajtis			y0_or_zero[2] = -y0;
73425c28e8Piotr Jasiukajtis			x1 = (x1 - t1) + y1;
73525c28e8Piotr Jasiukajtis			x2 = (x2 - t2) + y2;
73625c28e8Piotr Jasiukajtis			z0 = x0 * x0;
73725c28e8Piotr Jasiukajtis			z1 = x1 * x1;
73825c28e8Piotr Jasiukajtis			z2 = x2 * x2;
73925c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[j0] + z0 * poly4[j0]);
74025c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
74125c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
74225c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0));
74325c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
74425c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
74525c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
74625c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
74725c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
74825c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
74925c28e8Piotr Jasiukajtis			n1 ^= (xsb1 & ~(n1 << 1));
75025c28e8Piotr Jasiukajtis			n2 ^= (xsb2 & ~(n2 << 1));
75125c28e8Piotr Jasiukajtis			xsb1 |= 1;
75225c28e8Piotr Jasiukajtis			xsb2 |= 1;
75325c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+n1];
75425c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+n2];
75525c28e8Piotr Jasiukajtis			t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0);
75625c28e8Piotr Jasiukajtis			t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1];
75725c28e8Piotr Jasiukajtis			t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2];
75825c28e8Piotr Jasiukajtis			*py0 = t0;
75925c28e8Piotr Jasiukajtis			*py1 = ( a1 + t1 );
76025c28e8Piotr Jasiukajtis			*py2 = ( a2 + t2 );
76125c28e8Piotr Jasiukajtis			break;
76225c28e8Piotr Jasiukajtis
76325c28e8Piotr Jasiukajtis		case 2:
76425c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
76525c28e8Piotr Jasiukajtis			j1 = n1 & 1;
76625c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
76725c28e8Piotr Jasiukajtis			HI(&t0) = j0;
76825c28e8Piotr Jasiukajtis			HI(&t2) = j2;
76925c28e8Piotr Jasiukajtis			LO(&t0) = 0;
77025c28e8Piotr Jasiukajtis			LO(&t2) = 0;
77125c28e8Piotr Jasiukajtis			x1_or_one[0] = x1;
77225c28e8Piotr Jasiukajtis			x1_or_one[2] = -x1;
77325c28e8Piotr Jasiukajtis			x0 = (x0 - t0) + y0;
77425c28e8Piotr Jasiukajtis			y1_or_zero[0] = y1;
77525c28e8Piotr Jasiukajtis			y1_or_zero[2] = -y1;
77625c28e8Piotr Jasiukajtis			x2 = (x2 - t2) + y2;
77725c28e8Piotr Jasiukajtis			z0 = x0 * x0;
77825c28e8Piotr Jasiukajtis			z1 = x1 * x1;
77925c28e8Piotr Jasiukajtis			z2 = x2 * x2;
78025c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
78125c28e8Piotr Jasiukajtis			t1 = z1 * (poly3[j1] + z1 * poly4[j1]);
78225c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
78325c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
78425c28e8Piotr Jasiukajtis			t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1));
78525c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
78625c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
78725c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
78825c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
78925c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
79025c28e8Piotr Jasiukajtis			n0 ^= (xsb0 & ~(n0 << 1));
79125c28e8Piotr Jasiukajtis			n2 ^= (xsb2 & ~(n2 << 1));
79225c28e8Piotr Jasiukajtis			xsb0 |= 1;
79325c28e8Piotr Jasiukajtis			xsb2 |= 1;
79425c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+n0];
79525c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+n2];
79625c28e8Piotr Jasiukajtis			t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0];
79725c28e8Piotr Jasiukajtis			t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1);
79825c28e8Piotr Jasiukajtis			t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2];
79925c28e8Piotr Jasiukajtis			*py0 = ( a0 + t0 );
80025c28e8Piotr Jasiukajtis			*py1 = t1;
80125c28e8Piotr Jasiukajtis			*py2 = ( a2 + t2 );
80225c28e8Piotr Jasiukajtis			break;
80325c28e8Piotr Jasiukajtis
80425c28e8Piotr Jasiukajtis		case 3:
80525c28e8Piotr Jasiukajtis			j0 = n0 & 1;
80625c28e8Piotr Jasiukajtis			j1 = n1 & 1;
80725c28e8Piotr Jasiukajtis			j2 = (xsb2 + 0x4000) & 0xffff8000;
80825c28e8Piotr Jasiukajtis			HI(&t2) = j2;
80925c28e8Piotr Jasiukajtis			LO(&t2) = 0;
81025c28e8Piotr Jasiukajtis			x0_or_one[0] = x0;
81125c28e8Piotr Jasiukajtis			x0_or_one[2] = -x0;
81225c28e8Piotr Jasiukajtis			x1_or_one[0] = x1;
81325c28e8Piotr Jasiukajtis			x1_or_one[2] = -x1;
81425c28e8Piotr Jasiukajtis			y0_or_zero[0] = y0;
81525c28e8Piotr Jasiukajtis			y0_or_zero[2] = -y0;
81625c28e8Piotr Jasiukajtis			y1_or_zero[0] = y1;
81725c28e8Piotr Jasiukajtis			y1_or_zero[2] = -y1;
81825c28e8Piotr Jasiukajtis			x2 = (x2 - t2) + y2;
81925c28e8Piotr Jasiukajtis			z0 = x0 * x0;
82025c28e8Piotr Jasiukajtis			z1 = x1 * x1;
82125c28e8Piotr Jasiukajtis			z2 = x2 * x2;
82225c28e8Piotr Jasiukajtis			t0 = z0 * (poly3[j0] + z0 * poly4[j0]);
82325c28e8Piotr Jasiukajtis			t1 = z1 * (poly3[j1] + z1 * poly4[j1]);
82425c28e8Piotr Jasiukajtis			t2 = z2 * (qq1 + z2 * qq2);
82525c28e8Piotr Jasiukajtis			t0 = z0 * (poly1[j0] + z0 * (poly2[j0] + t0));
82625c28e8Piotr Jasiukajtis			t1 = z1 * (poly1[j1] + z1 * (poly2[j1] + t1));
82725c28e8Piotr Jasiukajtis			w2 = x2 * (one + z2 * (pp1 + z2 * pp2));
82825c28e8Piotr Jasiukajtis			j2 = (((j2 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
82925c28e8Piotr Jasiukajtis			xsb2 = (xsb2 >> 30) & 2;
83025c28e8Piotr Jasiukajtis			n2 ^= (xsb2 & ~(n2 << 1));
83125c28e8Piotr Jasiukajtis			xsb2 |= 1;
83225c28e8Piotr Jasiukajtis			a2 = __vlibm_TBL_sincos_hi[j2+n2];
83325c28e8Piotr Jasiukajtis			t0 = x0_or_one[n0] + (y0_or_zero[n0] + x0_or_one[n0] * t0);
83425c28e8Piotr Jasiukajtis			t1 = x1_or_one[n1] + (y1_or_zero[n1] + x1_or_one[n1] * t1);
83525c28e8Piotr Jasiukajtis			t2 = (__vlibm_TBL_sincos_hi[j2+((n2+xsb2)&3)] * w2 + a2 * t2) + __vlibm_TBL_sincos_lo[j2+n2];
83625c28e8Piotr Jasiukajtis			*py0 = t0;
83725c28e8Piotr Jasiukajtis			*py1 = t1;
83825c28e8Piotr Jasiukajtis			*py2 = ( a2 + t2 );
83925c28e8Piotr Jasiukajtis			break;
84025c28e8Piotr Jasiukajtis
84125c28e8Piotr Jasiukajtis		case 4:
84225c28e8Piotr Jasiukajtis			j0 = (xsb0 + 0x4000) & 0xffff8000;
84325c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
84425c28e8Piotr Jasiukajtis			j2 = n2 & 1;
84525c28e8Piotr Jasiukajtis			HI(&t0) = j0;
84625c28e8Piotr Jasiukajtis			HI(&t1) = j1;
84725c28e8Piotr Jasiukajtis			LO(&t0) = 0;
84825c28e8Piotr Jasiukajtis			LO(&t1) = 0;
84925c28e8Piotr Jasiukajtis			x2_or_one[0] = x2;
85025c28e8Piotr Jasiukajtis			x2_or_one[2] = -x2;
85125c28e8Piotr Jasiukajtis			x0 = (x0 - t0) + y0;
85225c28e8Piotr Jasiukajtis			x1 = (x1 - t1) + y1;
85325c28e8Piotr Jasiukajtis			y2_or_zero[0] = y2;
85425c28e8Piotr Jasiukajtis			y2_or_zero[2] = -y2;
85525c28e8Piotr Jasiukajtis			z0 = x0 * x0;
85625c28e8Piotr Jasiukajtis			z1 = x1 * x1;
85725c28e8Piotr Jasiukajtis			z2 = x2 * x2;
85825c28e8Piotr Jasiukajtis			t0 = z0 * (qq1 + z0 * qq2);
85925c28e8Piotr Jasiukajtis			t1 = z1 * (qq1 + z1 * qq2);
86025c28e8Piotr Jasiukajtis			t2 = z2 * (poly3[j2] + z2 * poly4[j2]);
86125c28e8Piotr Jasiukajtis			w0 = x0 * (one + z0 * (pp1 + z0 * pp2));
86225c28e8Piotr Jasiukajtis			w1 = x1 * (one + z1 * (pp1 + z1 * pp2));
86325c28e8Piotr Jasiukajtis			t2 = z2 * (poly1[j2] + z2 * (poly2[j2] + t2));
86425c28e8Piotr Jasiukajtis			j0 = (((j0 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
86525c28e8Piotr Jasiukajtis			j1 = (((j1 & ~0x80000000) - 0x3fc40000) >> 13) & ~0x3;
86625c28e8Piotr Jasiukajtis			xsb0 = (xsb0 >> 30) & 2;
86725c28e8Piotr Jasiukajtis			xsb1 = (xsb1 >> 30) & 2;
86825c28e8Piotr Jasiukajtis			n0 ^= (xsb0 & ~(n0 << 1));
86925c28e8Piotr Jasiukajtis			n1 ^= (xsb1 & ~(n1 << 1));
87025c28e8Piotr Jasiukajtis			xsb0 |= 1;
87125c28e8Piotr Jasiukajtis			xsb1 |= 1;
87225c28e8Piotr Jasiukajtis			a0 = __vlibm_TBL_sincos_hi[j0+n0];
87325c28e8Piotr Jasiukajtis			a1 = __vlibm_TBL_sincos_hi[j1+n1];
87425c28e8Piotr Jasiukajtis			t0 = (__vlibm_TBL_sincos_hi[j0+((n0+xsb0)&3)] * w0 + a0 * t0) + __vlibm_TBL_sincos_lo[j0+n0];
87525c28e8Piotr Jasiukajtis			t1 = (__vlibm_TBL_sincos_hi[j1+((n1+xsb1)&3)] * w1 + a1 * t1) + __vlibm_TBL_sincos_lo[j1+n1];
87625c28e8Piotr Jasiukajtis			t2 = x2_or_one[n2] + (y2_or_zero[n2] + x2_or_one[n2] * t2);
87725c28e8Piotr Jasiukajtis			*py0 = ( a0 + t0 );
87825c28e8Piotr Jasiukajtis			*py1 = ( a1 + t1 );
87925c28e8Piotr Jasiukajtis			*py2 = t2;
88025c28e8Piotr Jasiukajtis			break;
88125c28e8Piotr Jasiukajtis
88225c28e8Piotr Jasiukajtis		case 5:
88325c28e8Piotr Jasiukajtis			j0 = n0 & 1;
88425c28e8Piotr Jasiukajtis			j1 = (xsb1 + 0x4000) & 0xffff8000;
88525c28e8Piotr Jasiukajtis			j2 = n2 & 1;
88625c28e8Piotr Jasiukajtis			HI(&t1) = j1;
88725c28e8Piotr Jasiukajtis			LO(&t1) = 0;
88825c28e8Piotr Jasiukajtis			x0_or_one[0] = x0;
88925c28e8Piotr Jasiukajtis			x0_or_one[2] = -x0;
89025c28e8Piotr Jasiukajtis			x2_or_one[0] = x2;
89125c28e8Piotr Jasiukajtis			x2_or_one[2] = -x2;
89225c28e8Piotr Jasiukajtis			y0_or_zero[0] = y0;
89325c28e8Piotr Jasiukajtis			y0_or_zero[2] = -y0;
89425c28e8Piotr Jasiukajtis			x1 = (x1 - t1) + y1;
895