diff options
Diffstat (limited to 'gcc/real.c')
| -rw-r--r-- | gcc/real.c | 100 |
1 files changed, 100 insertions, 0 deletions
diff --git a/gcc/real.c b/gcc/real.c index d9a4b801e3d..1598254c993 100644 --- a/gcc/real.c +++ b/gcc/real.c @@ -23,6 +23,8 @@ #include "config.h" #include "system.h" +#include "coretypes.h" +#include "tm.h" #include "tree.h" #include "toplev.h" #include "real.h" @@ -2721,6 +2723,7 @@ const struct real_format ieee_single_format = 24, -125, 128, + 31, true, true, true, @@ -2914,6 +2917,7 @@ const struct real_format ieee_double_format = 53, -1021, 1024, + 63, true, true, true, @@ -3169,6 +3173,7 @@ const struct real_format ieee_extended_motorola_format = 64, -16382, 16384, + 95, true, true, true, @@ -3185,6 +3190,7 @@ const struct real_format ieee_extended_intel_96_format = 64, -16381, 16384, + 79, true, true, true, @@ -3201,6 +3207,7 @@ const struct real_format ieee_extended_intel_128_format = 64, -16381, 16384, + 79, true, true, true, @@ -3294,6 +3301,7 @@ const struct real_format ibm_extended_format = 53 + 53, -1021, 1024, + -1, true, true, true, @@ -3547,6 +3555,7 @@ const struct real_format ieee_quad_format = 113, -16381, 16384, + 127, true, true, true, @@ -3854,6 +3863,7 @@ const struct real_format vax_f_format = 24, -127, 127, + 15, false, false, false, @@ -3870,6 +3880,7 @@ const struct real_format vax_d_format = 56, -127, 127, + 15, false, false, false, @@ -3886,6 +3897,7 @@ const struct real_format vax_g_format = 53, -1023, 1023, + 15, false, false, false, @@ -4067,6 +4079,7 @@ const struct real_format i370_single_format = 6, -64, 63, + 31, false, false, false, /* ??? The encoding does allow for "unnormals". */ @@ -4083,6 +4096,7 @@ const struct real_format i370_double_format = 14, -64, 63, + 63, false, false, false, /* ??? The encoding does allow for "unnormals". */ @@ -4297,6 +4311,7 @@ const struct real_format c4x_single_format = 24, -126, 128, + -1, false, false, false, @@ -4313,6 +4328,7 @@ const struct real_format c4x_extended_format = 32, -126, 128, + -1, false, false, false, @@ -4358,6 +4374,7 @@ const struct real_format real_internal_format = SIGNIFICAND_BITS - 2, -MAX_EXP, MAX_EXP, + -1, true, true, false, @@ -4381,3 +4398,86 @@ const struct real_format *real_format_for_mode[TFmode - QFmode + 1] = NULL, /* XFmode */ &ieee_quad_format /* TFmode */ }; + + +/* Calculate the square root of X in mode MODE, and store the result + in R. For details see "High Precision Division and Square Root", + Alan H. Karp and Peter Markstein, HP Lab Report 93-93-42, June + 1993. http://www.hpl.hp.com/techreports/93/HPL-93-42.pdf. */ + +void +real_sqrt (r, mode, x) + REAL_VALUE_TYPE *r; + enum machine_mode mode; + const REAL_VALUE_TYPE *x; +{ + static REAL_VALUE_TYPE halfthree; + static REAL_VALUE_TYPE half; + static bool init = false; + REAL_VALUE_TYPE h, t, i; + int iter, exp; + + /* sqrt(-0.0) is -0.0. */ + if (real_isnegzero (x)) + { + *r = *x; + return; + } + + /* Negative arguments return NaN. */ + if (real_isneg (x)) + { + /* Mode is ignored for canonical NaN. */ + real_nan (r, "", 1, SFmode); + return; + } + + /* Infinity and NaN return themselves. */ + if (real_isinf (x) || real_isnan (x)) + { + *r = *x; + return; + } + + if (!init) + { + real_arithmetic (&half, RDIV_EXPR, &dconst1, &dconst2); + real_arithmetic (&halfthree, PLUS_EXPR, &dconst1, &half); + init = true; + } + + /* Initial guess for reciprocal sqrt, i. */ + exp = real_exponent (x); + real_ldexp (&i, &dconst1, -exp/2); + + /* Newton's iteration for reciprocal sqrt, i. */ + for (iter = 0; iter < 16; iter++) + { + /* i(n+1) = i(n) * (1.5 - 0.5*i(n)*i(n)*x). */ + real_arithmetic (&t, MULT_EXPR, x, &i); + real_arithmetic (&h, MULT_EXPR, &t, &i); + real_arithmetic (&t, MULT_EXPR, &h, &half); + real_arithmetic (&h, MINUS_EXPR, &halfthree, &t); + real_arithmetic (&t, MULT_EXPR, &i, &h); + + /* Check for early convergence. */ + if (iter >= 6 && real_identical (&i, &t)) + break; + + /* ??? Unroll loop to avoid copying. */ + i = t; + } + + /* Final iteration: r = i*x + 0.5*i*x*(1.0 - i*(i*x)). */ + real_arithmetic (&t, MULT_EXPR, x, &i); + real_arithmetic (&h, MULT_EXPR, &t, &i); + real_arithmetic (&i, MINUS_EXPR, &dconst1, &h); + real_arithmetic (&h, MULT_EXPR, &t, &i); + real_arithmetic (&i, MULT_EXPR, &half, &h); + real_arithmetic (&h, PLUS_EXPR, &t, &i); + + /* ??? We need a Tuckerman test to get the last bit. */ + + real_convert (r, mode, &h); +} + |