Daily bump.
[official-gcc.git] / libquadmath / math / sqrtq.c
blob56ea5d3243c06df2276689849b448eebacae6367
1 #include "quadmath-imp.h"
2 #include <math.h>
3 #include <float.h>
5 __float128
6 sqrtq (const __float128 x)
8 __float128 y;
9 int exp;
11 if (isnanq (x) || (isinfq (x) && x > 0))
12 return x;
14 if (x == 0)
15 return x;
17 if (x < 0)
19 /* Return NaN with invalid signal. */
20 return (x - x) / (x - x);
23 if (x <= DBL_MAX && x >= DBL_MIN)
25 /* Use double result as starting point. */
26 y = sqrt ((double) x);
28 /* Two Newton iterations. */
29 y -= 0.5q * (y - x / y);
30 y -= 0.5q * (y - x / y);
31 return y;
34 #ifdef HAVE_SQRTL
36 long double xl = (long double) x;
37 if (xl <= LDBL_MAX && xl >= LDBL_MIN)
39 /* Use long double result as starting point. */
40 y = (__float128) sqrtl (xl);
42 /* One Newton iteration. */
43 y -= 0.5q * (y - x / y);
44 return y;
47 #endif
49 /* If we're outside of the range of C types, we have to compute
50 the initial guess the hard way. */
51 y = frexpq (x, &exp);
52 if (exp % 2)
53 y *= 2, exp--;
55 y = sqrt (y);
56 y = scalbnq (y, exp / 2);
58 /* Two Newton iterations. */
59 y -= 0.5q * (y - x / y);
60 y -= 0.5q * (y - x / y);
61 return y;