From 1eba0867060b1643d71d4185fbc23995248092bf Mon Sep 17 00:00:00 2001
From: Jakub Jelinek <jakub@redhat.com>
Date: Wed, 19 Jul 2017 15:12:58 +0200
Subject: [PATCH] re PR libquadmath/65757 (gfortran gives incorrect result for
 anint with real*16 argument)

	PR libquadmath/65757
	* quadmath-imp.h (math_opt_barrier, math_force_eval,
	math_narrow_eval, math_check_force_underflow,
	math_check_force_underflow_nonneg): Define.
	* math/ceilq.c: Backport changes from upstream glibc
	between 2012-11-01 and 2017-07-13.
	* math/remquoq.c: Likewise.
	* math/expq.c: Likewise.
	* math/llroundq.c: Likewise.
	* math/logq.c: Likewise.
	* math/atanq.c: Likewise.
	* math/nearbyintq.c: Likewise.
	* math/scalblnq.c: Likewise.
	* math/finiteq.c: Likewise.
	* math/atanhq.c: Likewise.
	* math/expm1q.c: Likewise.
	* math/sinhq.c: Likewise.
	* math/log10q.c: Likewise.
	* math/rintq.c: Likewise.
	* math/roundq.c: Likewise.
	* math/fmaq.c: Likewise.
	* math/erfq.c: Likewise.
	* math/log2q.c: Likewise.
	* math/lroundq.c: Likewise.
	* math/j1q.c: Likewise.
	* math/scalbnq.c: Likewise.
	* math/truncq.c: Likewise.
	* math/frexpq.c: Likewise.
	* math/sincosq.c: Likewise.
	* math/tanhq.c: Likewise.
	* math/asinq.c: Likewise.
	* math/coshq.c: Likewise.
	* math/j0q.c: Likewise.
	* math/asinhq.c: Likewise.
	* math/floorq.c: Likewise.
	* math/sinq_kernel.c: Likewise.
	* math/powq.c: Likewise.
	* math/hypotq.c: Likewise.
	* math/sincos_table.c: Likewise.
	* math/rem_pio2q.c: Likewise.
	* math/nextafterq.c: Likewise.
	* math/log1pq.c: Likewise.
	* math/sincosq_kernel.c: Likewise.
	* math/tanq.c: Likewise.
	* math/acosq.c: Likewise.
	* math/lrintq.c: Likewise.
	* math/llrintq.c: Likewise.

From-SVN: r250343
---
 libquadmath/ChangeLog             |  50 +++++
 libquadmath/math/acosq.c          |   2 +-
 libquadmath/math/asinhq.c         |   1 +
 libquadmath/math/asinq.c          |   6 +-
 libquadmath/math/atanhq.c         |   6 +-
 libquadmath/math/atanq.c          |   5 +-
 libquadmath/math/ceilq.c          |  33 ++-
 libquadmath/math/coshq.c          |   4 +-
 libquadmath/math/erfq.c           |  33 ++-
 libquadmath/math/expm1q.c         |  33 ++-
 libquadmath/math/expq.c           |  14 +-
 libquadmath/math/finiteq.c        |   2 +-
 libquadmath/math/floorq.c         |  35 ++--
 libquadmath/math/fmaq.c           |  55 ++---
 libquadmath/math/frexpq.c         |   2 +-
 libquadmath/math/hypotq.c         |  15 +-
 libquadmath/math/j0q.c            |  87 ++++----
 libquadmath/math/j1q.c            | 103 ++++++----
 libquadmath/math/llrintq.c        |  43 +++-
 libquadmath/math/llroundq.c       |  29 ++-
 libquadmath/math/log10q.c         |   5 +-
 libquadmath/math/log1pq.c         |  14 +-
 libquadmath/math/log2q.c          |   5 +-
 libquadmath/math/logq.c           |   7 +-
 libquadmath/math/lrintq.c         |  84 ++++++--
 libquadmath/math/lroundq.c        |  58 ++++--
 libquadmath/math/nearbyintq.c     |   4 +-
 libquadmath/math/nextafterq.c     |  11 +-
 libquadmath/math/powq.c           |  32 +--
 libquadmath/math/rem_pio2q.c      | 330 +++++++++++++++---------------
 libquadmath/math/remquoq.c        |  11 +-
 libquadmath/math/rintq.c          |   2 +-
 libquadmath/math/roundq.c         |  40 ++--
 libquadmath/math/scalblnq.c       |   2 +-
 libquadmath/math/scalbnq.c        |   2 +-
 libquadmath/math/sincos_table.c   |   4 +-
 libquadmath/math/sincosq.c        |   5 +-
 libquadmath/math/sincosq_kernel.c |  17 +-
 libquadmath/math/sinhq.c          |   7 +-
 libquadmath/math/sinq_kernel.c    |   5 +-
 libquadmath/math/tanhq.c          |   5 +-
 libquadmath/math/tanq.c           |  13 +-
 libquadmath/math/truncq.c         |   4 +-
 libquadmath/quadmath-imp.h        |  41 ++++
 44 files changed, 789 insertions(+), 477 deletions(-)

diff --git a/libquadmath/ChangeLog b/libquadmath/ChangeLog
index 4bd166b0ff9..73d0e3077c5 100644
--- a/libquadmath/ChangeLog
+++ b/libquadmath/ChangeLog
@@ -1,3 +1,53 @@
+2017-07-19  Jakub Jelinek  <jakub@redhat.com>
+
+	PR libquadmath/65757
+	* quadmath-imp.h (math_opt_barrier, math_force_eval,
+	math_narrow_eval, math_check_force_underflow,
+	math_check_force_underflow_nonneg): Define.
+	* math/ceilq.c: Backport changes from upstream glibc
+	between 2012-11-01 and 2017-07-13.
+	* math/remquoq.c: Likewise.
+	* math/expq.c: Likewise.
+	* math/llroundq.c: Likewise.
+	* math/logq.c: Likewise.
+	* math/atanq.c: Likewise.
+	* math/nearbyintq.c: Likewise.
+	* math/scalblnq.c: Likewise.
+	* math/finiteq.c: Likewise.
+	* math/atanhq.c: Likewise.
+	* math/expm1q.c: Likewise.
+	* math/sinhq.c: Likewise.
+	* math/log10q.c: Likewise.
+	* math/rintq.c: Likewise.
+	* math/roundq.c: Likewise.
+	* math/fmaq.c: Likewise.
+	* math/erfq.c: Likewise.
+	* math/log2q.c: Likewise.
+	* math/lroundq.c: Likewise.
+	* math/j1q.c: Likewise.
+	* math/scalbnq.c: Likewise.
+	* math/truncq.c: Likewise.
+	* math/frexpq.c: Likewise.
+	* math/sincosq.c: Likewise.
+	* math/tanhq.c: Likewise.
+	* math/asinq.c: Likewise.
+	* math/coshq.c: Likewise.
+	* math/j0q.c: Likewise.
+	* math/asinhq.c: Likewise.
+	* math/floorq.c: Likewise.
+	* math/sinq_kernel.c: Likewise.
+	* math/powq.c: Likewise.
+	* math/hypotq.c: Likewise.
+	* math/sincos_table.c: Likewise.
+	* math/rem_pio2q.c: Likewise.
+	* math/nextafterq.c: Likewise.
+	* math/log1pq.c: Likewise.
+	* math/sincosq_kernel.c: Likewise.
+	* math/tanq.c: Likewise.
+	* math/acosq.c: Likewise.
+	* math/lrintq.c: Likewise.
+	* math/llrintq.c: Likewise.
+
 2017-02-09  Gerald Pfeifer  <gerald@pfeifer.com>
 
 	* configure.ac (ACX_BUGURL): Update.
