* arith.c (gfc_arith_init_1): Fix off by one problem;

  (gfc_check_integer_range): Chop extra bits in subnormal numbers.

git-svn-id: svn+ssh://gcc.gnu.org/svn/gcc/trunk@100299 138bc75d-0d04-0410-961f-82ee72b054a4

diff --git a/gcc/fortran/ChangeLog b/gcc/fortran/ChangeLog
index ed4f712..097ef92 100644 (file)
--- a/gcc/fortran/ChangeLog
+++ b/gcc/fortran/ChangeLog
@@ -1,3 +1,8 @@
+2005-05-28  Steven G. Kargl  <kargls@comcast.net>
+
+       * arith.c (gfc_arith_init_1): Fix off by one problem;
+       (gfc_check_integer_range): Chop extra bits in subnormal numbers.
+
 2005-05-28  Jerry DeLisle   <jvdelisle@verizon.net>
             Steven G. Kargl  <kargls@comcast.net>
 

diff --git a/gcc/fortran/arith.c b/gcc/fortran/arith.c
index ef19217..88b6c36 100644 (file)
--- a/gcc/fortran/arith.c
+++ b/gcc/fortran/arith.c
@@ -259,9 +259,9 @@ gfc_arith_init_1 (void)
       mpfr_init (real_info->tiny);
       mpfr_set (real_info->tiny, b, GFC_RND_MODE);
 
-      /* subnormal (x) = b**(emin - digit + 1) */
+      /* subnormal (x) = b**(emin - digit) */
       mpfr_set_ui (b, real_info->radix, GFC_RND_MODE);
-      mpfr_pow_si (b, b, real_info->min_exponent - real_info->digits + 1,
+      mpfr_pow_si (b, b, real_info->min_exponent - real_info->digits,
                   GFC_RND_MODE);
 
       mpfr_init (real_info->subnormal);
@@ -381,9 +381,42 @@ gfc_check_real_range (mpfr_t p, int kind)
   if (mpfr_sgn (q) == 0)
     retval = ARITH_OK;
   else if (mpfr_cmp (q, gfc_real_kinds[i].huge) > 0)
-      retval = ARITH_OVERFLOW;
+    retval = ARITH_OVERFLOW;
   else if (mpfr_cmp (q, gfc_real_kinds[i].subnormal) < 0)
     retval = ARITH_UNDERFLOW;
+  else if (mpfr_cmp (q, gfc_real_kinds[i].tiny) < 0)
+    {
+      /* MPFR operates on a numbers with a given precision and enormous
+       exponential range.  To represent subnormal numbers the exponent is
+       allowed to become smaller than emin, but always retains the full
+       precision.  This function resets unused bits to 0 to alleviate
+       rounding problems.  Note, a future version of MPFR will have a
+       mpfr_subnormalize() function, which handles this truncation in a
+       more efficient and robust way.  */
+
+      int j, k;
+      char *bin, *s;
+      mp_exp_t e;
+
+      bin = mpfr_get_str (NULL, &e, gfc_real_kinds[i].radix, 0, q, GMP_RNDN);
+      k = gfc_real_kinds[i].digits - (gfc_real_kinds[i].min_exponent - e);
+      for (j = k; j < gfc_real_kinds[i].digits; j++)
+       bin[j] = '0';
+      /* Need space for '0.', bin, 'E', and e */
+      s = (char *) gfc_getmem (strlen(bin)+10);
+      sprintf (s, "0.%sE%d", bin, (int) e);
+      mpfr_set_str (q, s, gfc_real_kinds[i].radix, GMP_RNDN);
+
+      if (mpfr_sgn (p) < 0)
+       mpfr_neg (p, q, GMP_RNDN);
+      else
+       mpfr_set (p, q, GMP_RNDN);
+
+      gfc_free (s);
+      gfc_free (bin);
+
+      retval = ARITH_OK;
+    }
   else
     retval = ARITH_OK;