1 files changed, 290 insertions, 309 deletions

diff --git a/arch/x86/math-emu/poly_sin.c b/arch/x86/math-emu/poly_sin.c
index a36313fb06f1..b862039c728e 100644
--- a/arch/x86/math-emu/poly_sin.c
+++ b/arch/x86/math-emu/poly_sin.c
@@ -11,7 +11,6 @@
  |                                                                           |
  +---------------------------------------------------------------------------*/
 
-
 #include "exception.h"
 #include "reg_constant.h"
 #include "fpu_emu.h"
@@ -19,379 +18,361 @@
 #include "control_w.h"
 #include "poly.h"
 
-
 #define	N_COEFF_P	4
 #define	N_COEFF_N	4
 
-static const unsigned long long pos_terms_l[N_COEFF_P] =
-{
-  0xaaaaaaaaaaaaaaabLL,
-  0x00d00d00d00cf906LL,
-  0x000006b99159a8bbLL,
-  0x000000000d7392e6LL
+static const unsigned long long pos_terms_l[N_COEFF_P] = {
+	0xaaaaaaaaaaaaaaabLL,
+	0x00d00d00d00cf906LL,
+	0x000006b99159a8bbLL,
+	0x000000000d7392e6LL
 };
 
-static const unsigned long long neg_terms_l[N_COEFF_N] =
-{
-  0x2222222222222167LL,
-  0x0002e3bc74aab624LL,
-  0x0000000b09229062LL,
-  0x00000000000c7973LL
+static const unsigned long long neg_terms_l[N_COEFF_N] = {
+	0x2222222222222167LL,
+	0x0002e3bc74aab624LL,
+	0x0000000b09229062LL,
+	0x00000000000c7973LL
 };
 
-
-
 #define	N_COEFF_PH	4
 #define	N_COEFF_NH	4
-static const unsigned long long pos_terms_h[N_COEFF_PH] =
-{
-  0x0000000000000000LL,
-  0x05b05b05b05b0406LL,
-  0x000049f93edd91a9LL,
-  0x00000000c9c9ed62LL
+static const unsigned long long pos_terms_h[N_COEFF_PH] = {
+	0x0000000000000000LL,
+	0x05b05b05b05b0406LL,
+	0x000049f93edd91a9LL,
+	0x00000000c9c9ed62LL
 };
 
-static const unsigned long long neg_terms_h[N_COEFF_NH] =
-{
-  0xaaaaaaaaaaaaaa98LL,
-  0x001a01a01a019064LL,
-  0x0000008f76c68a77LL,
-  0x0000000000d58f5eLL
+static const unsigned long long neg_terms_h[N_COEFF_NH] = {
+	0xaaaaaaaaaaaaaa98LL,
+	0x001a01a01a019064LL,
+	0x0000008f76c68a77LL,
+	0x0000000000d58f5eLL
 };
 
-
 /*--- poly_sine() -----------------------------------------------------------+
  |                                                                           |
  +---------------------------------------------------------------------------*/
-void	poly_sine(FPU_REG *st0_ptr)
+void poly_sine(FPU_REG *st0_ptr)
 {
-  int                 exponent, echange;
-  Xsig                accumulator, argSqrd, argTo4;
-  unsigned long       fix_up, adj;
-  unsigned long long  fixed_arg;
-  FPU_REG	      result;
+	int exponent, echange;
+	Xsig accumulator, argSqrd, argTo4;
+	unsigned long fix_up, adj;
+	unsigned long long fixed_arg;
+	FPU_REG result;
 
-  exponent = exponent(st0_ptr);
+	exponent = exponent(st0_ptr);
 
-  accumulator.lsw = accumulator.midw = accumulator.msw = 0;
+	accumulator.lsw = accumulator.midw = accumulator.msw = 0;
 
-  /* Split into two ranges, for arguments below and above 1.0 */
-  /* The boundary between upper and lower is approx 0.88309101259 */
-  if ( (exponent < -1) || ((exponent == -1) && (st0_ptr->sigh <= 0xe21240aa)) )
-    {
-      /* The argument is <= 0.88309101259 */
+	/* Split into two ranges, for arguments below and above 1.0 */
+	/* The boundary between upper and lower is approx 0.88309101259 */
+	if ((exponent < -1)
+	    || ((exponent == -1) && (st0_ptr->sigh <= 0xe21240aa))) {
+		/* The argument is <= 0.88309101259 */
+
+		argSqrd.msw = st0_ptr->sigh;
+		argSqrd.midw = st0_ptr->sigl;
+		argSqrd.lsw = 0;
+		mul64_Xsig(&argSqrd, &significand(st0_ptr));
+		shr_Xsig(&argSqrd, 2 * (-1 - exponent));
+		argTo4.msw = argSqrd.msw;
+		argTo4.midw = argSqrd.midw;
+		argTo4.lsw = argSqrd.lsw;
+		mul_Xsig_Xsig(&argTo4, &argTo4);
 
-      argSqrd.msw = st0_ptr->sigh; argSqrd.midw = st0_ptr->sigl; argSqrd.lsw = 0;
-      mul64_Xsig(&argSqrd, &significand(st0_ptr));
-      shr_Xsig(&argSqrd, 2*(-1-exponent));
-      argTo4.msw = argSqrd.msw; argTo4.midw = argSqrd.midw;
-      argTo4.lsw = argSqrd.lsw;
-      mul_Xsig_Xsig(&argTo4, &argTo4);
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_l,
+				N_COEFF_N - 1);
+		mul_Xsig_Xsig(&accumulator, &argSqrd);
+		negate_Xsig(&accumulator);
 
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_l,
-		      N_COEFF_N-1);
-      mul_Xsig_Xsig(&accumulator, &argSqrd);
-      negate_Xsig(&accumulator);
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_l,
+				N_COEFF_P - 1);
 
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_l,
-		      N_COEFF_P-1);
+		shr_Xsig(&accumulator, 2);	/* Divide by four */
+		accumulator.msw |= 0x80000000;	/* Add 1.0 */
 
-      shr_Xsig(&accumulator, 2);    /* Divide by four */
-      accumulator.msw |= 0x80000000;  /* Add 1.0 */
+		mul64_Xsig(&accumulator, &significand(st0_ptr));
+		mul64_Xsig(&accumulator, &significand(st0_ptr));
+		mul64_Xsig(&accumulator, &significand(st0_ptr));
 
-      mul64_Xsig(&accumulator, &significand(st0_ptr));
-      mul64_Xsig(&accumulator, &significand(st0_ptr));
-      mul64_Xsig(&accumulator, &significand(st0_ptr));
+		/* Divide by four, FPU_REG compatible, etc */
+		exponent = 3 * exponent;
 
-      /* Divide by four, FPU_REG compatible, etc */
-      exponent = 3*exponent;
+		/* The minimum exponent difference is 3 */
+		shr_Xsig(&accumulator, exponent(st0_ptr) - exponent);
 
-      /* The minimum exponent difference is 3 */
-      shr_Xsig(&accumulator, exponent(st0_ptr) - exponent);
+		negate_Xsig(&accumulator);
+		XSIG_LL(accumulator) += significand(st0_ptr);
 
-      negate_Xsig(&accumulator);
-      XSIG_LL(accumulator) += significand(st0_ptr);
+		echange = round_Xsig(&accumulator);
 
-      echange = round_Xsig(&accumulator);
+		setexponentpos(&result, exponent(st0_ptr) + echange);
+	} else {
+		/* The argument is > 0.88309101259 */
+		/* We use sin(st(0)) = cos(pi/2-st(0)) */
 
-      setexponentpos(&result, exponent(st0_ptr) + echange);
-    }
-  else
-    {
-      /* The argument is > 0.88309101259 */
-      /* We use sin(st(0)) = cos(pi/2-st(0)) */
+		fixed_arg = significand(st0_ptr);
 
-      fixed_arg = significand(st0_ptr);
+		if (exponent == 0) {
+			/* The argument is >= 1.0 */
 
-      if ( exponent == 0 )
-	{
-	  /* The argument is >= 1.0 */
+			/* Put the binary point at the left. */
+			fixed_arg <<= 1;
+		}
+		/* pi/2 in hex is: 1.921fb54442d18469 898CC51701B839A2 52049C1 */
+		fixed_arg = 0x921fb54442d18469LL - fixed_arg;
+		/* There is a special case which arises due to rounding, to fix here. */
+		if (fixed_arg == 0xffffffffffffffffLL)
+			fixed_arg = 0;
 
-	  /* Put the binary point at the left. */
-	  fixed_arg <<= 1;
-	}
-      /* pi/2 in hex is: 1.921fb54442d18469 898CC51701B839A2 52049C1 */
-      fixed_arg = 0x921fb54442d18469LL - fixed_arg;
-      /* There is a special case which arises due to rounding, to fix here. */
-      if ( fixed_arg == 0xffffffffffffffffLL )
-	fixed_arg = 0;
+		XSIG_LL(argSqrd) = fixed_arg;
+		argSqrd.lsw = 0;
+		mul64_Xsig(&argSqrd, &fixed_arg);
 
-      XSIG_LL(argSqrd) = fixed_arg; argSqrd.lsw = 0;
-      mul64_Xsig(&argSqrd, &fixed_arg);
+		XSIG_LL(argTo4) = XSIG_LL(argSqrd);
+		argTo4.lsw = argSqrd.lsw;
+		mul_Xsig_Xsig(&argTo4, &argTo4);
 
-      XSIG_LL(argTo4) = XSIG_LL(argSqrd); argTo4.lsw = argSqrd.lsw;
-      mul_Xsig_Xsig(&argTo4, &argTo4);
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_h,
+				N_COEFF_NH - 1);
+		mul_Xsig_Xsig(&accumulator, &argSqrd);
+		negate_Xsig(&accumulator);
 
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_h,
-		      N_COEFF_NH-1);
-      mul_Xsig_Xsig(&accumulator, &argSqrd);
-      negate_Xsig(&accumulator);
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_h,
+				N_COEFF_PH - 1);
+		negate_Xsig(&accumulator);
 
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_h,
-		      N_COEFF_PH-1);
-      negate_Xsig(&accumulator);
+		mul64_Xsig(&accumulator, &fixed_arg);
+		mul64_Xsig(&accumulator, &fixed_arg);
 