diff --git a/libquadmath/math/acosq.c b/libquadmath/math/acosq.c
index 7ef79474651..7f2ed2725d1 100644
--- a/libquadmath/math/acosq.c
+++ b/libquadmath/math/acosq.c
@@ -172,7 +172,7 @@ acosq (__float128 x)
     }
   else if (ix < 0x3ffe0000)	/* |x| < 0.5 */
     {
-      if (ix < 0x3fc60000)	/* |x| < 2**-57 */
+      if (ix < 0x3f8e0000)	/* |x| < 2**-113 */
 	return pio2_hi + pio2_lo;
       if (ix < 0x3ffde000)	/* |x| < .4375 */
 	{
diff --git a/libquadmath/math/asinhq.c b/libquadmath/math/asinhq.c
index 9be0aa1f053..f3644939b3a 100644
--- a/libquadmath/math/asinhq.c
+++ b/libquadmath/math/asinhq.c
@@ -46,6 +46,7 @@ asinhq (__float128 x)
     return x + x;		/* x is inf or NaN */
   if (ix < 0x3fc70000)
     {				/* |x| < 2^ -56 */
+      math_check_force_underflow (x);
       if (huge + x > one)
 	return x;		/* return x inexact except 0 */
     }
diff --git a/libquadmath/math/asinq.c b/libquadmath/math/asinq.c
index 7bd4d768c97..5dff2817694 100644
--- a/libquadmath/math/asinq.c
+++ b/libquadmath/math/asinq.c
@@ -151,8 +151,10 @@ asinq (__float128 x)
     {
       if (ix < 0x3fc60000) /* |x| < 2**-57 */
 	{
-	  if (huge + x > one)
-	    return x;		/* return x with inexact if x!=0 */
+	  math_check_force_underflow (x);
+	  __float128 force_inexact = huge + x;
+	  math_force_eval (force_inexact);
+	  return x;		/* return x with inexact if x!=0 */
 	}
       else
 	{
diff --git a/libquadmath/math/atanhq.c b/libquadmath/math/atanhq.c
index b34036715d7..652138d495c 100644
--- a/libquadmath/math/atanhq.c
+++ b/libquadmath/math/atanhq.c
@@ -55,7 +55,11 @@ atanhq (__float128 x)
 	    else
 	      return (x-x)/(x-x);
 	  }
-	if(ix<0x3fc60000 && (huge+x)>zero) return x;	/* x < 2^-57 */
+	if(ix<0x3fc60000 && (huge+x)>zero)	/* x < 2^-57 */
+	  {
+	    math_check_force_underflow (x);
+	    return x;
+	  }
 
 	if(ix<0x3ffe0000) {		/* x < 0.5 */
 	    t = u.value+u.value;
diff --git a/libquadmath/math/atanq.c b/libquadmath/math/atanq.c
index 8eccdc3317d..01fd9d69e57 100644
--- a/libquadmath/math/atanq.c
+++ b/libquadmath/math/atanq.c
@@ -42,7 +42,7 @@
  *
  */
 
-/* Copyright 2001 by Stephen L. Moshier <moshier@na-net.ornl.gov> 
+/* Copyright 2001 by Stephen L. Moshier <moshier@na-net.ornl.gov>
 
     This library is free software; you can redistribute it and/or
     modify it under the terms of the GNU Lesser General Public
@@ -167,7 +167,7 @@ static const __float128
   q4 = 2.173623741810414221251136181221172551416E1Q;
   /* q5 = 1.000000000000000000000000000000000000000E0 */
 
-static const long double huge = 1.0e4930Q;
+static const __float128 huge = 1.0e4930Q;
 
 __float128
 atanq (__float128 x)
@@ -200,6 +200,7 @@ atanq (__float128 x)
 
   if (k <= 0x3fc50000) /* |x| < 2**-58 */
     {
+      math_check_force_underflow (x);
       /* Raise inexact. */
       if (huge + x > 0.0)
 	return x;
diff --git a/libquadmath/math/ceilq.c b/libquadmath/math/ceilq.c
index 0d9bb8b87f3..1adc1e1b9f0 100644
--- a/libquadmath/math/ceilq.c
+++ b/libquadmath/math/ceilq.c
@@ -15,8 +15,6 @@
 
 #include "quadmath-imp.h"
 
-static const __float128 huge = 1.0e4930Q;
-
 __float128
 ceilq (__float128 x)
 {
@@ -25,18 +23,15 @@ ceilq (__float128 x)
 	GET_FLT128_WORDS64(i0,i1,x);
 	j0 = ((i0>>48)&0x7fff)-0x3fff;
 	if(j0<48) {
-	    if(j0<0) { 	/* raise inexact if x != 0 */
-		if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */
-		    if(i0<0) {i0=0x8000000000000000ULL;i1=0;}
-		    else if((i0|i1)!=0) { i0=0x3fff000000000000ULL;i1=0;}
-		}
+	    if(j0<0) {
+		/* return 0*sign(x) if |x|<1 */
+		if(i0<0) {i0=0x8000000000000000ULL;i1=0;}
+		else if((i0|i1)!=0) { i0=0x3fff000000000000ULL;i1=0;}
 	    } else {
 		i = (0x0000ffffffffffffULL)>>j0;
 		if(((i0&i)|i1)==0) return x; /* x is integral */
-		if(huge+x>0.0) {	/* raise inexact flag */
-		    if(i0>0) i0 += (0x0001000000000000LL)>>j0;
-		    i0 &= (~i); i1=0;
-		}
+		if(i0>0) i0 += (0x0001000000000000LL)>>j0;
+		i0 &= (~i); i1=0;
 	    }
 	} else if (j0>111) {
 	    if(j0==0x4000) return x+x;	/* inf or NaN */
@@ -44,17 +39,15 @@ ceilq (__float128 x)
 	} else {
 	    i = -1ULL>>(j0-48);
 	    if((i1&i)==0) return x;	/* x is integral */
-	    if(huge+x>0.0) { 		/* raise inexact flag */
-		if(i0>0) {
-		    if(j0==48) i0+=1;
-		    else {
-			j = i1+(1LL<<(112-j0));
-			if(j<i1) i0 +=1 ; 	/* got a carry */
-			i1=j;
-		    }
+	    if(i0>0) {
+		if(j0==48) i0+=1;
+		else {
+		    j = i1+(1LL<<(112-j0));
+		    if(j<i1) i0 +=1 ; 	/* got a carry */
+		    i1=j;
 		}
-		i1 &= (~i);
 	    }
+	    i1 &= (~i);
 	}
 	SET_FLT128_WORDS64(x,i0,i1);
 	return x;
diff --git a/libquadmath/math/coshq.c b/libquadmath/math/coshq.c
index 77f4b98338b..61397509527 100644
--- a/libquadmath/math/coshq.c
+++ b/libquadmath/math/coshq.c
@@ -76,10 +76,10 @@ coshq (__float128 x)
   /* |x| in [0,0.5*ln2], return 1+expm1l(|x|)^2/(2*expq(|x|)) */
   if (ex < 0x3ffd62e4) /* 0.3465728759765625 */
     {
+      if (ex < 0x3fb80000) /* |x| < 2^-116 */
+	return one;		/* cosh(tiny) = 1 */
       t = expm1q (u.value);
       w = one + t;
-      if (ex < 0x3fb80000) /* |x| < 2^-116 */
-	return w;		/* cosh(tiny) = 1 */
 
       return one + (t * t) / (w + w);
     }
diff --git a/libquadmath/math/erfq.c b/libquadmath/math/erfq.c
index 8d383e9ca70..45a8c014e5c 100644
--- a/libquadmath/math/erfq.c
+++ b/libquadmath/math/erfq.c
@@ -11,9 +11,9 @@
 
 /* Modifications and expansions for 128-bit long double are
    Copyright (C) 2001 Stephen L. Moshier <moshier@na-net.ornl.gov>
-   and are incorporated herein by permission of the author.  The author 
+   and are incorporated herein by permission of the author.  The author
    reserves the right to distribute this material elsewhere under different
-   copying permissions.  These modifications are distributed here under 
+   copying permissions.  These modifications are distributed here under
    the following terms:
 
     This library is free software; you can redistribute it and/or
@@ -96,6 +96,7 @@
  *		erfc/erf(NaN) is NaN
  */
 
+#include <errno.h>
 #include "quadmath-imp.h"
 
 
@@ -142,13 +143,10 @@ deval (__float128 x, const __float128 *p, int n)
 
 static const __float128
 tiny = 1e-4931Q,
-  half = 0.5Q,
   one = 1.0Q,
   two = 2.0Q,
   /* 2/sqrt(pi) - 1 */
-  efx = 1.2837916709551257389615890312154517168810E-1Q,
-  /* 8 * (2/sqrt(pi) - 1) */
-  efx8 = 1.0270333367641005911692712249723613735048E0Q;
+  efx = 1.2837916709551257389615890312154517168810E-1Q;
 
 
 /* erf(x)  = x  + x R(x^2)
@@ -773,6 +771,8 @@ erfq (__float128 x)
 
   if (ix >= 0x3fff0000) /* |x| >= 1.0 */
     {
+      if (ix >= 0x40030000 && sign > 0)
+	return one; /* x >= 16, avoid spurious underflow from erfc.  */
       y = erfcq (x);
       return (one - y);
       /*    return (one - erfcq (x)); */
@@ -785,7 +785,12 @@ erfq (__float128 x)
       if (ix < 0x3fc60000) /* |x|<2**-57 */
 	{
 	  if (ix < 0x00080000)
-	    return 0.125 * (8.0 * x + efx8 * x);	/*avoid underflow */
+	    {
+	      /* Avoid spurious underflow.  */
+	      __float128 ret =  0.0625 * (16.0 * x + (16.0 * efx) * x);
+	      math_check_force_underflow (ret);
+	      return ret;
+	    }
 	  return x + efx * x;
 	}
       y = a + a * neval (z, TN1, NTN1) / deval (z, TD1, NTD1);
@@ -867,7 +872,7 @@ erfcq (__float128 x)
 	  y = C19b + z * neval (z, RNr19, NRNr19) / deval (z, RDr19, NRDr19);
 	  y += C19a;
 	  break;
-	case 9:
+	default: /* i == 9.  */
 	  z = x - 1.125Q;
 	  y = C20b + z * neval (z, RNr20, NRNr20) / deval (z, RDr20, NRDr20);
 	  y += C20a;
@@ -921,14 +926,22 @@ erfcq (__float128 x)
       z = u.value;
       r = expq (-z * z - 0.5625) * expq ((z - x) * (z + x) + p);
       if ((sign & 0x80000000) == 0)
-	return r / x;
+	{
+	  __float128 ret = r / x;
+	  if (ret == 0)
+	    errno = ERANGE;
+	  return ret;
+	}
       else
 	return two - r / x;
     }
   else
     {
       if ((sign & 0x80000000) == 0)
-	return tiny * tiny;
+	{
+	  errno = ERANGE;
+	  return tiny * tiny;
+	}
       else
 	return two - tiny;
     }
diff --git a/libquadmath/math/expm1q.c b/libquadmath/math/expm1q.c
index 8cfdd8eec94..9060d480858 100644
--- a/libquadmath/math/expm1q.c
+++ b/libquadmath/math/expm1q.c
@@ -35,7 +35,7 @@
  *
  */
 
-/* Copyright 2001 by Stephen L. Moshier 
+/* Copyright 2001 by Stephen L. Moshier
 
     This library is free software; you can redistribute it and/or
     modify it under the terms of the GNU Lesser General Public
@@ -82,8 +82,6 @@ static const __float128
 
   C1 = 6.93145751953125E-1Q,
   C2 = 1.428606820309417232121458176568075500134E-6Q,
-/* ln (2^16384 * (1 - 2^-113)) */
-  maxlog = 1.1356523406294143949491931077970764891253E4Q,
 /* ln 2^-114 */
   minarg = -7.9018778583833765273564461846232128760607E1Q;
 
@@ -108,33 +106,30 @@ expm1q (__float128 x)
     }
   if (ix >= 0x7fff0000)
     {
-      /* Infinity. */
+      /* Infinity (which must be negative infinity). */
       if (((ix & 0xffff) | u.words32.w1 | u.words32.w2 | u.words32.w3) == 0)
-	{
-	  if (sign)
-	    return -1.0Q;
-	  else
-	    return x;
-	}
-      /* NaN. No invalid exception. */
-      return x;
+	return -1.0Q;
+      /* NaN.  Invalid exception if signaling.  */
+      return x + x;
     }
 
   /* expm1(+- 0) = +- 0.  */
   if ((ix == 0) && (u.words32.w1 | u.words32.w2 | u.words32.w3) == 0)
     return x;
 
-  /* Overflow.  */
-  if (x > maxlog)
-    {
-      errno = ERANGE;
-      return (HUGE_VALQ * HUGE_VALQ);
-    }
-
   /* Minimum value.  */
   if (x < minarg)
     return (4.0/HUGE_VALQ - 1.0Q);
 
+  /* Avoid internal underflow when result does not underflow, while
+     ensuring underflow (without returning a zero of the wrong sign)
+     when the result does underflow.  */
+  if (fabsq (x) < 0x1p-113Q)
+    {
+      math_check_force_underflow (x);
+      return x;
+    }
+
   /* Express x = ln 2 (k + remainder), remainder not exceeding 1/2. */
   xx = C1 + C2;			/* ln 2. */
   px = floorq (0.5 + x / xx);
diff --git a/libquadmath/math/expq.c b/libquadmath/math/expq.c
index 70c638d7a08..5df6cd8e192 100644
--- a/libquadmath/math/expq.c
+++ b/libquadmath/math/expq.c
@@ -1,5 +1,5 @@
 /* Quad-precision floating point e^x.
-   Copyright (C) 1999 Free Software Foundation, Inc.
+   Copyright (C) 1999-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Jakub Jelinek <jj@ultra.linux.cz>
    Partly based on double-precision code
@@ -1075,7 +1075,7 @@ static const __float128 C[] = {
 #define TWO15 C[11]
  32768.0Q,
 
-/* Chebyshev polynom coeficients for (exp(x)-1)/x */
+/* Chebyshev polynom coefficients for (exp(x)-1)/x */
 #define P1 C[12]
 #define P2 C[13]
 #define P3 C[14]
@@ -1142,7 +1142,7 @@ expq (__float128 x)
 		* __expq_table[T_EXPL_RES2 + tval2];
       n_i = (int)n;
       /* 'unsafe' is 1 iff n_1 != 0.  */
-      unsafe = abs(n_i) >= -FLT128_MIN_EXP - 1;
+      unsafe = abs(n_i) >= 15000;
       ex2_u.ieee.exponent += n_i >> unsafe;
 
       /* Compute scale = 2^n_1.  */
@@ -1179,7 +1179,7 @@ expq (__float128 x)
 	  ex3_u.d = (result - ex2_u.d) - x22 * ex2_u.d;
 	  ex2_u.d = result;
 	  ex3_u.ieee.exponent += LDBL_MANT_DIG + 15 + IEEE854_LONG_DOUBLE_BIAS
-	  			 - ex2_u.ieee.exponent;
+				 - ex2_u.ieee.exponent;
 	  n_i = abs (ex3_u.d);
 	  n_i = (n_i + 1) / 2;
 #ifdef USE_FENV_H
@@ -1196,7 +1196,11 @@ expq (__float128 x)
       if (!unsafe)
 	return result;
       else
-	return result * scale_u.value;
+	{
+	  result *= scale_u.value;
+	  math_check_force_underflow_nonneg (result);
+	  return result;
+	}
     }
   /* Exceptional cases:  */
   else if (__builtin_isless (x, himark))
diff --git a/libquadmath/math/finiteq.c b/libquadmath/math/finiteq.c
index c5326d2194b..e6703fb2261 100644
--- a/libquadmath/math/finiteq.c
+++ b/libquadmath/math/finiteq.c
@@ -25,6 +25,6 @@ finiteq (const __float128 x)
 {
   int64_t hx;
   GET_FLT128_MSW64(hx,x);
-  return (int)((uint64_t)((hx&0x7fffffffffffffffLL)
+  return (int)((uint64_t)((hx&0x7fff000000000000LL)
 			   -0x7fff000000000000LL)>>63);
 }
diff --git a/libquadmath/math/floorq.c b/libquadmath/math/floorq.c
index 0d74f94e27d..41b993fa7a0 100644
--- a/libquadmath/math/floorq.c
+++ b/libquadmath/math/floorq.c
@@ -15,8 +15,6 @@
 
 #include "quadmath-imp.h"
 
-static const __float128 huge = 1.0e4930Q;
-
 __float128
 floorq (__float128 x)
 {
@@ -25,19 +23,16 @@ floorq (__float128 x)
 	GET_FLT128_WORDS64(i0,i1,x);
 	j0 = ((i0>>48)&0x7fff)-0x3fff;
 	if(j0<48) {
-	    if(j0<0) { 	/* raise inexact if x != 0 */
-		if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */
-		    if(i0>=0) {i0=i1=0;}
-		    else if(((i0&0x7fffffffffffffffLL)|i1)!=0)
-			{ i0=0xbfff000000000000ULL;i1=0;}
-		}
+	    if(j0<0) {
+		/* return 0*sign(x) if |x|<1 */
+		if(i0>=0) {i0=i1=0;}
+		else if(((i0&0x7fffffffffffffffLL)|i1)!=0)
+		    { i0=0xbfff000000000000ULL;i1=0;}
 	    } else {
 		i = (0x0000ffffffffffffULL)>>j0;
 		if(((i0&i)|i1)==0) return x; /* x is integral */
-		if(huge+x>0.0) {	/* raise inexact flag */
-		    if(i0<0) i0 += (0x0001000000000000LL)>>j0;
-		    i0 &= (~i); i1=0;
-		}
+		if(i0<0) i0 += (0x0001000000000000LL)>>j0;
+		i0 &= (~i); i1=0;
 	    }
 	} else if (j0>111) {
 	    if(j0==0x4000) return x+x;	/* inf or NaN */
@@ -45,17 +40,15 @@ floorq (__float128 x)
 	} else {
 	    i = -1ULL>>(j0-48);
 	    if((i1&i)==0) return x;	/* x is integral */
-	    if(huge+x>0.0) { 		/* raise inexact flag */
-		if(i0<0) {
-		    if(j0==48) i0+=1;
-		    else {
-			j = i1+(1LL<<(112-j0));
-			if(j<i1) i0 +=1 ; 	/* got a carry */
-			i1=j;
-		    }
+	    if(i0<0) {
+		if(j0==48) i0+=1;
+		else {
+		    j = i1+(1LL<<(112-j0));
+		    if(j<i1) i0 +=1 ; 	/* got a carry */
+		    i1=j;
 		}
-		i1 &= (~i);
 	    }
+	    i1 &= (~i);
 	}
 	SET_FLT128_WORDS64(x,i0,i1);
 	return x;
diff --git a/libquadmath/math/fmaq.c b/libquadmath/math/fmaq.c
index d3c5fb3901a..68a63cf8fd6 100644
--- a/libquadmath/math/fmaq.c
+++ b/libquadmath/math/fmaq.c
@@ -1,5 +1,5 @@
 /* Compute x * y + z as ternary operation.
-   Copyright (C) 2010-2012 Free Software Foundation, Inc.
+   Copyright (C) 2010-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Jakub Jelinek <jakub@redhat.com>, 2010.
 
@@ -97,8 +97,8 @@ fmaq (__float128 x, __float128 y, __float128 z)
 		     && w.ieee.mant_low == 0
 		     && w.ieee.mant_high == 0)))
 	    {
-	      volatile __float128 force_underflow = x * y;
-	      (void) force_underflow;
+	      __float128 force_underflow = x * y;
+	      math_force_eval (force_underflow);
 	    }
 	  return v.value * 0x1p-114Q;
 	}
@@ -161,15 +161,15 @@ fmaq (__float128 x, __float128 y, __float128 z)
 		  <= IEEE854_FLOAT128_BIAS + FLT128_MANT_DIG) */
 	{
 	  if (u.ieee.exponent > v.ieee.exponent)
-	    u.ieee.exponent += 2 * FLT128_MANT_DIG;
+	    u.ieee.exponent += 2 * FLT128_MANT_DIG + 2;
 	  else
-	    v.ieee.exponent += 2 * FLT128_MANT_DIG;
-	  if (w.ieee.exponent <= 4 * FLT128_MANT_DIG + 4)
+	    v.ieee.exponent += 2 * FLT128_MANT_DIG + 2;
+	  if (w.ieee.exponent <= 4 * FLT128_MANT_DIG + 6)
 	    {
 	      if (w.ieee.exponent)
-		w.ieee.exponent += 2 * FLT128_MANT_DIG;
+		w.ieee.exponent += 2 * FLT128_MANT_DIG + 2;
 	      else
-		w.value *= 0x1p226Q;
+		w.value *= 0x1p228Q;
 	      adjust = -1;
 	    }
 	  /* Otherwise x * y should just affect inexact
@@ -182,7 +182,10 @@ fmaq (__float128 x, __float128 y, __float128 z)
 
   /* Ensure correct sign of exact 0 + 0.  */
   if (__builtin_expect ((x == 0 || y == 0) && z == 0, 0))
-    return x * y + z;
+    {
+      x = math_opt_barrier (x);
+      return x * y + z;
+    }
 
 #ifdef USE_FENV_H
   fenv_t env;
@@ -208,24 +211,24 @@ fmaq (__float128 x, __float128 y, __float128 z)
   t1 = m1 - t1;
   t2 = z - t2;
   __float128 a2 = t1 + t2;
+  /* Ensure the arithmetic is not scheduled after feclearexcept call.  */
+  math_force_eval (m2);
+  math_force_eval (a2);
 #ifdef USE_FENV_H
   feclearexcept (FE_INEXACT);
 #endif
 
-  /* If the result is an exact zero, ensure it has the correct
-     sign.  */
+  /* If the result is an exact zero, ensure it has the correct sign.  */
   if (a1 == 0 && m2 == 0)
     {
 #ifdef USE_FENV_H
       feupdateenv (&env);
 #endif
-      /* Ensure that round-to-nearest value of z + m1 is not
-	 reused.  */
-      asm volatile ("" : "=m" (z) : "m" (z));
+      /* Ensure that round-to-nearest value of z + m1 is not reused.  */
+      z = math_opt_barrier (z);
       return z + m1;
     }
 
-
 #ifdef USE_FENV_H
   fesetround (FE_TOWARDZERO);
 #endif
@@ -273,19 +276,19 @@ fmaq (__float128 x, __float128 y, __float128 z)
       /* If a1 + u.value is exact, the only rounding happens during
 	 scaling down.  */
       if (j == 0)
-	return v.value * 0x1p-226Q;
+	return v.value * 0x1p-228Q;
       /* If result rounded to zero is not subnormal, no double
 	 rounding will occur.  */
-      if (v.ieee.exponent > 226)
-	return (a1 + u.value) * 0x1p-226Q;
-      /* If v.value * 0x1p-226Q with round to zero is a subnormal above
-	 or equal to FLT128_MIN / 2, then v.value * 0x1p-226Q shifts mantissa
+      if (v.ieee.exponent > 228)
+	return (a1 + u.value) * 0x1p-228Q;
+      /* If v.value * 0x1p-228Q with round to zero is a subnormal above
+	 or equal to FLT128_MIN / 2, then v.value * 0x1p-228Q shifts mantissa
 	 down just by 1 bit, which means v.ieee.mant_low |= j would
 	 change the round bit, not sticky or guard bit.
-	 v.value * 0x1p-226Q never normalizes by shifting up,
+	 v.value * 0x1p-228Q never normalizes by shifting up,
 	 so round bit plus sticky bit should be already enough
 	 for proper rounding.  */
-      if (v.ieee.exponent == 226)
+      if (v.ieee.exponent == 228)
 	{
 	  /* If the exponent would be in the normal range when
 	     rounding to normal precision with unbounded exponent
@@ -295,8 +298,8 @@ fmaq (__float128 x, __float128 y, __float128 z)
 	  if (TININESS_AFTER_ROUNDING)
 	    {
 	      w.value = a1 + u.value;
-	      if (w.ieee.exponent == 227)
-		return w.value * 0x1p-226Q;
+	      if (w.ieee.exponent == 229)
+		return w.value * 0x1p-228Q;
 	    }
 	  /* v.ieee.mant_low & 2 is LSB bit of the result before rounding,
 	     v.ieee.mant_low & 1 is the round bit and j is our sticky
@@ -305,11 +308,11 @@ fmaq (__float128 x, __float128 y, __float128 z)
 	  w.ieee.mant_low = ((v.ieee.mant_low & 3) << 1) | j;
 	  w.ieee.negative = v.ieee.negative;
 	  v.ieee.mant_low &= ~3U;
-	  v.value *= 0x1p-226Q;
+	  v.value *= 0x1p-228Q;
 	  w.value *= 0x1p-2Q;
 	  return v.value + w.value;
 	}
       v.ieee.mant_low |= j;
-      return v.value * 0x1p-226Q;
+      return v.value * 0x1p-228Q;
     }
 }
diff --git a/libquadmath/math/frexpq.c b/libquadmath/math/frexpq.c
index b0b305af574..2bd77829bf2 100644
--- a/libquadmath/math/frexpq.c
+++ b/libquadmath/math/frexpq.c
@@ -35,7 +35,7 @@ frexpq (__float128 x, int *eptr)
   GET_FLT128_WORDS64(hx,lx,x);
   ix = 0x7fffffffffffffffULL&hx;
   *eptr = 0;
-  if(ix>=0x7fff000000000000ULL||((ix|lx)==0)) return x;	/* 0,inf,nan */
+  if(ix>=0x7fff000000000000ULL||((ix|lx)==0)) return x + x;/* 0,inf,nan */
   if (ix<0x0001000000000000ULL) {		/* subnormal */
     x *= two114;
     GET_FLT128_MSW64(hx,x);
diff --git a/libquadmath/math/hypotq.c b/libquadmath/math/hypotq.c
index 2df317f3681..057901073dc 100644
--- a/libquadmath/math/hypotq.c
+++ b/libquadmath/math/hypotq.c
@@ -89,6 +89,17 @@ hypotq (__float128 x, __float128 y)
   	b *= t1;
   	a *= t1;
   	k -= 16382;
+	GET_FLT128_MSW64 (ha, a);
+	GET_FLT128_MSW64 (hb, b);
+	if (hb > ha)
+	  {
+	    t1 = a;
+	    a = b;
+	    b = t1;
+	    j = ha;
+	    ha = hb;
+	    hb = j;
+	  }
       } else {		/* scale a and b by 2^9600 */
           ha += 0x2580000000000000LL; 	/* a *= 2^9600 */
   	hb += 0x2580000000000000LL;	/* b *= 2^9600 */
@@ -119,6 +130,8 @@ hypotq (__float128 x, __float128 y)
       t1 = 1.0Q;
       GET_FLT128_MSW64(high,t1);
       SET_FLT128_MSW64(t1,high+(k<<48));
-      return t1*w;
+      w *= t1;
+      math_check_force_underflow_nonneg (w);
+      return w;
   } else return w;
 }
diff --git a/libquadmath/math/j0q.c b/libquadmath/math/j0q.c
index 8c6f811125e..c6e482b1c51 100644
--- a/libquadmath/math/j0q.c
+++ b/libquadmath/math/j0q.c
@@ -681,7 +681,7 @@ j0q (__float128 x)
   if (! finiteq (x))
     {
       if (x != x)
-	return x;
+	return x + x;
       else
 	return 0.0Q;
     }
@@ -691,6 +691,8 @@ j0q (__float128 x)
   xx = fabsq (x);
   if (xx <= 2.0Q)
     {
+      if (xx < 0x1p-57Q)
+	return 1.0Q;
       /* 0 <= x <= 2 */
       z = xx * xx;
       p = z * z * neval (z, J0_2N, NJ0_2N) / deval (z, J0_2D, NJ0_2D);
@@ -699,6 +701,28 @@ j0q (__float128 x)
       return p;
     }
 
+  /* X = x - pi/4
+     cos(X) = cos(x) cos(pi/4) + sin(x) sin(pi/4)
+     = 1/sqrt(2) * (cos(x) + sin(x))
+     sin(X) = sin(x) cos(pi/4) - cos(x) sin(pi/4)
+     = 1/sqrt(2) * (sin(x) - cos(x))
+     sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
+     cf. Fdlibm.  */
+  sincosq (xx, &s, &c);
+  ss = s - c;
+  cc = s + c;
+  if (xx <= FLT128_MAX / 2.0Q)
+    {
+      z = -cosq (xx + xx);
+      if ((s * c) < 0)
+	cc = z / ss;
+      else
+	ss = z / cc;
+    }
+
+  if (xx > 0x1p256Q)
+    return ONEOSQPI * cc / sqrtq (xx);
+
   xinv = 1.0Q / xx;
   z = xinv * xinv;
   if (xinv <= 0.25)
@@ -760,21 +784,6 @@ j0q (__float128 x)
   p = 1.0Q + z * p;
   q = z * xinv * q;
   q = q - 0.125Q * xinv;
-  /* X = x - pi/4
-     cos(X) = cos(x) cos(pi/4) + sin(x) sin(pi/4)
-     = 1/sqrt(2) * (cos(x) + sin(x))
-     sin(X) = sin(x) cos(pi/4) - cos(x) sin(pi/4)
-     = 1/sqrt(2) * (sin(x) - cos(x))
-     sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
-     cf. Fdlibm.  */
-  sincosq (xx, &s, &c);
-  ss = s - c;
-  cc = s + c;
-  z = - cosq (xx + xx);
-  if ((s * c) < 0)
-    cc = z / ss;
-  else
-    ss = z / cc;
   z = ONEOSQPI * (p * cc - q * ss) / sqrtq (xx);
   return z;
 }
@@ -817,17 +826,12 @@ y0q (__float128 x)
   __float128 xx, xinv, z, p, q, c, s, cc, ss;
 
   if (! finiteq (x))
-    {
-      if (x != x)
-	return x;
-      else
-	return 0.0Q;
-    }
+    return 1 / (x + x * x);
   if (x <= 0.0Q)
     {
       if (x < 0.0Q)
 	return (zero / (zero * x));
-      return -HUGE_VALQ + x;
+      return -1 / zero; /* -inf and divide by zero exception.  */
     }
   xx = fabsq (x);
   if (xx <= 0x1p-57)
@@ -841,6 +845,28 @@ y0q (__float128 x)
       return p;
     }
 
+  /* X = x - pi/4
+     cos(X) = cos(x) cos(pi/4) + sin(x) sin(pi/4)
+     = 1/sqrt(2) * (cos(x) + sin(x))
+     sin(X) = sin(x) cos(pi/4) - cos(x) sin(pi/4)
+     = 1/sqrt(2) * (sin(x) - cos(x))
+     sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
+     cf. Fdlibm.  */
+  sincosq (x, &s, &c);
+  ss = s - c;
+  cc = s + c;
+  if (xx <= FLT128_MAX / 2.0Q)
+    {
+      z = -cosq (x + x);
+      if ((s * c) < 0)
+	cc = z / ss;
+      else
+	ss = z / cc;
+    }
+
+  if (xx > 0x1p256Q)
+    return ONEOSQPI * ss / sqrtq (x);
+
   xinv = 1.0Q / xx;
   z = xinv * xinv;
   if (xinv <= 0.25)
@@ -902,21 +928,6 @@ y0q (__float128 x)
   p = 1.0Q + z * p;
   q = z * xinv * q;
   q = q - 0.125Q * xinv;
-  /* X = x - pi/4
-     cos(X) = cos(x) cos(pi/4) + sin(x) sin(pi/4)
-     = 1/sqrt(2) * (cos(x) + sin(x))
-     sin(X) = sin(x) cos(pi/4) - cos(x) sin(pi/4)
-     = 1/sqrt(2) * (sin(x) - cos(x))
-     sin(x) +- cos(x) = -cos(2x)/(sin(x) -+ cos(x))
-     cf. Fdlibm.  */
-  sincosq (x, &s, &c);
-  ss = s - c;
-  cc = s + c;
-  z = - cosq (x + x);
-  if ((s * c) < 0)
-    cc = z / ss;
-  else
-    ss = z / cc;
   z = ONEOSQPI * (p * ss + q * cc) / sqrtq (x);
   return z;
 }
diff --git a/libquadmath/math/j1q.c b/libquadmath/math/j1q.c
index eb599c949a9..5eb705084e2 100644
--- a/libquadmath/math/j1q.c
+++ b/libquadmath/math/j1q.c
@@ -95,6 +95,7 @@
     License along with this library; if not, write to the Free Software
     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307  USA */
 
+#include <errno.h>
 #include "quadmath-imp.h"
 
 /* 1 / sqrt(pi) */
@@ -687,13 +688,21 @@ j1q (__float128 x)
   if (! finiteq (x))
     {
       if (x != x)
-	return x;
+	return x + x;
       else
 	return 0.0Q;
     }
   if (x == 0.0Q)
     return x;
   xx = fabsq (x);
+  if (xx <= 0x1p-58Q)
+    {
+      __float128 ret = x * 0.5Q;
+      math_check_force_underflow (ret);
+      if (ret == 0)
+	errno = ERANGE;
+      return ret;
+    }
   if (xx <= 2.0Q)
     {
       /* 0 <= x <= 2 */
@@ -705,6 +714,32 @@ j1q (__float128 x)
       return p;
     }
 
+  /* X = x - 3 pi/4
+     cos(X) = cos(x) cos(3 pi/4) + sin(x) sin(3 pi/4)
+     = 1/sqrt(2) * (-cos(x) + sin(x))
+     sin(X) = sin(x) cos(3 pi/4) - cos(x) sin(3 pi/4)
+     = -1/sqrt(2) * (sin(x) + cos(x))
+     cf. Fdlibm.  */
+  sincosq (xx, &s, &c);
+  ss = -s - c;
+  cc = s - c;
+  if (xx <= FLT128_MAX / 2.0Q)
+    {
+      z = cosq (xx + xx);
+      if ((s * c) > 0)
+	cc = z / ss;
+      else
+	ss = z / cc;
+    }
+
+  if (xx > 0x1p256Q)
+    {
+      z = ONEOSQPI * cc / sqrtq (xx);
+      if (x < 0)
+	z = -z;
+      return z;
+    }
+
   xinv = 1.0Q / xx;
   z = xinv * xinv;
   if (xinv <= 0.25)
@@ -766,20 +801,6 @@ j1q (__float128 x)
   p = 1.0Q + z * p;
   q = z * q;
   q = q * xinv + 0.375Q * xinv;
-  /* X = x - 3 pi/4
-     cos(X) = cos(x) cos(3 pi/4) + sin(x) sin(3 pi/4)
-     = 1/sqrt(2) * (-cos(x) + sin(x))
-     sin(X) = sin(x) cos(3 pi/4) - cos(x) sin(3 pi/4)
-     = -1/sqrt(2) * (sin(x) + cos(x))
-     cf. Fdlibm.  */
-  sincosq (xx, &s, &c);
-  ss = -s - c;
-  cc = s - c;
-  z = cosq (xx + xx);
-  if ((s * c) > 0)
-    cc = z / ss;
-  else
-    ss = z / cc;
   z = ONEOSQPI * (p * cc - q * ss) / sqrtq (xx);
   if (x < 0)
     z = -z;
@@ -823,24 +844,25 @@ y1q (__float128 x)
   __float128 xx, xinv, z, p, q, c, s, cc, ss;
 
   if (! finiteq (x))
-    {
-      if (x != x)
-	return x;
-      else
-	return 0.0Q;
-    }
+    return 1 / (x + x * x);
   if (x <= 0.0Q)
     {
       if (x < 0.0Q)
 	return (zero / (zero * x));
-      return -HUGE_VALQ + x;
+      return -1 / zero; /* -inf and divide by zero exception.  */
     }
   xx = fabsq (x);
   if (xx <= 0x1p-114)
-    return -TWOOPI / x;
+    {
+      z = -TWOOPI / x;
+      if (isinfq (z))
+	errno = ERANGE;
+      return z;
+    }
   if (xx <= 2.0Q)
    {
       /* 0 <= x <= 2 */
+      /* FIXME:  SET_RESTORE_ROUNDL (FE_TONEAREST);  */
       z = xx * xx;
       p = xx * neval (z, Y0_2N, NY0_2N) / deval (z, Y0_2D, NY0_2D);
       p = -TWOOPI / xx + p;
@@ -848,6 +870,27 @@ y1q (__float128 x)
       return p;
     }
 
+  /* X = x - 3 pi/4
+     cos(X) = cos(x) cos(3 pi/4) + sin(x) sin(3 pi/4)
+     = 1/sqrt(2) * (-cos(x) + sin(x))
+     sin(X) = sin(x) cos(3 pi/4) - cos(x) sin(3 pi/4)
+     = -1/sqrt(2) * (sin(x) + cos(x))
+     cf. Fdlibm.  */
+  sincosq (xx, &s, &c);
+  ss = -s - c;
+  cc = s - c;
+  if (xx <= FLT128_MAX / 2.0Q)
+    {
+      z = cosq (xx + xx);
+      if ((s * c) > 0)
+	cc = z / ss;
+      else
+	ss = z / cc;
+    }
+
+  if (xx > 0x1p256Q)
+    return ONEOSQPI * ss / sqrtq (xx);
+
   xinv = 1.0Q / xx;
   z = xinv * xinv;
   if (xinv <= 0.25)
@@ -909,20 +952,6 @@ y1q (__float128 x)
   p = 1.0Q + z * p;
   q = z * q;
   q = q * xinv + 0.375Q * xinv;
-  /* X = x - 3 pi/4
-     cos(X) = cos(x) cos(3 pi/4) + sin(x) sin(3 pi/4)
-     = 1/sqrt(2) * (-cos(x) + sin(x))
-     sin(X) = sin(x) cos(3 pi/4) - cos(x) sin(3 pi/4)
-     = -1/sqrt(2) * (sin(x) + cos(x))
-     cf. Fdlibm.  */
-  sincosq (xx, &s, &c);
-  ss = -s - c;
-  cc = s - c;
-  z = cosq (xx + xx);
-  if ((s * c) > 0)
-    cc = z / ss;
-  else
-    ss = z / cc;
   z = ONEOSQPI * (p * ss + q * cc) / sqrtq (xx);
   return z;
 }
diff --git a/libquadmath/math/llrintq.c b/libquadmath/math/llrintq.c
index eef31d823b6..a6a0ae64bd3 100644
--- a/libquadmath/math/llrintq.c
+++ b/libquadmath/math/llrintq.c
@@ -1,9 +1,9 @@
 /* Round argument to nearest integral value according to current rounding
    direction.
-   Copyright (C) 1997, 1999, 2006 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
-   		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
+		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
 
    The GNU C Library is free software; you can redistribute it and/or
    modify it under the terms of the GNU Lesser General Public
@@ -33,7 +33,7 @@ llrintq (__float128 x)
 {
   int32_t j0;
   uint64_t i0,i1;
-  volatile __float128 w;
+  __float128 w;
   __float128 t;
   long long int result;
   int sx;
@@ -46,8 +46,23 @@ llrintq (__float128 x)
 
   if (j0 < (int32_t) (8 * sizeof (long long int)) - 1)
     {
-      w = two112[sx] + x;
-      t = w - two112[sx];
+#if defined FE_INVALID || defined FE_INEXACT
+      /* X < LLONG_MAX + 1 implied by J0 < 63.  */
+      if (x > (__float128) LLONG_MAX)
+	{
+	  /* In the event of overflow we must raise the "invalid"
+	     exception, but not "inexact".  */
+	  t = nearbyintq (x);
+#ifdef USE_FENV_H
+	  feraiseexcept (t == LLONG_MAX ? FE_INEXACT : FE_INVALID);
+#endif
+	}
+      else
+#endif
+	{
+	  w = two112[sx] + x;
+	  t = w - two112[sx];
+	}
       GET_FLT128_WORDS64 (i0, i1, t);
       j0 = ((i0 >> 48) & 0x7fff) - 0x3fff;
       i0 &= 0x0000ffffffffffffLL;
@@ -62,6 +77,24 @@ llrintq (__float128 x)
     }
   else
     {
+      /* The number is too large.  Unless it rounds to LLONG_MIN,
+	 FE_INVALID must be raised and the return value is
+	 unspecified.  */
+#if defined FE_INVALID || defined FE_INEXACT
+      if (x < (__float128) LLONG_MIN
+	  && x > (__float128) LLONG_MIN - 1.0Q)
+	{
+	  /* If truncation produces LLONG_MIN, the cast will not raise
+	     the exception, but may raise "inexact".  */
+	  t = nearbyintq (x);
+#ifdef USE_FENV_H
+	  feraiseexcept (t == LLONG_MIN ? FE_INEXACT : FE_INVALID);
+#endif
+	  return LLONG_MIN;
+	}
+
+#endif
+
       /* The number is too large.  It is left implementation defined
 	 what happens.  */
       return (long long int) x;
diff --git a/libquadmath/math/llroundq.c b/libquadmath/math/llroundq.c
index d22180d6bba..098fb9ef72b 100644
--- a/libquadmath/math/llroundq.c
+++ b/libquadmath/math/llroundq.c
@@ -1,8 +1,8 @@
 /* Round __float128 value to long long int.
-   Copyright (C) 1997, 1999, 2004 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
-   		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
+		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
 
    The GNU C Library is free software; you can redistribute it and/or
    modify it under the terms of the GNU Lesser General Public
@@ -59,13 +59,32 @@ llroundq (__float128 x)
 	  if (j0 == 48)
 	    result = (long long int) i0;
 	  else
-	    result = ((long long int) i0 << (j0 - 48)) | (j >> (112 - j0));
+	    {
+	      result = ((long long int) i0 << (j0 - 48)) | (j >> (112 - j0));
+#if defined FE_INVALID && defined USE_FENV_H
+	      if (sign == 1 && result == LLONG_MIN)
+		/* Rounding brought the value out of range.  */
+		feraiseexcept (FE_INVALID);
+#endif
+	    }
 	}
     }
   else
     {
-      /* The number is too large.  It is left implementation defined
-	 what happens.  */
+      /* The number is too large.  Unless it rounds to LLONG_MIN,
+	 FE_INVALID must be raised and the return value is
+	 unspecified.  */
+#ifdef FE_INVALID
+      if (x <= (__float128) LLONG_MIN - 0.5Q)
+	{
+	  /* If truncation produces LLONG_MIN, the cast will not raise
+	     the exception, but may raise "inexact".  */
+#ifdef USE_FENV_H
+	  feraiseexcept (FE_INVALID);
+#endif
+	  return LLONG_MIN;
+	}
+#endif
       return (long long int) x;
     }
 
diff --git a/libquadmath/math/log10q.c b/libquadmath/math/log10q.c
index 9eeb9ae3fc4..3afb1121267 100644
--- a/libquadmath/math/log10q.c
+++ b/libquadmath/math/log10q.c
@@ -188,12 +188,15 @@ log10q (__float128 x)
 /* Test for domain */
   GET_FLT128_WORDS64 (hx, lx, x);
   if (((hx & 0x7fffffffffffffffLL) | lx) == 0)
-    return (-1.0Q / (x - x));
+    return (-1.0Q / fabsq (x));		/* log10l(+-0)=-inf  */
   if (hx < 0)
     return (x - x) / (x - x);
   if (hx >= 0x7fff000000000000LL)
     return (x + x);
 
+  if (x == 1.0Q)
+    return 0.0Q;
+
 /* separate mantissa from exponent */
 
 /* Note, frexp is used so that denormal numbers
diff --git a/libquadmath/math/log1pq.c b/libquadmath/math/log1pq.c
index d8bff405dff..c59ceef5c9e 100644
--- a/libquadmath/math/log1pq.c
+++ b/libquadmath/math/log1pq.c
@@ -36,7 +36,7 @@
  *    IEEE      -1, 8       100000      1.9e-34     4.3e-35
  */
 
-/* Copyright 2001 by Stephen L. Moshier 
+/* Copyright 2001 by Stephen L. Moshier
 
     This library is free software; you can redistribute it and/or
     modify it under the terms of the GNU Lesser General Public
@@ -128,8 +128,8 @@ log1pq (__float128 xm1)
   /* Test for NaN or infinity input. */
   u.value = xm1;
   hx = u.words32.w0;
-  if (hx >= 0x7fff0000)
-    return xm1;
+  if ((hx & 0x7fffffff) >= 0x7fff0000)
+    return xm1 + fabsq (xm1);
 
   /* log1p(+- 0) = +- 0.  */
   if (((hx & 0x7fffffff) == 0)
@@ -138,17 +138,21 @@ log1pq (__float128 xm1)
 
   if ((hx & 0x7fffffff) < 0x3f8e0000)
     {
+      math_check_force_underflow (xm1);
       if ((int) xm1 == 0)
        return xm1;
     }
 
-  x = xm1 + 1.0Q;
+  if (xm1 >= 0x1p113Q)
+    x = xm1;
+  else
+    x = xm1 + 1.0Q;
 
   /* log1p(-1) = -inf */
   if (x <= 0.0Q)
     {
       if (x == 0.0Q)
-	return (-1.0Q / (x - x));
+	return (-1.0Q / zero);	/* log1p(-1) = -inf */
       else
 	return (zero / (x - x));
     }
diff --git a/libquadmath/math/log2q.c b/libquadmath/math/log2q.c
index f8275369b37..865f341f03d 100644
--- a/libquadmath/math/log2q.c
+++ b/libquadmath/math/log2q.c
@@ -181,12 +181,15 @@ log2q (__float128 x)
 /* Test for domain */
   GET_FLT128_WORDS64 (hx, lx, x);
   if (((hx & 0x7fffffffffffffffLL) | lx) == 0)
-    return (-1.0Q / (x - x));
+    return (-1.0Q / fabsq (x));		/* log2l(+-0)=-inf  */
   if (hx < 0)
     return (x - x) / (x - x);
   if (hx >= 0x7fff000000000000LL)
     return (x + x);
 
+  if (x == 1.0Q)
+    return 0.0Q;
+
 /* separate mantissa from exponent */
 
 /* Note, frexp is used so that denormal numbers
diff --git a/libquadmath/math/logq.c b/libquadmath/math/logq.c
index 7aae9b101ad..43249bb498a 100644
--- a/libquadmath/math/logq.c
+++ b/libquadmath/math/logq.c
@@ -212,9 +212,8 @@ logq (__float128 x)
     }
 
   /* Extract exponent and reduce domain to 0.703125 <= u < 1.40625  */
-  e = (int) (m >> 16) - (int) 0x3ffe;
-  m &= 0xffff;
-  u.words32.w0 = m | 0x3ffe0000;
+  u.value = frexpq (x, &e);
+  m = u.words32.w0 & 0xffff;
   m |= 0x10000;
   /* Find lookup table index k from high order bits of the significand. */
   if (m < 0x16800)
@@ -241,6 +240,8 @@ logq (__float128 x)
   /* On this interval the table is not used due to cancellation error.  */
   if ((x <= 1.0078125Q) && (x >= 0.9921875Q))
     {
+      if (x == 1.0Q)
+	return 0.0Q;
       z = x - 1.0Q;
       k = 64;
       t.value  = 1.0Q;
diff --git a/libquadmath/math/lrintq.c b/libquadmath/math/lrintq.c
index d1497ae3882..50a300554b9 100644
--- a/libquadmath/math/lrintq.c
+++ b/libquadmath/math/lrintq.c
@@ -1,9 +1,9 @@
 /* Round argument to nearest integral value according to current rounding
    direction.
-   Copyright (C) 1997, 1999, 2004, 2006 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
-   		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
+		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
 
    The GNU C Library is free software; you can redistribute it and/or
    modify it under the terms of the GNU Lesser General Public
@@ -33,7 +33,7 @@ lrintq (__float128 x)
 {
   int32_t j0;
   uint64_t i0,i1;
-  volatile __float128 w;
+  __float128 w;
   __float128 t;
   long int result;
   int sx;
@@ -44,25 +44,57 @@ lrintq (__float128 x)
   i0 &= 0x0000ffffffffffffLL;
   i0 |= 0x0001000000000000LL;
 
-  if (j0 < 48)
+  if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
     {
-      w = two112[sx] + x;
-      t = w - two112[sx];
-      GET_FLT128_WORDS64 (i0, i1, t);
-      j0 = ((i0 >> 48) & 0x7fff) - 0x3fff;
-      i0 &= 0x0000ffffffffffffLL;
-      i0 |= 0x0001000000000000LL;
+      if (j0 < 48)
+	{
+#if defined FE_INVALID || defined FE_INEXACT
+	  /* X < LONG_MAX + 1 implied by J0 < 31.  */
+	  if (sizeof (long int) == 4
+	      && x > (__float128) LONG_MAX)
+	    {
+	      /* In the event of overflow we must raise the "invalid"
+		 exception, but not "inexact".  */
+	      t = nearbyintq (x);
+#ifdef USE_FENV_H
+	      feraiseexcept (t == LONG_MAX ? FE_INEXACT : FE_INVALID);
+#endif
+	    }
+	  else
+#endif
+	    {
+	      w = two112[sx] + x;
+	      t = w - two112[sx];
+	    }
+	  GET_FLT128_WORDS64 (i0, i1, t);
+	  j0 = ((i0 >> 48) & 0x7fff) - 0x3fff;
+	  i0 &= 0x0000ffffffffffffLL;
+	  i0 |= 0x0001000000000000LL;
 
-      result = (j0 < 0 ? 0 : i0 >> (48 - j0));
-    }
-  else if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
-    {
-      if (j0 >= 112)
+	  result = (j0 < 0 ? 0 : i0 >> (48 - j0));
+	}
+      else if (j0 >= 112)
 	result = ((long int) i0 << (j0 - 48)) | (i1 << (j0 - 112));
       else
 	{
-	  w = two112[sx] + x;
-	  t = w - two112[sx];
+#if defined FE_INVALID || defined FE_INEXACT
+	  /* X < LONG_MAX + 1 implied by J0 < 63.  */
+	  if (sizeof (long int) == 8
+	      && x > (__float128) LONG_MAX)
+	    {
+	      /* In the event of overflow we must raise the "invalid"
+		 exception, but not "inexact".  */
+	      t = nearbyintq (x);
+#ifdef USE_FENV_H
+	      feraiseexcept (t == LONG_MAX ? FE_INEXACT : FE_INVALID);
+#endif
+	    }
+	  else
+#endif
+	    {
+	      w = two112[sx] + x;
+	      t = w - two112[sx];
+	    }
 	  GET_FLT128_WORDS64 (i0, i1, t);
 	  j0 = ((i0 >> 48) & 0x7fff) - 0x3fff;
 	  i0 &= 0x0000ffffffffffffLL;
@@ -76,8 +108,22 @@ lrintq (__float128 x)
     }
   else
     {
-      /* The number is too large.  It is left implementation defined
-	 what happens.  */
+      /* The number is too large.  Unless it rounds to LONG_MIN,
+	 FE_INVALID must be raised and the return value is
+	 unspecified.  */
+#if defined FE_INVALID || defined FE_INEXACT
+      if (x < (__float128) LONG_MIN
+	  && x > (__float128) LONG_MIN - 1.0Q)
+	{
+	  /* If truncation produces LONG_MIN, the cast will not raise
+	     the exception, but may raise "inexact".  */
+	  t = nearbyintq (x);
+#ifdef USE_FENV_H
+	  feraiseexcept (t == LONG_MIN ? FE_INEXACT : FE_INVALID);
+#endif
+	  return LONG_MIN;
+	}
+#endif
       return (long int) x;
     }
 
diff --git a/libquadmath/math/lroundq.c b/libquadmath/math/lroundq.c
index 59c883a1464..55285034cec 100644
--- a/libquadmath/math/lroundq.c
+++ b/libquadmath/math/lroundq.c
@@ -1,8 +1,8 @@
 /* Round __float128 value to long int.
-   Copyright (C) 1997, 1999, 2004 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
-   		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
+		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
 
    The GNU C Library is free software; you can redistribute it and/or
    modify it under the terms of the GNU Lesser General Public
@@ -36,19 +36,26 @@ lroundq (__float128 x)
   i0 &= 0x0000ffffffffffffLL;
   i0 |= 0x0001000000000000LL;
 
-  if (j0 < 48)
+  if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
     {
-      if (j0 < 0)
-	return j0 < -1 ? 0 : sign;
-      else
+      if (j0 < 48)
 	{
-	  i0 += 0x0000800000000000LL >> j0;
-	  result = i0 >> (48 - j0);
+	  if (j0 < 0)
+	    return j0 < -1 ? 0 : sign;
+	  else
+	    {
+	      i0 += 0x0000800000000000LL >> j0;
+	      result = i0 >> (48 - j0);
+#if defined FE_INVALID && defined USE_FENV_H
+	      if (sizeof (long int) == 4
+		  && sign == 1
+		  && result == LONG_MIN)
+		/* Rounding brought the value out of range.  */
+		feraiseexcept (FE_INVALID);
+#endif
+	    }
 	}
-    }
-  else if (j0 < (int32_t) (8 * sizeof (long int)) - 1)
-    {
-      if (j0 >= 112)
+      else if (j0 >= 112)
 	result = ((long int) i0 << (j0 - 48)) | (i1 << (j0 - 112));
       else
 	{
@@ -59,13 +66,34 @@ lroundq (__float128 x)
 	  if (j0 == 48)
 	    result = (long int) i0;
 	  else
-	    result = ((long int) i0 << (j0 - 48)) | (j >> (112 - j0));
+	    {
+	      result = ((long int) i0 << (j0 - 48)) | (j >> (112 - j0));
+#if defined FE_INVALID && defined USE_FENV_H
+	      if (sizeof (long int) == 8
+		  && sign == 1
+		  && result == LONG_MIN)
+		/* Rounding brought the value out of range.  */
+		feraiseexcept (FE_INVALID);
+#endif
+	    }
 	}
     }
   else
     {
-      /* The number is too large.  It is left implementation defined
-	 what happens.  */
+      /* The number is too large.  Unless it rounds to LONG_MIN,
+	 FE_INVALID must be raised and the return value is
+	 unspecified.  */
+#ifdef FE_INVALID
+      if (x <= (__float128) LONG_MIN - 0.5Q)
+	{
+	  /* If truncation produces LONG_MIN, the cast will not raise
+	     the exception, but may raise "inexact".  */
+#ifdef USE_FENV_H
+	  feraiseexcept (FE_INVALID);
+#endif
+	  return LONG_MIN;
+	}
+#endif
       return (long int) x;
     }
 
diff --git a/libquadmath/math/nearbyintq.c b/libquadmath/math/nearbyintq.c
index 207124808e7..b250927ea2a 100644
--- a/libquadmath/math/nearbyintq.c
+++ b/libquadmath/math/nearbyintq.c
@@ -44,7 +44,7 @@ nearbyintq(__float128 x)
 	fenv_t env;
 #endif
 	int64_t i0,j0,sx;
-	uint64_t i1;
+	uint64_t i1 __attribute__ ((unused));
 	__float128 w,t;
 	GET_FLT128_WORDS64(i0,i1,x);
 	sx = (((uint64_t)i0)>>63);
@@ -56,6 +56,7 @@ nearbyintq(__float128 x)
 #endif
 	        w = TWO112[sx]+x;
 	        t = w-TWO112[sx];
+		math_force_eval (t);
 #ifdef USE_FENV_H
 	        fesetenv (&env);
 #endif
@@ -72,6 +73,7 @@ nearbyintq(__float128 x)
 #endif
 	w = TWO112[sx]+x;
 	t = w-TWO112[sx];
+	math_force_eval (t);
 #ifdef USE_FENV_H	
 	fesetenv (&env);
 #endif
diff --git a/libquadmath/math/nextafterq.c b/libquadmath/math/nextafterq.c
index 04d63deb8e2..a030e9c6444 100644
--- a/libquadmath/math/nextafterq.c
+++ b/libquadmath/math/nextafterq.c
@@ -13,6 +13,7 @@
  * ====================================================
  */
 
+#include <errno.h>
 #include "quadmath-imp.h"
 
 __float128
@@ -54,9 +55,15 @@ nextafterq (__float128 x, __float128 y)
 	    }
 	}
 	hy = hx&0x7fff000000000000LL;
-	if(hy==0x7fff000000000000LL) return x+x;/* overflow  */
+	if(hy==0x7fff000000000000LL) {
+	    __float128 u = x + x;		/* overflow  */
+	    math_force_eval (u);
+	    errno = ERANGE;
+	}
 	if(hy==0) {
-	    /* here we should raise an underflow flag */
+	    __float128 u = x*x;			/* underflow */
+	    math_force_eval (u);		/* raise underflow flag */
+	    errno = ERANGE;
 	}
 	SET_FLT128_WORDS64(x,hx,lx);
 	return x;
diff --git a/libquadmath/math/powq.c b/libquadmath/math/powq.c
index dd44b7c175a..acad2057516 100644
--- a/libquadmath/math/powq.c
+++ b/libquadmath/math/powq.c
@@ -147,7 +147,7 @@ __float128
 powq (__float128 x, __float128 y)
 {
   __float128 z, ax, z_h, z_l, p_h, p_l;
-  __float128 y1, t1, t2, r, s, t, u, v, w;
+  __float128 y1, t1, t2, r, s, sgn, t, u, v, w;
   __float128 s2, s_h, s_l, t_h, t_l, ay;
   int32_t i, j, k, yisint, n;
   uint32_t ix, iy;
@@ -261,6 +261,11 @@ powq (__float128 x, __float128 y)
   if (((((uint32_t) hx >> 31) - 1) | yisint) == 0)
     return (x - x) / (x - x);
 
+  /* sgn (sign of result -ve**odd) = -1 else = 1 */
+  sgn = one;
+  if (((((u_int32_t) hx >> 31) - 1) | (yisint - 1)) == 0)
+    sgn = -one;				/* (-ve)**(odd int) */
+
   /* |y| is huge.
      2^-16495 = 1/2 of smallest representable value.
      If (1 - 1/131072)^y underflows, y > 1.4986e9 */
@@ -276,9 +281,9 @@ powq (__float128 x, __float128 y)
 	}
       /* over/underflow if x is not close to one */
       if (ix < 0x3ffeffff)
-	return (hy < 0) ? huge * huge : tiny * tiny;
+	return (hy < 0) ? sgn * huge * huge : sgn * tiny * tiny;
       if (ix > 0x3fff0000)
-	return (hy > 0) ? huge * huge : tiny * tiny;
+	return (hy > 0) ? sgn * huge * huge : sgn * tiny * tiny;
     }
 
   ay = y > 0 ? y : -y;
@@ -365,11 +370,6 @@ powq (__float128 x, __float128 y)
   t1 = o.value;
   t2 = z_l - (((t1 - t) - dp_h[k]) - z_h);
 
-  /* s (sign of result -ve**odd) = -1 else = 1 */
-  s = one;
-  if (((((uint32_t) hx >> 31) - 1) | (yisint - 1)) == 0)
-    s = -one;			/* (-ve)**(odd int) */
-
   /* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
   y1 = y;
   o.value = y1;
@@ -385,11 +385,11 @@ powq (__float128 x, __float128 y)
     {
       /* if z > 16384 */
       if (((j - 0x400d0000) | o.words32.w1 | o.words32.w2 | o.words32.w3) != 0)
-	return s * huge * huge;	/* overflow */
+	return sgn * huge * huge;	/* overflow */
       else
 	{
 	  if (p_l + ovt > z - p_h)
-	    return s * huge * huge;	/* overflow */
+	    return sgn * huge * huge;	/* overflow */
 	}
     }
   else if ((j & 0x7fffffff) >= 0x400d01b9)	/* z <= -16495 */
@@ -397,11 +397,11 @@ powq (__float128 x, __float128 y)
       /* z < -16495 */
       if (((j - 0xc00d01bc) | o.words32.w1 | o.words32.w2 | o.words32.w3)
 	  != 0)