-      mul64_Xsig(&accumulator, &fixed_arg);
-      mul64_Xsig(&accumulator, &fixed_arg);
+		shr_Xsig(&accumulator, 3);
+		negate_Xsig(&accumulator);
 
-      shr_Xsig(&accumulator, 3);
-      negate_Xsig(&accumulator);
+		add_Xsig_Xsig(&accumulator, &argSqrd);
 
-      add_Xsig_Xsig(&accumulator, &argSqrd);
+		shr_Xsig(&accumulator, 1);
 
-      shr_Xsig(&accumulator, 1);
+		accumulator.lsw |= 1;	/* A zero accumulator here would cause problems */
+		negate_Xsig(&accumulator);
 
-      accumulator.lsw |= 1;  /* A zero accumulator here would cause problems */
-      negate_Xsig(&accumulator);
+		/* The basic computation is complete. Now fix the answer to
+		   compensate for the error due to the approximation used for
+		   pi/2
+		 */
 
-      /* The basic computation is complete. Now fix the answer to
-	 compensate for the error due to the approximation used for
-	 pi/2
-	 */
+		/* This has an exponent of -65 */
+		fix_up = 0x898cc517;
+		/* The fix-up needs to be improved for larger args */
+		if (argSqrd.msw & 0xffc00000) {
+			/* Get about 32 bit precision in these: */
+			fix_up -= mul_32_32(0x898cc517, argSqrd.msw) / 6;
+		}
+		fix_up = mul_32_32(fix_up, LL_MSW(fixed_arg));
 