-	return s * tiny * tiny;	/* underflow */
+	return sgn * tiny * tiny;	/* underflow */
       else
 	{
 	  if (p_l <= z - p_h)
-	    return s * tiny * tiny;	/* underflow */
+	    return sgn * tiny * tiny;	/* underflow */
 	}
     }
   /* compute 2**(p_h+p_l) */
@@ -434,11 +434,15 @@ powq (__float128 x, __float128 y)
   j = o.words32.w0;
   j += (n << 16);
   if ((j >> 16) <= 0)
-    z = scalbnq (z, n);	/* subnormal output */
+    {
+      z = scalbnq (z, n);	/* subnormal output */
+      __float128 force_underflow = z * z;
+      math_force_eval (force_underflow);
+    }
   else
     {
       o.words32.w0 = j;
       z = o.value;
     }
-  return s * z;
+  return sgn * z;
 }
diff --git a/libquadmath/math/rem_pio2q.c b/libquadmath/math/rem_pio2q.c
index 60653c8d1d3..3308b218473 100644
--- a/libquadmath/math/rem_pio2q.c
+++ b/libquadmath/math/rem_pio2q.c
@@ -312,7 +312,7 @@ recompute:
 
 
 /* Quad-precision floating point argument reduction.
-   Copyright (C) 1999 Free Software Foundation, Inc.
+   Copyright (C) 1999-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Jakub Jelinek <jj@ultra.linux.cz>
 
@@ -332,176 +332,176 @@ recompute:
    02111-1307 USA.  */
 
 /*
- * Table of constants for 2/pi, 5628 hexadecimal digits of 2/pi 
+ * Table of constants for 2/pi, 5628 hexadecimal digits of 2/pi
  */
 static const int32_t two_over_pi[] = {
-0xa2f983, 0x6e4e44, 0x1529fc, 0x2757d1, 0xf534dd, 0xc0db62, 
-0x95993c, 0x439041, 0xfe5163, 0xabdebb, 0xc561b7, 0x246e3a, 
-0x424dd2, 0xe00649, 0x2eea09, 0xd1921c, 0xfe1deb, 0x1cb129, 
-0xa73ee8, 0x8235f5, 0x2ebb44, 0x84e99c, 0x7026b4, 0x5f7e41, 
-0x3991d6, 0x398353, 0x39f49c, 0x845f8b, 0xbdf928, 0x3b1ff8, 
-0x97ffde, 0x05980f, 0xef2f11, 0x8b5a0a, 0x6d1f6d, 0x367ecf, 
-0x27cb09, 0xb74f46, 0x3f669e, 0x5fea2d, 0x7527ba, 0xc7ebe5, 
-0xf17b3d, 0x0739f7, 0x8a5292, 0xea6bfb, 0x5fb11f, 0x8d5d08, 
-0x560330, 0x46fc7b, 0x6babf0, 0xcfbc20, 0x9af436, 0x1da9e3, 
-0x91615e, 0xe61b08, 0x659985, 0x5f14a0, 0x68408d, 0xffd880, 
-0x4d7327, 0x310606, 0x1556ca, 0x73a8c9, 0x60e27b, 0xc08c6b, 
-0x47c419, 0xc367cd, 0xdce809, 0x2a8359, 0xc4768b, 0x961ca6, 
-0xddaf44, 0xd15719, 0x053ea5, 0xff0705, 0x3f7e33, 0xe832c2, 
-0xde4f98, 0x327dbb, 0xc33d26, 0xef6b1e, 0x5ef89f, 0x3a1f35, 
-0xcaf27f, 0x1d87f1, 0x21907c, 0x7c246a, 0xfa6ed5, 0x772d30, 
-0x433b15, 0xc614b5, 0x9d19c3, 0xc2c4ad, 0x414d2c, 0x5d000c, 
-0x467d86, 0x2d71e3, 0x9ac69b, 0x006233, 0x7cd2b4, 0x97a7b4, 
-0xd55537, 0xf63ed7, 0x1810a3, 0xfc764d, 0x2a9d64, 0xabd770, 
-0xf87c63, 0x57b07a, 0xe71517, 0x5649c0, 0xd9d63b, 0x3884a7, 
-0xcb2324, 0x778ad6, 0x23545a, 0xb91f00, 0x1b0af1, 0xdfce19, 
-0xff319f, 0x6a1e66, 0x615799, 0x47fbac, 0xd87f7e, 0xb76522, 
-0x89e832, 0x60bfe6, 0xcdc4ef, 0x09366c, 0xd43f5d, 0xd7de16, 
-0xde3b58, 0x929bde, 0x2822d2, 0xe88628, 0x4d58e2, 0x32cac6, 
-0x16e308, 0xcb7de0, 0x50c017, 0xa71df3, 0x5be018, 0x34132e, 
-0x621283, 0x014883, 0x5b8ef5, 0x7fb0ad, 0xf2e91e, 0x434a48, 
-0xd36710, 0xd8ddaa, 0x425fae, 0xce616a, 0xa4280a, 0xb499d3, 
-0xf2a606, 0x7f775c, 0x83c2a3, 0x883c61, 0x78738a, 0x5a8caf, 
-0xbdd76f, 0x63a62d, 0xcbbff4, 0xef818d, 0x67c126, 0x45ca55, 
-0x36d9ca, 0xd2a828, 0x8d61c2, 0x77c912, 0x142604, 0x9b4612, 
-0xc459c4, 0x44c5c8, 0x91b24d, 0xf31700, 0xad43d4, 0xe54929, 
-0x10d5fd, 0xfcbe00, 0xcc941e, 0xeece70, 0xf53e13, 0x80f1ec, 
-0xc3e7b3, 0x28f8c7, 0x940593, 0x3e71c1, 0xb3092e, 0xf3450b, 
-0x9c1288, 0x7b20ab, 0x9fb52e, 0xc29247, 0x2f327b, 0x6d550c, 
-0x90a772, 0x1fe76b, 0x96cb31, 0x4a1679, 0xe27941, 0x89dff4, 
-0x9794e8, 0x84e6e2, 0x973199, 0x6bed88, 0x365f5f, 0x0efdbb, 
-0xb49a48, 0x6ca467, 0x427271, 0x325d8d, 0xb8159f, 0x09e5bc, 
-0x25318d, 0x3974f7, 0x1c0530, 0x010c0d, 0x68084b, 0x58ee2c, 
-0x90aa47, 0x02e774, 0x24d6bd, 0xa67df7, 0x72486e, 0xef169f, 
-0xa6948e, 0xf691b4, 0x5153d1, 0xf20acf, 0x339820, 0x7e4bf5, 
-0x6863b2, 0x5f3edd, 0x035d40, 0x7f8985, 0x295255, 0xc06437, 
-0x10d86d, 0x324832, 0x754c5b, 0xd4714e, 0x6e5445, 0xc1090b, 
-0x69f52a, 0xd56614, 0x9d0727, 0x50045d, 0xdb3bb4, 0xc576ea, 
-0x17f987, 0x7d6b49, 0xba271d, 0x296996, 0xacccc6, 0x5414ad, 
-0x6ae290, 0x89d988, 0x50722c, 0xbea404, 0x940777, 0x7030f3, 
-0x27fc00, 0xa871ea, 0x49c266, 0x3de064, 0x83dd97, 0x973fa3, 
-0xfd9443, 0x8c860d, 0xde4131, 0x9d3992, 0x8c70dd, 0xe7b717, 
-0x3bdf08, 0x2b3715, 0xa0805c, 0x93805a, 0x921110, 0xd8e80f, 
-0xaf806c, 0x4bffdb, 0x0f9038, 0x761859, 0x15a562, 0xbbcb61, 
-0xb989c7, 0xbd4010, 0x04f2d2, 0x277549, 0xf6b6eb, 0xbb22db, 
-0xaa140a, 0x2f2689, 0x768364, 0x333b09, 0x1a940e, 0xaa3a51, 
-0xc2a31d, 0xaeedaf, 0x12265c, 0x4dc26d, 0x9c7a2d, 0x9756c0, 
-0x833f03, 0xf6f009, 0x8c402b, 0x99316d, 0x07b439, 0x15200c, 
-0x5bc3d8, 0xc492f5, 0x4badc6, 0xa5ca4e, 0xcd37a7, 0x36a9e6, 
-0x9492ab, 0x6842dd, 0xde6319, 0xef8c76, 0x528b68, 0x37dbfc, 
-0xaba1ae, 0x3115df, 0xa1ae00, 0xdafb0c, 0x664d64, 0xb705ed, 
-0x306529, 0xbf5657, 0x3aff47, 0xb9f96a, 0xf3be75, 0xdf9328, 
-0x3080ab, 0xf68c66, 0x15cb04, 0x0622fa, 0x1de4d9, 0xa4b33d, 
-0x8f1b57, 0x09cd36, 0xe9424e, 0xa4be13, 0xb52333, 0x1aaaf0, 
-0xa8654f, 0xa5c1d2, 0x0f3f0b, 0xcd785b, 0x76f923, 0x048b7b, 
-0x721789, 0x53a6c6, 0xe26e6f, 0x00ebef, 0x584a9b, 0xb7dac4, 
-0xba66aa, 0xcfcf76, 0x1d02d1, 0x2df1b1, 0xc1998c, 0x77adc3, 
-0xda4886, 0xa05df7, 0xf480c6, 0x2ff0ac, 0x9aecdd, 0xbc5c3f, 
-0x6dded0, 0x1fc790, 0xb6db2a, 0x3a25a3, 0x9aaf00, 0x9353ad, 
-0x0457b6, 0xb42d29, 0x7e804b, 0xa707da, 0x0eaa76, 0xa1597b, 
-0x2a1216, 0x2db7dc, 0xfde5fa, 0xfedb89, 0xfdbe89, 0x6c76e4, 
-0xfca906, 0x70803e, 0x156e85, 0xff87fd, 0x073e28, 0x336761, 
-0x86182a, 0xeabd4d, 0xafe7b3, 0x6e6d8f, 0x396795, 0x5bbf31, 
-0x48d784, 0x16df30, 0x432dc7, 0x356125, 0xce70c9, 0xb8cb30, 
-0xfd6cbf, 0xa200a4, 0xe46c05, 0xa0dd5a, 0x476f21, 0xd21262, 
-0x845cb9, 0x496170, 0xe0566b, 0x015299, 0x375550, 0xb7d51e, 
-0xc4f133, 0x5f6e13, 0xe4305d, 0xa92e85, 0xc3b21d, 0x3632a1, 
-0xa4b708, 0xd4b1ea, 0x21f716, 0xe4698f, 0x77ff27, 0x80030c, 
-0x2d408d, 0xa0cd4f, 0x99a520, 0xd3a2b3, 0x0a5d2f, 0x42f9b4, 
-0xcbda11, 0xd0be7d, 0xc1db9b, 0xbd17ab, 0x81a2ca, 0x5c6a08, 
-0x17552e, 0x550027, 0xf0147f, 0x8607e1, 0x640b14, 0x8d4196, 
-0xdebe87, 0x2afdda, 0xb6256b, 0x34897b, 0xfef305, 0x9ebfb9, 
-0x4f6a68, 0xa82a4a, 0x5ac44f, 0xbcf82d, 0x985ad7, 0x95c7f4, 
-0x8d4d0d, 0xa63a20, 0x5f57a4, 0xb13f14, 0x953880, 0x0120cc, 
-0x86dd71, 0xb6dec9, 0xf560bf, 0x11654d, 0x6b0701, 0xacb08c, 
-0xd0c0b2, 0x485551, 0x0efb1e, 0xc37295, 0x3b06a3, 0x3540c0, 
-0x7bdc06, 0xcc45e0, 0xfa294e, 0xc8cad6, 0x41f3e8, 0xde647c, 
-0xd8649b, 0x31bed9, 0xc397a4, 0xd45877, 0xc5e369, 0x13daf0, 
-0x3c3aba, 0x461846, 0x5f7555, 0xf5bdd2, 0xc6926e, 0x5d2eac, 
-0xed440e, 0x423e1c, 0x87c461, 0xe9fd29, 0xf3d6e7, 0xca7c22, 
-0x35916f, 0xc5e008, 0x8dd7ff, 0xe26a6e, 0xc6fdb0, 0xc10893, 
-0x745d7c, 0xb2ad6b, 0x9d6ecd, 0x7b723e, 0x6a11c6, 0xa9cff7, 
-0xdf7329, 0xbac9b5, 0x5100b7, 0x0db2e2, 0x24ba74, 0x607de5, 
-0x8ad874, 0x2c150d, 0x0c1881, 0x94667e, 0x162901, 0x767a9f, 
-0xbefdfd, 0xef4556, 0x367ed9, 0x13d9ec, 0xb9ba8b, 0xfc97c4, 
-0x27a831, 0xc36ef1, 0x36c594, 0x56a8d8, 0xb5a8b4, 0x0ecccf, 
-0x2d8912, 0x34576f, 0x89562c, 0xe3ce99, 0xb920d6, 0xaa5e6b, 
-0x9c2a3e, 0xcc5f11, 0x4a0bfd, 0xfbf4e1, 0x6d3b8e, 0x2c86e2, 
-0x84d4e9, 0xa9b4fc, 0xd1eeef, 0xc9352e, 0x61392f, 0x442138, 
-0xc8d91b, 0x0afc81, 0x6a4afb, 0xd81c2f, 0x84b453, 0x8c994e, 
-0xcc2254, 0xdc552a, 0xd6c6c0, 0x96190b, 0xb8701a, 0x649569, 
-0x605a26, 0xee523f, 0x0f117f, 0x11b5f4, 0xf5cbfc, 0x2dbc34, 
-0xeebc34, 0xcc5de8, 0x605edd, 0x9b8e67, 0xef3392, 0xb817c9, 
-0x9b5861, 0xbc57e1, 0xc68351, 0x103ed8, 0x4871dd, 0xdd1c2d, 
-0xa118af, 0x462c21, 0xd7f359, 0x987ad9, 0xc0549e, 0xfa864f, 
-0xfc0656, 0xae79e5, 0x362289, 0x22ad38, 0xdc9367, 0xaae855, 
-0x382682, 0x9be7ca, 0xa40d51, 0xb13399, 0x0ed7a9, 0x480569, 
-0xf0b265, 0xa7887f, 0x974c88, 0x36d1f9, 0xb39221, 0x4a827b, 
-0x21cf98, 0xdc9f40, 0x5547dc, 0x3a74e1, 0x42eb67, 0xdf9dfe, 
-0x5fd45e, 0xa4677b, 0x7aacba, 0xa2f655, 0x23882b, 0x55ba41, 
-0x086e59, 0x862a21, 0x834739, 0xe6e389, 0xd49ee5, 0x40fb49, 
-0xe956ff, 0xca0f1c, 0x8a59c5, 0x2bfa94, 0xc5c1d3, 0xcfc50f, 
-0xae5adb, 0x86c547, 0x624385, 0x3b8621, 0x94792c, 0x876110, 
-0x7b4c2a, 0x1a2c80, 0x12bf43, 0x902688, 0x893c78, 0xe4c4a8, 
-0x7bdbe5, 0xc23ac4, 0xeaf426, 0x8a67f7, 0xbf920d, 0x2ba365, 
-0xb1933d, 0x0b7cbd, 0xdc51a4, 0x63dd27, 0xdde169, 0x19949a, 
-0x9529a8, 0x28ce68, 0xb4ed09, 0x209f44, 0xca984e, 0x638270, 
-0x237c7e, 0x32b90f, 0x8ef5a7, 0xe75614, 0x08f121, 0x2a9db5, 
-0x4d7e6f, 0x5119a5, 0xabf9b5, 0xd6df82, 0x61dd96, 0x023616, 
-0x9f3ac4, 0xa1a283, 0x6ded72, 0x7a8d39, 0xa9b882, 0x5c326b, 
-0x5b2746, 0xed3400, 0x7700d2, 0x55f4fc, 0x4d5901, 0x8071e0, 
-0xe13f89, 0xb295f3, 0x64a8f1, 0xaea74b, 0x38fc4c, 0xeab2bb, 
-0x47270b, 0xabc3a7, 0x34ba60, 0x52dd34, 0xf8563a, 0xeb7e8a, 
-0x31bb36, 0x5895b7, 0x47f7a9, 0x94c3aa, 0xd39225, 0x1e7f3e, 
-0xd8974e, 0xbba94f, 0xd8ae01, 0xe661b4, 0x393d8e, 0xa523aa, 
-0x33068e, 0x1633b5, 0x3bb188, 0x1d3a9d, 0x4013d0, 0xcc1be5, 
-0xf862e7, 0x3bf28f, 0x39b5bf, 0x0bc235, 0x22747e, 0xa247c0, 
-0xd52d1f, 0x19add3, 0x9094df, 0x9311d0, 0xb42b25, 0x496db2, 
-0xe264b2, 0x5ef135, 0x3bc6a4, 0x1a4ad0, 0xaac92e, 0x64e886, 
-0x573091, 0x982cfb, 0x311b1a, 0x08728b, 0xbdcee1, 0x60e142, 
-0xeb641d, 0xd0bba3, 0xe559d4, 0x597b8c, 0x2a4483, 0xf332ba, 
-0xf84867, 0x2c8d1b, 0x2fa9b0, 0x50f3dd, 0xf9f573, 0xdb61b4, 
-0xfe233e, 0x6c41a6, 0xeea318, 0x775a26, 0xbc5e5c, 0xcea708, 
-0x94dc57, 0xe20196, 0xf1e839, 0xbe4851, 0x5d2d2f, 0x4e9555, 
-0xd96ec2, 0xe7d755, 0x6304e0, 0xc02e0e, 0xfc40a0, 0xbbf9b3, 
-0x7125a7, 0x222dfb, 0xf619d8, 0x838c1c, 0x6619e6, 0xb20d55, 
-0xbb5137, 0x79e809, 0xaf9149, 0x0d73de, 0x0b0da5, 0xce7f58, 
-0xac1934, 0x724667, 0x7a1a13, 0x9e26bc, 0x4555e7, 0x585cb5, 
-0x711d14, 0x486991, 0x480d60, 0x56adab, 0xd62f64, 0x96ee0c, 
-0x212ff3, 0x5d6d88, 0xa67684, 0x95651e, 0xab9e0a, 0x4ddefe, 
-0x571010, 0x836a39, 0xf8ea31, 0x9e381d, 0xeac8b1, 0xcac96b, 
-0x37f21e, 0xd505e9, 0x984743, 0x9fc56c, 0x0331b7, 0x3b8bf8, 
-0x86e56a, 0x8dc343, 0x6230e7, 0x93cfd5, 0x6a8f2d, 0x733005, 
-0x1af021, 0xa09fcb, 0x7415a1, 0xd56b23, 0x6ff725, 0x2f4bc7, 
-0xb8a591, 0x7fac59, 0x5c55de, 0x212c38, 0xb13296, 0x5cff50, 
-0x366262, 0xfa7b16, 0xf4d9a6, 0x2acfe7, 0xf07403, 0xd4d604, 
-0x6fd916, 0x31b1bf, 0xcbb450, 0x5bd7c8, 0x0ce194, 0x6bd643, 
-0x4fd91c, 0xdf4543, 0x5f3453, 0xe2b5aa, 0xc9aec8, 0x131485, 
-0xf9d2bf, 0xbadb9e, 0x76f5b9, 0xaf15cf, 0xca3182, 0x14b56d, 
-0xe9fe4d, 0x50fc35, 0xf5aed5, 0xa2d0c1, 0xc96057, 0x192eb6, 
-0xe91d92, 0x07d144, 0xaea3c6, 0x343566, 0x26d5b4, 0x3161e2, 
-0x37f1a2, 0x209eff, 0x958e23, 0x493798, 0x35f4a6, 0x4bdc02, 
-0xc2be13, 0xbe80a0, 0x0b72a3, 0x115c5f, 0x1e1bd1, 0x0db4d3, 
-0x869e85, 0x96976b, 0x2ac91f, 0x8a26c2, 0x3070f0, 0x041412, 
-0xfc9fa5, 0xf72a38, 0x9c6878, 0xe2aa76, 0x50cfe1, 0x559274, 
-0x934e38, 0x0a92f7, 0x5533f0, 0xa63db4, 0x399971, 0xe2b755, 
-0xa98a7c, 0x008f19, 0xac54d2, 0x2ea0b4, 0xf5f3e0, 0x60c849, 
-0xffd269, 0xae52ce, 0x7a5fdd, 0xe9ce06, 0xfb0ae8, 0xa50cce, 
-0xea9d3e, 0x3766dd, 0xb834f5, 0x0da090, 0x846f88, 0x4ae3d5, 
-0x099a03, 0x2eae2d, 0xfcb40a, 0xfb9b33, 0xe281dd, 0x1b16ba, 
-0xd8c0af, 0xd96b97, 0xb52dc9, 0x9c277f, 0x5951d5, 0x21ccd6, 
-0xb6496b, 0x584562, 0xb3baf2, 0xa1a5c4, 0x7ca2cf, 0xa9b93d, 
-0x7b7b89, 0x483d38, 
+0xa2f983, 0x6e4e44, 0x1529fc, 0x2757d1, 0xf534dd, 0xc0db62,
+0x95993c, 0x439041, 0xfe5163, 0xabdebb, 0xc561b7, 0x246e3a,
+0x424dd2, 0xe00649, 0x2eea09, 0xd1921c, 0xfe1deb, 0x1cb129,
+0xa73ee8, 0x8235f5, 0x2ebb44, 0x84e99c, 0x7026b4, 0x5f7e41,
+0x3991d6, 0x398353, 0x39f49c, 0x845f8b, 0xbdf928, 0x3b1ff8,
+0x97ffde, 0x05980f, 0xef2f11, 0x8b5a0a, 0x6d1f6d, 0x367ecf,
+0x27cb09, 0xb74f46, 0x3f669e, 0x5fea2d, 0x7527ba, 0xc7ebe5,
+0xf17b3d, 0x0739f7, 0x8a5292, 0xea6bfb, 0x5fb11f, 0x8d5d08,
+0x560330, 0x46fc7b, 0x6babf0, 0xcfbc20, 0x9af436, 0x1da9e3,
+0x91615e, 0xe61b08, 0x659985, 0x5f14a0, 0x68408d, 0xffd880,
+0x4d7327, 0x310606, 0x1556ca, 0x73a8c9, 0x60e27b, 0xc08c6b,
+0x47c419, 0xc367cd, 0xdce809, 0x2a8359, 0xc4768b, 0x961ca6,
+0xddaf44, 0xd15719, 0x053ea5, 0xff0705, 0x3f7e33, 0xe832c2,
+0xde4f98, 0x327dbb, 0xc33d26, 0xef6b1e, 0x5ef89f, 0x3a1f35,
+0xcaf27f, 0x1d87f1, 0x21907c, 0x7c246a, 0xfa6ed5, 0x772d30,
+0x433b15, 0xc614b5, 0x9d19c3, 0xc2c4ad, 0x414d2c, 0x5d000c,
+0x467d86, 0x2d71e3, 0x9ac69b, 0x006233, 0x7cd2b4, 0x97a7b4,
+0xd55537, 0xf63ed7, 0x1810a3, 0xfc764d, 0x2a9d64, 0xabd770,
+0xf87c63, 0x57b07a, 0xe71517, 0x5649c0, 0xd9d63b, 0x3884a7,
+0xcb2324, 0x778ad6, 0x23545a, 0xb91f00, 0x1b0af1, 0xdfce19,
+0xff319f, 0x6a1e66, 0x615799, 0x47fbac, 0xd87f7e, 0xb76522,
+0x89e832, 0x60bfe6, 0xcdc4ef, 0x09366c, 0xd43f5d, 0xd7de16,
+0xde3b58, 0x929bde, 0x2822d2, 0xe88628, 0x4d58e2, 0x32cac6,
+0x16e308, 0xcb7de0, 0x50c017, 0xa71df3, 0x5be018, 0x34132e,
+0x621283, 0x014883, 0x5b8ef5, 0x7fb0ad, 0xf2e91e, 0x434a48,
+0xd36710, 0xd8ddaa, 0x425fae, 0xce616a, 0xa4280a, 0xb499d3,
+0xf2a606, 0x7f775c, 0x83c2a3, 0x883c61, 0x78738a, 0x5a8caf,
+0xbdd76f, 0x63a62d, 0xcbbff4, 0xef818d, 0x67c126, 0x45ca55,
+0x36d9ca, 0xd2a828, 0x8d61c2, 0x77c912, 0x142604, 0x9b4612,
+0xc459c4, 0x44c5c8, 0x91b24d, 0xf31700, 0xad43d4, 0xe54929,
+0x10d5fd, 0xfcbe00, 0xcc941e, 0xeece70, 0xf53e13, 0x80f1ec,
+0xc3e7b3, 0x28f8c7, 0x940593, 0x3e71c1, 0xb3092e, 0xf3450b,
+0x9c1288, 0x7b20ab, 0x9fb52e, 0xc29247, 0x2f327b, 0x6d550c,
+0x90a772, 0x1fe76b, 0x96cb31, 0x4a1679, 0xe27941, 0x89dff4,
+0x9794e8, 0x84e6e2, 0x973199, 0x6bed88, 0x365f5f, 0x0efdbb,
+0xb49a48, 0x6ca467, 0x427271, 0x325d8d, 0xb8159f, 0x09e5bc,
+0x25318d, 0x3974f7, 0x1c0530, 0x010c0d, 0x68084b, 0x58ee2c,
+0x90aa47, 0x02e774, 0x24d6bd, 0xa67df7, 0x72486e, 0xef169f,
+0xa6948e, 0xf691b4, 0x5153d1, 0xf20acf, 0x339820, 0x7e4bf5,
+0x6863b2, 0x5f3edd, 0x035d40, 0x7f8985, 0x295255, 0xc06437,
+0x10d86d, 0x324832, 0x754c5b, 0xd4714e, 0x6e5445, 0xc1090b,
+0x69f52a, 0xd56614, 0x9d0727, 0x50045d, 0xdb3bb4, 0xc576ea,
+0x17f987, 0x7d6b49, 0xba271d, 0x296996, 0xacccc6, 0x5414ad,
+0x6ae290, 0x89d988, 0x50722c, 0xbea404, 0x940777, 0x7030f3,
+0x27fc00, 0xa871ea, 0x49c266, 0x3de064, 0x83dd97, 0x973fa3,
+0xfd9443, 0x8c860d, 0xde4131, 0x9d3992, 0x8c70dd, 0xe7b717,
+0x3bdf08, 0x2b3715, 0xa0805c, 0x93805a, 0x921110, 0xd8e80f,
+0xaf806c, 0x4bffdb, 0x0f9038, 0x761859, 0x15a562, 0xbbcb61,
+0xb989c7, 0xbd4010, 0x04f2d2, 0x277549, 0xf6b6eb, 0xbb22db,
+0xaa140a, 0x2f2689, 0x768364, 0x333b09, 0x1a940e, 0xaa3a51,
+0xc2a31d, 0xaeedaf, 0x12265c, 0x4dc26d, 0x9c7a2d, 0x9756c0,
+0x833f03, 0xf6f009, 0x8c402b, 0x99316d, 0x07b439, 0x15200c,
+0x5bc3d8, 0xc492f5, 0x4badc6, 0xa5ca4e, 0xcd37a7, 0x36a9e6,
+0x9492ab, 0x6842dd, 0xde6319, 0xef8c76, 0x528b68, 0x37dbfc,
+0xaba1ae, 0x3115df, 0xa1ae00, 0xdafb0c, 0x664d64, 0xb705ed,
+0x306529, 0xbf5657, 0x3aff47, 0xb9f96a, 0xf3be75, 0xdf9328,
+0x3080ab, 0xf68c66, 0x15cb04, 0x0622fa, 0x1de4d9, 0xa4b33d,
+0x8f1b57, 0x09cd36, 0xe9424e, 0xa4be13, 0xb52333, 0x1aaaf0,
+0xa8654f, 0xa5c1d2, 0x0f3f0b, 0xcd785b, 0x76f923, 0x048b7b,
+0x721789, 0x53a6c6, 0xe26e6f, 0x00ebef, 0x584a9b, 0xb7dac4,
+0xba66aa, 0xcfcf76, 0x1d02d1, 0x2df1b1, 0xc1998c, 0x77adc3,
+0xda4886, 0xa05df7, 0xf480c6, 0x2ff0ac, 0x9aecdd, 0xbc5c3f,
+0x6dded0, 0x1fc790, 0xb6db2a, 0x3a25a3, 0x9aaf00, 0x9353ad,
+0x0457b6, 0xb42d29, 0x7e804b, 0xa707da, 0x0eaa76, 0xa1597b,
+0x2a1216, 0x2db7dc, 0xfde5fa, 0xfedb89, 0xfdbe89, 0x6c76e4,
+0xfca906, 0x70803e, 0x156e85, 0xff87fd, 0x073e28, 0x336761,
+0x86182a, 0xeabd4d, 0xafe7b3, 0x6e6d8f, 0x396795, 0x5bbf31,
+0x48d784, 0x16df30, 0x432dc7, 0x356125, 0xce70c9, 0xb8cb30,
+0xfd6cbf, 0xa200a4, 0xe46c05, 0xa0dd5a, 0x476f21, 0xd21262,
+0x845cb9, 0x496170, 0xe0566b, 0x015299, 0x375550, 0xb7d51e,
+0xc4f133, 0x5f6e13, 0xe4305d, 0xa92e85, 0xc3b21d, 0x3632a1,
+0xa4b708, 0xd4b1ea, 0x21f716, 0xe4698f, 0x77ff27, 0x80030c,
+0x2d408d, 0xa0cd4f, 0x99a520, 0xd3a2b3, 0x0a5d2f, 0x42f9b4,
+0xcbda11, 0xd0be7d, 0xc1db9b, 0xbd17ab, 0x81a2ca, 0x5c6a08,
+0x17552e, 0x550027, 0xf0147f, 0x8607e1, 0x640b14, 0x8d4196,
+0xdebe87, 0x2afdda, 0xb6256b, 0x34897b, 0xfef305, 0x9ebfb9,
+0x4f6a68, 0xa82a4a, 0x5ac44f, 0xbcf82d, 0x985ad7, 0x95c7f4,
+0x8d4d0d, 0xa63a20, 0x5f57a4, 0xb13f14, 0x953880, 0x0120cc,
+0x86dd71, 0xb6dec9, 0xf560bf, 0x11654d, 0x6b0701, 0xacb08c,
+0xd0c0b2, 0x485551, 0x0efb1e, 0xc37295, 0x3b06a3, 0x3540c0,
+0x7bdc06, 0xcc45e0, 0xfa294e, 0xc8cad6, 0x41f3e8, 0xde647c,
+0xd8649b, 0x31bed9, 0xc397a4, 0xd45877, 0xc5e369, 0x13daf0,
+0x3c3aba, 0x461846, 0x5f7555, 0xf5bdd2, 0xc6926e, 0x5d2eac,
+0xed440e, 0x423e1c, 0x87c461, 0xe9fd29, 0xf3d6e7, 0xca7c22,
+0x35916f, 0xc5e008, 0x8dd7ff, 0xe26a6e, 0xc6fdb0, 0xc10893,
+0x745d7c, 0xb2ad6b, 0x9d6ecd, 0x7b723e, 0x6a11c6, 0xa9cff7,
+0xdf7329, 0xbac9b5, 0x5100b7, 0x0db2e2, 0x24ba74, 0x607de5,
+0x8ad874, 0x2c150d, 0x0c1881, 0x94667e, 0x162901, 0x767a9f,
+0xbefdfd, 0xef4556, 0x367ed9, 0x13d9ec, 0xb9ba8b, 0xfc97c4,
+0x27a831, 0xc36ef1, 0x36c594, 0x56a8d8, 0xb5a8b4, 0x0ecccf,
+0x2d8912, 0x34576f, 0x89562c, 0xe3ce99, 0xb920d6, 0xaa5e6b,
+0x9c2a3e, 0xcc5f11, 0x4a0bfd, 0xfbf4e1, 0x6d3b8e, 0x2c86e2,
+0x84d4e9, 0xa9b4fc, 0xd1eeef, 0xc9352e, 0x61392f, 0x442138,
+0xc8d91b, 0x0afc81, 0x6a4afb, 0xd81c2f, 0x84b453, 0x8c994e,
+0xcc2254, 0xdc552a, 0xd6c6c0, 0x96190b, 0xb8701a, 0x649569,
+0x605a26, 0xee523f, 0x0f117f, 0x11b5f4, 0xf5cbfc, 0x2dbc34,
+0xeebc34, 0xcc5de8, 0x605edd, 0x9b8e67, 0xef3392, 0xb817c9,
+0x9b5861, 0xbc57e1, 0xc68351, 0x103ed8, 0x4871dd, 0xdd1c2d,
+0xa118af, 0x462c21, 0xd7f359, 0x987ad9, 0xc0549e, 0xfa864f,
+0xfc0656, 0xae79e5, 0x362289, 0x22ad38, 0xdc9367, 0xaae855,
+0x382682, 0x9be7ca, 0xa40d51, 0xb13399, 0x0ed7a9, 0x480569,
+0xf0b265, 0xa7887f, 0x974c88, 0x36d1f9, 0xb39221, 0x4a827b,
+0x21cf98, 0xdc9f40, 0x5547dc, 0x3a74e1, 0x42eb67, 0xdf9dfe,
+0x5fd45e, 0xa4677b, 0x7aacba, 0xa2f655, 0x23882b, 0x55ba41,
+0x086e59, 0x862a21, 0x834739, 0xe6e389, 0xd49ee5, 0x40fb49,
+0xe956ff, 0xca0f1c, 0x8a59c5, 0x2bfa94, 0xc5c1d3, 0xcfc50f,
+0xae5adb, 0x86c547, 0x624385, 0x3b8621, 0x94792c, 0x876110,
+0x7b4c2a, 0x1a2c80, 0x12bf43, 0x902688, 0x893c78, 0xe4c4a8,
+0x7bdbe5, 0xc23ac4, 0xeaf426, 0x8a67f7, 0xbf920d, 0x2ba365,
+0xb1933d, 0x0b7cbd, 0xdc51a4, 0x63dd27, 0xdde169, 0x19949a,
+0x9529a8, 0x28ce68, 0xb4ed09, 0x209f44, 0xca984e, 0x638270,
+0x237c7e, 0x32b90f, 0x8ef5a7, 0xe75614, 0x08f121, 0x2a9db5,
+0x4d7e6f, 0x5119a5, 0xabf9b5, 0xd6df82, 0x61dd96, 0x023616,
+0x9f3ac4, 0xa1a283, 0x6ded72, 0x7a8d39, 0xa9b882, 0x5c326b,
+0x5b2746, 0xed3400, 0x7700d2, 0x55f4fc, 0x4d5901, 0x8071e0,
+0xe13f89, 0xb295f3, 0x64a8f1, 0xaea74b, 0x38fc4c, 0xeab2bb,
+0x47270b, 0xabc3a7, 0x34ba60, 0x52dd34, 0xf8563a, 0xeb7e8a,
+0x31bb36, 0x5895b7, 0x47f7a9, 0x94c3aa, 0xd39225, 0x1e7f3e,
+0xd8974e, 0xbba94f, 0xd8ae01, 0xe661b4, 0x393d8e, 0xa523aa,
+0x33068e, 0x1633b5, 0x3bb188, 0x1d3a9d, 0x4013d0, 0xcc1be5,
+0xf862e7, 0x3bf28f, 0x39b5bf, 0x0bc235, 0x22747e, 0xa247c0,
+0xd52d1f, 0x19add3, 0x9094df, 0x9311d0, 0xb42b25, 0x496db2,
+0xe264b2, 0x5ef135, 0x3bc6a4, 0x1a4ad0, 0xaac92e, 0x64e886,
+0x573091, 0x982cfb, 0x311b1a, 0x08728b, 0xbdcee1, 0x60e142,
+0xeb641d, 0xd0bba3, 0xe559d4, 0x597b8c, 0x2a4483, 0xf332ba,
+0xf84867, 0x2c8d1b, 0x2fa9b0, 0x50f3dd, 0xf9f573, 0xdb61b4,
+0xfe233e, 0x6c41a6, 0xeea318, 0x775a26, 0xbc5e5c, 0xcea708,
+0x94dc57, 0xe20196, 0xf1e839, 0xbe4851, 0x5d2d2f, 0x4e9555,
+0xd96ec2, 0xe7d755, 0x6304e0, 0xc02e0e, 0xfc40a0, 0xbbf9b3,
+0x7125a7, 0x222dfb, 0xf619d8, 0x838c1c, 0x6619e6, 0xb20d55,
+0xbb5137, 0x79e809, 0xaf9149, 0x0d73de, 0x0b0da5, 0xce7f58,
+0xac1934, 0x724667, 0x7a1a13, 0x9e26bc, 0x4555e7, 0x585cb5,
+0x711d14, 0x486991, 0x480d60, 0x56adab, 0xd62f64, 0x96ee0c,
+0x212ff3, 0x5d6d88, 0xa67684, 0x95651e, 0xab9e0a, 0x4ddefe,
+0x571010, 0x836a39, 0xf8ea31, 0x9e381d, 0xeac8b1, 0xcac96b,
+0x37f21e, 0xd505e9, 0x984743, 0x9fc56c, 0x0331b7, 0x3b8bf8,
+0x86e56a, 0x8dc343, 0x6230e7, 0x93cfd5, 0x6a8f2d, 0x733005,
+0x1af021, 0xa09fcb, 0x7415a1, 0xd56b23, 0x6ff725, 0x2f4bc7,
+0xb8a591, 0x7fac59, 0x5c55de, 0x212c38, 0xb13296, 0x5cff50,
+0x366262, 0xfa7b16, 0xf4d9a6, 0x2acfe7, 0xf07403, 0xd4d604,
+0x6fd916, 0x31b1bf, 0xcbb450, 0x5bd7c8, 0x0ce194, 0x6bd643,
+0x4fd91c, 0xdf4543, 0x5f3453, 0xe2b5aa, 0xc9aec8, 0x131485,
+0xf9d2bf, 0xbadb9e, 0x76f5b9, 0xaf15cf, 0xca3182, 0x14b56d,
+0xe9fe4d, 0x50fc35, 0xf5aed5, 0xa2d0c1, 0xc96057, 0x192eb6,
+0xe91d92, 0x07d144, 0xaea3c6, 0x343566, 0x26d5b4, 0x3161e2,
+0x37f1a2, 0x209eff, 0x958e23, 0x493798, 0x35f4a6, 0x4bdc02,
+0xc2be13, 0xbe80a0, 0x0b72a3, 0x115c5f, 0x1e1bd1, 0x0db4d3,
+0x869e85, 0x96976b, 0x2ac91f, 0x8a26c2, 0x3070f0, 0x041412,
+0xfc9fa5, 0xf72a38, 0x9c6878, 0xe2aa76, 0x50cfe1, 0x559274,
+0x934e38, 0x0a92f7, 0x5533f0, 0xa63db4, 0x399971, 0xe2b755,
+0xa98a7c, 0x008f19, 0xac54d2, 0x2ea0b4, 0xf5f3e0, 0x60c849,
+0xffd269, 0xae52ce, 0x7a5fdd, 0xe9ce06, 0xfb0ae8, 0xa50cce,
+0xea9d3e, 0x3766dd, 0xb834f5, 0x0da090, 0x846f88, 0x4ae3d5,
+0x099a03, 0x2eae2d, 0xfcb40a, 0xfb9b33, 0xe281dd, 0x1b16ba,
+0xd8c0af, 0xd96b97, 0xb52dc9, 0x9c277f, 0x5951d5, 0x21ccd6,
+0xb6496b, 0x584562, 0xb3baf2, 0xa1a5c4, 0x7ca2cf, 0xa9b93d,
+0x7b7b89, 0x483d38,
 };
 
 static const __float128 c[] = {
-/* 93 bits of pi/2 */
+/* 113 bits of pi/2 */
 #define PI_2_1 c[0]
- 1.57079632679489661923132169155131424e+00Q, /* 3fff921fb54442d18469898cc5100000 */
+ 0x1.921fb54442d18469898cc51701b8p+0Q,
 
 /* pi/2 - PI_2_1 */
 #define PI_2_1t c[1]
- 8.84372056613570112025531863263659260e-29Q, /* 3fa1c06e0e68948127044533e63a0106 */
+ 0x3.9a252049c1114cf98e804177d4c8p-116Q,
 };
 
 