-      /* This has an exponent of -65 */
-      fix_up = 0x898cc517;
-      /* The fix-up needs to be improved for larger args */
-      if ( argSqrd.msw & 0xffc00000 )
-	{
-	  /* Get about 32 bit precision in these: */
-	  fix_up -= mul_32_32(0x898cc517, argSqrd.msw) / 6;
-	}
-      fix_up = mul_32_32(fix_up, LL_MSW(fixed_arg));
+		adj = accumulator.lsw;	/* temp save */
+		accumulator.lsw -= fix_up;
+		if (accumulator.lsw > adj)
+			XSIG_LL(accumulator)--;
 
-      adj = accumulator.lsw;    /* temp save */
-      accumulator.lsw -= fix_up;
-      if ( accumulator.lsw > adj )
-	XSIG_LL(accumulator) --;
+		echange = round_Xsig(&accumulator);
 
-      echange = round_Xsig(&accumulator);
-
-      setexponentpos(&result, echange - 1);
-    }
+		setexponentpos(&result, echange - 1);
+	}
 
-  significand(&result) = XSIG_LL(accumulator);
-  setsign(&result, getsign(st0_ptr));
-  FPU_copy_to_reg0(&result, TAG_Valid);
+	significand(&result) = XSIG_LL(accumulator);
+	setsign(&result, getsign(st0_ptr));
+	FPU_copy_to_reg0(&result, TAG_Valid);
 
 #ifdef PARANOID
-  if ( (exponent(&result) >= 0)
-      && (significand(&result) > 0x8000000000000000LL) )
-    {
-      EXCEPTION(EX_INTERNAL|0x150);
-    }
+	if ((exponent(&result) >= 0)
+	    && (significand(&result) > 0x8000000000000000LL)) {
+		EXCEPTION(EX_INTERNAL | 0x150);
+	}
 #endif /* PARANOID */
 
 }
 
-
-
 /*--- poly_cos() ------------------------------------------------------------+
  |                                                                           |
  +---------------------------------------------------------------------------*/
-void	poly_cos(FPU_REG *st0_ptr)
+void poly_cos(FPU_REG *st0_ptr)
 {
-  FPU_REG	      result;
-  long int            exponent, exp2, echange;
-  Xsig                accumulator, argSqrd, fix_up, argTo4;
-  unsigned long long  fixed_arg;
+	FPU_REG result;
+	long int exponent, exp2, echange;
+	Xsig accumulator, argSqrd, fix_up, argTo4;
+	unsigned long long fixed_arg;
 
 #ifdef PARANOID
-  if ( (exponent(st0_ptr) > 0)
-      || ((exponent(st0_ptr) == 0)
-	  && (significand(st0_ptr) > 0xc90fdaa22168c234LL)) )
-    {
-      EXCEPTION(EX_Invalid);
-      FPU_copy_to_reg0(&CONST_QNaN, TAG_Special);
-      return;
-    }
-#endif /* PARANOID */
-
-  exponent = exponent(st0_ptr);
-
-  accumulator.lsw = accumulator.midw = accumulator.msw = 0;
-
-  if ( (exponent < -1) || ((exponent == -1) && (st0_ptr->sigh <= 0xb00d6f54)) )
-    {
-      /* arg is < 0.687705 */
-
-      argSqrd.msw = st0_ptr->sigh; argSqrd.midw = st0_ptr->sigl;
-      argSqrd.lsw = 0;
-      mul64_Xsig(&argSqrd, &significand(st0_ptr));
-
-      if ( exponent < -1 )
-	{
-	  /* shift the argument right by the required places */
-	  shr_Xsig(&argSqrd, 2*(-1-exponent));
-	}
-
-      argTo4.msw = argSqrd.msw; argTo4.midw = argSqrd.midw;
-      argTo4.lsw = argSqrd.lsw;
-      mul_Xsig_Xsig(&argTo4, &argTo4);
-
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_h,
-		      N_COEFF_NH-1);
-      mul_Xsig_Xsig(&accumulator, &argSqrd);
-      negate_Xsig(&accumulator);
-
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_h,
-		      N_COEFF_PH-1);
-      negate_Xsig(&accumulator);
-
-      mul64_Xsig(&accumulator, &significand(st0_ptr));
-      mul64_Xsig(&accumulator, &significand(st0_ptr));
-      shr_Xsig(&accumulator, -2*(1+exponent));
-
-      shr_Xsig(&accumulator, 3);
-      negate_Xsig(&accumulator);
-
-      add_Xsig_Xsig(&accumulator, &argSqrd);
-
-      shr_Xsig(&accumulator, 1);
-
-      /* It doesn't matter if accumulator is all zero here, the
-	 following code will work ok */
-      negate_Xsig(&accumulator);
-
-      if ( accumulator.lsw & 0x80000000 )
-	XSIG_LL(accumulator) ++;
-      if ( accumulator.msw == 0 )
-	{
-	  /* The result is 1.0 */
-	  FPU_copy_to_reg0(&CONST_1, TAG_Valid);
-	  return;
-	}
-      else
-	{
-	  significand(&result) = XSIG_LL(accumulator);
-      
-	  /* will be a valid positive nr with expon = -1 */
-	  setexponentpos(&result, -1);
-	}
-    }
-  else
-    {
-      fixed_arg = significand(st0_ptr);
-
-      if ( exponent == 0 )
-	{
-	  /* The argument is >= 1.0 */
-
-	  /* Put the binary point at the left. */
-	  fixed_arg <<= 1;
-	}
-      /* pi/2 in hex is: 1.921fb54442d18469 898CC51701B839A2 52049C1 */
-      fixed_arg = 0x921fb54442d18469LL - fixed_arg;
-      /* There is a special case which arises due to rounding, to fix here. */
-      if ( fixed_arg == 0xffffffffffffffffLL )
-	fixed_arg = 0;
-
-      exponent = -1;
-      exp2 = -1;
-
-      /* A shift is needed here only for a narrow range of arguments,
-	 i.e. for fixed_arg approx 2^-32, but we pick up more... */
-      if ( !(LL_MSW(fixed_arg) & 0xffff0000) )
-	{
-	  fixed_arg <<= 16;
-	  exponent -= 16;
-	  exp2 -= 16;
+	if ((exponent(st0_ptr) > 0)
+	    || ((exponent(st0_ptr) == 0)
+		&& (significand(st0_ptr) > 0xc90fdaa22168c234LL))) {
+		EXCEPTION(EX_Invalid);
+		FPU_copy_to_reg0(&CONST_QNaN, TAG_Special);
+		return;
 	}
+#endif /* PARANOID */
 
-      XSIG_LL(argSqrd) = fixed_arg; argSqrd.lsw = 0;
-      mul64_Xsig(&argSqrd, &fixed_arg);
-
-      if ( exponent < -1 )
-	{
-	  /* shift the argument right by the required places */
-	  shr_Xsig(&argSqrd, 2*(-1-exponent));
-	}
-
-      argTo4.msw = argSqrd.msw; argTo4.midw = argSqrd.midw;
-      argTo4.lsw = argSqrd.lsw;
-      mul_Xsig_Xsig(&argTo4, &argTo4);
-
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_l,
-		      N_COEFF_N-1);
-      mul_Xsig_Xsig(&accumulator, &argSqrd);
-      negate_Xsig(&accumulator);
-
-      polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_l,
-		      N_COEFF_P-1);
-
-      shr_Xsig(&accumulator, 2);    /* Divide by four */
-      accumulator.msw |= 0x80000000;  /* Add 1.0 */
-
-      mul64_Xsig(&accumulator, &fixed_arg);
-      mul64_Xsig(&accumulator, &fixed_arg);
-      mul64_Xsig(&accumulator, &fixed_arg);
-
-      /* Divide by four, FPU_REG compatible, etc */
-      exponent = 3*exponent;
-
-      /* The minimum exponent difference is 3 */
-      shr_Xsig(&accumulator, exp2 - exponent);
-
-      negate_Xsig(&accumulator);
-      XSIG_LL(accumulator) += fixed_arg;
-
-      /* The basic computation is complete. Now fix the answer to
-	 compensate for the error due to the approximation used for
-	 pi/2
-	 */
-
-      /* This has an exponent of -65 */
-      XSIG_LL(fix_up) = 0x898cc51701b839a2ll;
-      fix_up.lsw = 0;
-
-      /* The fix-up needs to be improved for larger args */
-      if ( argSqrd.msw & 0xffc00000 )
-	{
-	  /* Get about 32 bit precision in these: */
-	  fix_up.msw -= mul_32_32(0x898cc517, argSqrd.msw) / 2;
-	  fix_up.msw += mul_32_32(0x898cc517, argTo4.msw) / 24;
+	exponent = exponent(st0_ptr);
+
+	accumulator.lsw = accumulator.midw = accumulator.msw = 0;
+
+	if ((exponent < -1)
+	    || ((exponent == -1) && (st0_ptr->sigh <= 0xb00d6f54))) {
+		/* arg is < 0.687705 */
+
+		argSqrd.msw = st0_ptr->sigh;
+		argSqrd.midw = st0_ptr->sigl;
+		argSqrd.lsw = 0;
+		mul64_Xsig(&argSqrd, &significand(st0_ptr));
+
+		if (exponent < -1) {
+			/* shift the argument right by the required places */
+			shr_Xsig(&argSqrd, 2 * (-1 - exponent));
+		}
+
+		argTo4.msw = argSqrd.msw;
+		argTo4.midw = argSqrd.midw;
+		argTo4.lsw = argSqrd.lsw;
+		mul_Xsig_Xsig(&argTo4, &argTo4);
+
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_h,
+				N_COEFF_NH - 1);
+		mul_Xsig_Xsig(&accumulator, &argSqrd);
+		negate_Xsig(&accumulator);
+
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_h,
+				N_COEFF_PH - 1);
+		negate_Xsig(&accumulator);
+
+		mul64_Xsig(&accumulator, &significand(st0_ptr));
+		mul64_Xsig(&accumulator, &significand(st0_ptr));
+		shr_Xsig(&accumulator, -2 * (1 + exponent));
+
+		shr_Xsig(&accumulator, 3);
+		negate_Xsig(&accumulator);
+
+		add_Xsig_Xsig(&accumulator, &argSqrd);
+
+		shr_Xsig(&accumulator, 1);
+
+		/* It doesn't matter if accumulator is all zero here, the
+		   following code will work ok */
+		negate_Xsig(&accumulator);
+
+		if (accumulator.lsw & 0x80000000)
+			XSIG_LL(accumulator)++;
+		if (accumulator.msw == 0) {
+			/* The result is 1.0 */
+			FPU_copy_to_reg0(&CONST_1, TAG_Valid);
+			return;
+		} else {
+			significand(&result) = XSIG_LL(accumulator);
+
+			/* will be a valid positive nr with expon = -1 */
+			setexponentpos(&result, -1);
+		}
+	} else {
+		fixed_arg = significand(st0_ptr);
+
+		if (exponent == 0) {
+			/* The argument is >= 1.0 */
+
+			/* Put the binary point at the left. */
+			fixed_arg <<= 1;
+		}
+		/* pi/2 in hex is: 1.921fb54442d18469 898CC51701B839A2 52049C1 */
+		fixed_arg = 0x921fb54442d18469LL - fixed_arg;
+		/* There is a special case which arises due to rounding, to fix here. */
+		if (fixed_arg == 0xffffffffffffffffLL)
+			fixed_arg = 0;
+
+		exponent = -1;
+		exp2 = -1;
+
+		/* A shift is needed here only for a narrow range of arguments,
+		   i.e. for fixed_arg approx 2^-32, but we pick up more... */
+		if (!(LL_MSW(fixed_arg) & 0xffff0000)) {
+			fixed_arg <<= 16;
+			exponent -= 16;
+			exp2 -= 16;
+		}
+
+		XSIG_LL(argSqrd) = fixed_arg;
+		argSqrd.lsw = 0;
+		mul64_Xsig(&argSqrd, &fixed_arg);
+
+		if (exponent < -1) {
+			/* shift the argument right by the required places */
+			shr_Xsig(&argSqrd, 2 * (-1 - exponent));
+		}
+
+		argTo4.msw = argSqrd.msw;
+		argTo4.midw = argSqrd.midw;
+		argTo4.lsw = argSqrd.lsw;
+		mul_Xsig_Xsig(&argTo4, &argTo4);
+
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), neg_terms_l,
+				N_COEFF_N - 1);
+		mul_Xsig_Xsig(&accumulator, &argSqrd);
+		negate_Xsig(&accumulator);
+
+		polynomial_Xsig(&accumulator, &XSIG_LL(argTo4), pos_terms_l,
+				N_COEFF_P - 1);
+
+		shr_Xsig(&accumulator, 2);	/* Divide by four */
+		accumulator.msw |= 0x80000000;	/* Add 1.0 */
+
+		mul64_Xsig(&accumulator, &fixed_arg);
+		mul64_Xsig(&accumulator, &fixed_arg);
+		mul64_Xsig(&accumulator, &fixed_arg);
+
+		/* Divide by four, FPU_REG compatible, etc */
+		exponent = 3 * exponent;
+
+		/* The minimum exponent difference is 3 */
+		shr_Xsig(&accumulator, exp2 - exponent);
+
+		negate_Xsig(&accumulator);
+		XSIG_LL(accumulator) += fixed_arg;
+
+		/* The basic computation is complete. Now fix the answer to
+		   compensate for the error due to the approximation used for
+		   pi/2
+		 */
+
+		/* This has an exponent of -65 */
+		XSIG_LL(fix_up) = 0x898cc51701b839a2ll;
+		fix_up.lsw = 0;
+
+		/* The fix-up needs to be improved for larger args */
+		if (argSqrd.msw & 0xffc00000) {
+			/* Get about 32 bit precision in these: */
+			fix_up.msw -= mul_32_32(0x898cc517, argSqrd.msw) / 2;
+			fix_up.msw += mul_32_32(0x898cc517, argTo4.msw) / 24;
+		}
+
+		exp2 += norm_Xsig(&accumulator);
+		shr_Xsig(&accumulator, 1);	/* Prevent overflow */
+		exp2++;
+		shr_Xsig(&fix_up, 65 + exp2);
+
+		add_Xsig_Xsig(&accumulator, &fix_up);
+
+		echange = round_Xsig(&accumulator);
+
+		setexponentpos(&result, exp2 + echange);
+		significand(&result) = XSIG_LL(accumulator);
 	}
 
-      exp2 += norm_Xsig(&accumulator);
-      shr_Xsig(&accumulator, 1); /* Prevent overflow */
-      exp2++;
-      shr_Xsig(&fix_up, 65 + exp2);
-
-      add_Xsig_Xsig(&accumulator, &fix_up);
-
-      echange = round_Xsig(&accumulator);
-
-      setexponentpos(&result, exp2 + echange);
-      significand(&result) = XSIG_LL(accumulator);
-    }
-
-  FPU_copy_to_reg0(&result, TAG_Valid);
+	FPU_copy_to_reg0(&result, TAG_Valid);
 
 #ifdef PARANOID
-  if ( (exponent(&result) >= 0)
-      && (significand(&result) > 0x8000000000000000LL) )
-    {
-      EXCEPTION(EX_INTERNAL|0x151);
-    }
+	if ((exponent(&result) >= 0)
+	    && (significand(&result) > 0x8000000000000000LL)) {
+		EXCEPTION(EX_INTERNAL | 0x151);
+	}
 #endif /* PARANOID */
 
 }