@@ -525,8 +525,8 @@ __quadmath_rem_pio2q (__float128 x, __float128 *y)
   if (ix < 0x40002d97c7f3321dLL)	/* |x| in <pi/4, 3pi/4) */
     {
       if (hx > 0)
-	{ 
-	  /* 113 + 93 bit PI is ok */
+	{
+	  /* 113 + 113 bit PI is ok */
 	  z = x - PI_2_1;
 	  y[0] = z - PI_2_1t;
 	  y[1] = (z - y[0]) - PI_2_1t;
@@ -534,7 +534,7 @@ __quadmath_rem_pio2q (__float128 x, __float128 *y)
 	}
       else
         {
-	  /* 113 + 93 bit PI is ok */
+	  /* 113 + 113 bit PI is ok */
 	  z = x + PI_2_1;
 	  y[0] = z + PI_2_1t;
 	  y[1] = (z - y[0]) + PI_2_1t;
diff --git a/libquadmath/math/remquoq.c b/libquadmath/math/remquoq.c
index f7001afc3e5..fa85c5b2d50 100644
--- a/libquadmath/math/remquoq.c
+++ b/libquadmath/math/remquoq.c
@@ -1,5 +1,5 @@
 /* Compute remainder and a congruent to the quotient.
-   Copyright (C) 1997, 1999, 2002 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
 		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
@@ -49,7 +49,7 @@ remquoq (__float128 x, __float128 y, int *quo)
 
   if (hy <= 0x7ffbffffffffffffLL)
     x = fmodq (x, 8 * y);              /* now x < 8y */
-      
+
   if (((hx - hy) | (lx - ly)) == 0)
     {
       *quo = qs ? -1 : 1;
@@ -60,12 +60,12 @@ remquoq (__float128 x, __float128 y, int *quo)
   y  = fabsq (y);
   cquo = 0;
 
-  if (x >= 4 * y)
+  if (hy <= 0x7ffcffffffffffffLL && x >= 4 * y)
     {
       x -= 4 * y;
       cquo += 4;
     }
-  if (x >= 2 * y)
+  if (hy <= 0x7ffdffffffffffffLL && x >= 2 * y)
     {
       x -= 2 * y;
       cquo += 2;
@@ -101,6 +101,9 @@ remquoq (__float128 x, __float128 y, int *quo)
 
   *quo = qs ? -cquo : cquo;
 
+  /* Ensure correct sign of zero result in round-downward mode.  */
+  if (x == 0.0Q)
+    x = 0.0Q;
   if (sx)
     x = -x;
   return x;
diff --git a/libquadmath/math/rintq.c b/libquadmath/math/rintq.c
index 15d4c781ba4..7bc9684e5ea 100644
--- a/libquadmath/math/rintq.c
+++ b/libquadmath/math/rintq.c
@@ -35,7 +35,7 @@ __float128
 rintq (__float128 x)
 {
 	int64_t i0,j0,sx;
-	uint64_t i1;
+	uint64_t i1 __attribute__ ((unused));
 	__float128 w,t;
 	GET_FLT128_WORDS64(i0,i1,x);
 	sx = (((uint64_t)i0)>>63);
diff --git a/libquadmath/math/roundq.c b/libquadmath/math/roundq.c
index 7c9d640e933..b48366315f4 100644
--- a/libquadmath/math/roundq.c
+++ b/libquadmath/math/roundq.c
@@ -1,5 +1,5 @@
 /* Round __float128 to integer away from zero.
-   Copyright (C) 1997, 1999 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
 		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
@@ -21,9 +21,6 @@
 
 #include "quadmath-imp.h"
 
-static const __float128 huge = 1.0E4930Q;
-
-
 __float128
 roundq (__float128 x)
 {
@@ -32,17 +29,14 @@ roundq (__float128 x)
 
   GET_FLT128_WORDS64 (i0, i1, x);
   j0 = ((i0 >> 48) & 0x7fff) - 0x3fff;
-  if (j0 < 31)
+  if (j0 < 48)
     {
       if (j0 < 0)
 	{
-	  if (huge + x > 0.0)
-	    {
-	      i0 &= 0x8000000000000000ULL;
-	      if (j0 == -1)
-		i0 |= 0x3fff000000000000LL;
-	      i1 = 0;
-	    }
+	  i0 &= 0x8000000000000000ULL;
+	  if (j0 == -1)
+	    i0 |= 0x3fff000000000000LL;
+	  i1 = 0;
 	}
       else
 	{
@@ -50,13 +44,9 @@ roundq (__float128 x)
 	  if (((i0 & i) | i1) == 0)
 	    /* X is integral.  */
 	    return x;
-	  if (huge + x > 0.0)
-	    {
-	      /* Raise inexact if x != 0.  */
-	      i0 += 0x0000800000000000LL >> j0;
-	      i0 &= ~i;
-	      i1 = 0;
-	    }
+	  i0 += 0x0000800000000000LL >> j0;
+	  i0 &= ~i;
+	  i1 = 0;
 	}
     }
   else if (j0 > 111)
@@ -74,14 +64,10 @@ roundq (__float128 x)
 	/* X is integral.  */
 	return x;
 
-      if (huge + x > 0.0)
-	{
-	  /* Raise inexact if x != 0.  */
-	  uint64_t j = i1 + (1LL << (111 - j0));
-	  if (j < i1)
-	    i0 += 1;
-	  i1 = j;
-	}
+      uint64_t j = i1 + (1LL << (111 - j0));
+      if (j < i1)
+	i0 += 1;
+      i1 = j;
       i1 &= ~i;
     }
 
diff --git a/libquadmath/math/scalblnq.c b/libquadmath/math/scalblnq.c
index a414045b51f..6aac21125f8 100644
--- a/libquadmath/math/scalblnq.c
+++ b/libquadmath/math/scalblnq.c
@@ -1,7 +1,7 @@
 /* scalblnq.c -- __float128 version of s_scalbn.c.
  * Conversion to IEEE quad long double by Jakub Jelinek, jj@ultra.linux.cz.
  */
-   
+
 /*
  * ====================================================
  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
diff --git a/libquadmath/math/scalbnq.c b/libquadmath/math/scalbnq.c
index 9975a47154c..9ed5fe61799 100644
--- a/libquadmath/math/scalbnq.c
+++ b/libquadmath/math/scalbnq.c
@@ -1,7 +1,7 @@
 /* scalbnq.c -- __float128 version of s_scalbn.c.
  * Conversion to IEEE quad long double by Jakub Jelinek, jj@ultra.linux.cz.
  */
-   
+
 /*
  * ====================================================
  * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
diff --git a/libquadmath/math/sincos_table.c b/libquadmath/math/sincos_table.c
index b7e7c750314..362f517a6ef 100644
--- a/libquadmath/math/sincos_table.c
+++ b/libquadmath/math/sincos_table.c
@@ -1,5 +1,5 @@
 /* Quad-precision floating point sine and cosine tables.
-   Copyright (C) 1999 Free Software Foundation, Inc.
+   Copyright (C) 1999-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Jakub Jelinek <jj@ultra.linux.cz>
 
@@ -37,7 +37,7 @@ const __float128 __sincosq_table[] = {
 /* sin(x) = 0.25dc50bc95711d0d9787d108fd438cf5959ee0bfb7a1e36e8b1a112968f356657420e9cc9ea */
  1.47892995873409608580026675734609314e-01Q, /* 3ffc2ee285e4ab88e86cbc3e8847ea1c */
  9.74950446464233268291647449768590886e-36Q, /* 3f8a9eb2b3dc17f6f43c6dd16342252d */
-    
+
 /* x = 1.56250000000000000000000000000000000e-01 3ffc4000000000000000000000000000 */
 /* cos(x) = 0.fce1a053e621438b6d60c76e8c45bf0a9dc71aa16f922acc10e95144ec796a249813c9cb649 */
  9.87817783816471944100503034363211317e-01Q, /* 3ffef9c340a7cc428716dac18edd188b */
diff --git a/libquadmath/math/sincosq.c b/libquadmath/math/sincosq.c
index b7c221486d7..d83b1a6b757 100644
--- a/libquadmath/math/sincosq.c
+++ b/libquadmath/math/sincosq.c
@@ -1,5 +1,5 @@
 /* Compute sine and cosine of argument.
-   Copyright (C) 1997, 1999 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
 		  Jakub Jelinek <jj@ultra.linux.cz>.
@@ -19,6 +19,7 @@
    Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
    02111-1307 USA.  */
 
+#include <errno.h>
 #include "quadmath-imp.h"
 
 void
@@ -37,6 +38,8 @@ sincosq (__float128 x, __float128 *sinx, __float128 *cosx)
     {
       /* sin(Inf or NaN) is NaN */
       *sinx = *cosx = x - x;
+      if (isinfq (x))
+	errno = EDOM;
     }
   else
     {
diff --git a/libquadmath/math/sincosq_kernel.c b/libquadmath/math/sincosq_kernel.c
index f6341a4d948..171ed6fc117 100644
--- a/libquadmath/math/sincosq_kernel.c
+++ b/libquadmath/math/sincosq_kernel.c
@@ -1,5 +1,5 @@
 /* Quad-precision floating point sine and cosine on <-pi/4,pi/4>.
-   Copyright (C) 1999 Free Software Foundation, Inc.
+   Copyright (C) 1999-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Jakub Jelinek <jj@ultra.linux.cz>
 
@@ -110,12 +110,15 @@ __quadmath_kernel_sincosq(__float128 x, __float128 y, __float128 *sinx,
       /* Argument is small enough to approximate it by a Chebyshev
 	 polynomial of degree 16(17).  */
       if (tix < 0x3fc60000)		/* |x| < 2^-57 */
-	if (!((int)x))			/* generate inexact */
-	  {
-	    *sinx = x;
-	    *cosx = ONE;
-	    return;
-	  }
+	{
+	  math_check_force_underflow (x);
+	  if (!((int)x))			/* generate inexact */
+	    {
+	      *sinx = x;
+	      *cosx = ONE;
+	      return;
+	    }
+	}
       z = x * x;
       *sinx = x + (x * (z*(SIN1+z*(SIN2+z*(SIN3+z*(SIN4+
 			z*(SIN5+z*(SIN6+z*(SIN7+z*SIN8)))))))));
diff --git a/libquadmath/math/sinhq.c b/libquadmath/math/sinhq.c
index eff8149b515..a4191c0fd3d 100644
--- a/libquadmath/math/sinhq.c
+++ b/libquadmath/math/sinhq.c
@@ -85,8 +85,11 @@ sinhq (__float128 x)
   if (ix <= 0x40044000)
     {
       if (ix < 0x3fc60000) /* |x| < 2^-57 */
-	if (shuge + x > one)
-	  return x;		/* sinh(tiny) = tiny with inexact */
+	{
+	  math_check_force_underflow (x);
+	  if (shuge + x > one)
+	    return x;		/* sinh(tiny) = tiny with inexact */
+	}
       t = expm1q (u.value);
       if (ix < 0x3fff0000)
 	return h * (2.0Q * t - t * t / (t + one));
diff --git a/libquadmath/math/sinq_kernel.c b/libquadmath/math/sinq_kernel.c
index 86034551d43..873341a7617 100644
--- a/libquadmath/math/sinq_kernel.c
+++ b/libquadmath/math/sinq_kernel.c
@@ -90,7 +90,10 @@ __quadmath_kernel_sinq (__float128 x, __float128 y, int iy)
       /* Argument is small enough to approximate it by a Chebyshev
 	 polynomial of degree 17.  */
       if (tix < 0x3fc60000)		/* |x| < 2^-57 */
-	if (!((int)x)) return x;	/* generate inexact */
+	{
+	  math_check_force_underflow (x);
+	  if (!((int)x)) return x;	/* generate inexact */
+	}
       z = x * x;
       return x + (x * (z*(SIN1+z*(SIN2+z*(SIN3+z*(SIN4+
 		       z*(SIN5+z*(SIN6+z*(SIN7+z*SIN8)))))))));
diff --git a/libquadmath/math/tanhq.c b/libquadmath/math/tanhq.c
index 4ef4fd02143..6b60dbf1b2f 100644
--- a/libquadmath/math/tanhq.c
+++ b/libquadmath/math/tanhq.c
@@ -72,7 +72,10 @@ tanhq (__float128 x)
       if (u.value == 0)
 	return x;		/* x == +- 0 */
       if (ix < 0x3fc60000)	/* |x| < 2^-57 */
-	return x * (one + tiny); /* tanh(small) = small */
+	{
+	  math_check_force_underflow (x);
+	  return x * (one + tiny); /* tanh(small) = small */
+	}
       u.words32.w0 = ix;	/* Absolute value of x.  */
       if (ix >= 0x3fff0000)
 	{			/* |x| >= 1  */
diff --git a/libquadmath/math/tanq.c b/libquadmath/math/tanq.c
index 690d94b782c..d2864f61c7d 100644
--- a/libquadmath/math/tanq.c
+++ b/libquadmath/math/tanq.c
@@ -12,9 +12,9 @@
 /*
   Long double expansions are
   Copyright (C) 2001 Stephen L. Moshier <moshier@na-net.ornl.gov>
-  and are incorporated herein by permission of the author.  The author 
+  and are incorporated herein by permission of the author.  The author
   reserves the right to distribute this material elsewhere under different
-  copying permissions.  These modifications are distributed here under 
+  copying permissions.  These modifications are distributed here under
   the following terms:
 
     This library is free software; you can redistribute it and/or
@@ -99,8 +99,13 @@ __quadmath_kernel_tanq (__float128 x, __float128 y, int iy)
 	  if ((ix | u.words32.w1 | u.words32.w2 | u.words32.w3
 	       | (iy + 1)) == 0)
 	    return one / fabsq (x);
+	  else if (iy == 1)
+	    {
+	      math_check_force_underflow (x);
+	      return x;
+	    }
 	  else
-	    return (iy == 1) ? x : -one / x;
+	    return -one / x;
 	}
     }
   if (ix >= 0x3ffe5942) /* |x| >= 0.6743316650390625 */
@@ -163,7 +168,7 @@ __quadmath_kernel_tanq (__float128 x, __float128 y, int iy)
 /* tanq.c -- __float128 version of s_tan.c.
  * Conversion to IEEE quad long double by Jakub Jelinek, jj@ultra.linux.cz.
  */
-  
+
 /* @(#)s_tan.c 5.1 93/09/24 */
 /*
  * ====================================================
diff --git a/libquadmath/math/truncq.c b/libquadmath/math/truncq.c
index 608871785d1..607c24985e3 100644
--- a/libquadmath/math/truncq.c
+++ b/libquadmath/math/truncq.c
@@ -1,8 +1,8 @@
 /* Truncate argument to nearest integral value not larger than the argument.
-   Copyright (C) 1997, 1999 Free Software Foundation, Inc.
+   Copyright (C) 1997-2017 Free Software Foundation, Inc.
    This file is part of the GNU C Library.
    Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997 and
-   		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
+		  Jakub Jelinek <jj@ultra.linux.cz>, 1999.
 
    The GNU C Library is free software; you can redistribute it and/or
    modify it under the terms of the GNU Lesser General Public
diff --git a/libquadmath/quadmath-imp.h b/libquadmath/quadmath-imp.h
index 40b346b6ff2..ed186144e75 100644
--- a/libquadmath/quadmath-imp.h
+++ b/libquadmath/quadmath-imp.h
@@ -186,4 +186,45 @@ do {                                   \
   __builtin_fpclassify (QUADFP_NAN, QUADFP_INFINITE, QUADFP_NORMAL, \
 			QUADFP_SUBNORMAL, QUADFP_ZERO, x)
 
+#ifndef math_opt_barrier
+# define math_opt_barrier(x) \
+({ __typeof (x) __x = (x); __asm ("" : "+m" (__x)); __x; })
+# define math_force_eval(x) \
+({ __typeof (x) __x = (x); __asm __volatile__ ("" : : "m" (__x)); })
+#endif
+
+/* math_narrow_eval reduces its floating-point argument to the range
+   and precision of its semantic type.  (The original evaluation may
+   still occur with excess range and precision, so the result may be
+   affected by double rounding.)  */
+#define math_narrow_eval(x) (x)
+
+/* If X (which is not a NaN) is subnormal, force an underflow
+   exception.  */
+#define math_check_force_underflow(x)				\
+  do								\
+    {								\
+      __float128 force_underflow_tmp = (x);			\
+      if (fabsq (force_underflow_tmp) < FLT128_MIN)		\
+	{							\
+	  __float128 force_underflow_tmp2			\
+	    = force_underflow_tmp * force_underflow_tmp;	\
+	  math_force_eval (force_underflow_tmp2);		\
+	}							\
+    }								\
+  while (0)
+/* Likewise, but X is also known to be nonnegative.  */
+#define math_check_force_underflow_nonneg(x)			\
+  do								\
+    {								\
+      __float128 force_underflow_tmp = (x);			\
+      if (force_underflow_tmp < FLT128_MIN)			\
+	{							\
+	  __float128 force_underflow_tmp2			\
+	    = force_underflow_tmp * force_underflow_tmp;	\
+	  math_force_eval (force_underflow_tmp2);		\
+	}							\
+    }								\
+  while (0)
+
 #endif
-- 
2.30.2