[Bug fortran/30981] a ** exp fails for integer exponents if exp is "-huge()-1" (endless loop)

[ray at ultramarine dot com]

------- Comment #4 from ray at ultramarine dot com  2007-02-27 19:52 -------
Subject: Re:  a ** exp fails for integer exponents if exp
 is "-huge()-1" (endless loop)

On Tue, 27 Feb 2007, burnus at gcc dot gnu dot org wrote:

> 
> 
> ------- Comment #3 from burnus at gcc dot gnu dot org  2007-02-27 19:35 -------
> Hi,
> 
> > > Could you post an example?
> 
> With example I mean something which actually compiles and runs. Here, I have
> two problems: 
>       include '$(where)/amos/include/essential.ecm'
> is missing as well as the main program which calls the routine. And as always:
> The smaller the better.


  I really try to do this, but I am an "old" fortran programer and 
  I do not know how to set +Infinity. Yes, of course, I have a case
  that computes it and then hangs, but my code is over 500,000 lines
  long and I did not think it is appropriate.

> In addition, please attach such long files, they clutter the bug report page
> quite a bit if they are put into the comment field.
> 
> >       if ( n_normal .ge. 10.-eps .and. n_normal .le. 10+eps ) then
> >
> >   PS. Notice the "eps" in the last segment of code. This is 
> >   really a bug also. It originally was 
> >       If ( n_normal .eq. 10. ) then
> >   But in some cases it missed the proper branch.
> 
> If you do:
>   n_normal = 10.0
> then  if(n_normal == 10.0)  should work. If you do
>   n_normal = 500.0/100.0*2.0
> then it might not work due to rounding errors and representation problems of a
> decimal number as binary number. Therefore, one should almost always use
> something like  if ( abs(n_normal - 10.0) < eps),  where eps is sensibly
> chosen, it could be e.g. eps = epsilon(n_normal), but this might be even too
> small. For irrational numbers such as 3.14159 it is more obvious that one
> should use, e.g. (abs(myNumber - 3.14159) < 1e-4) instead of  "if (myNumber ==
> 3.14159)".

  This is not as silly as it sounds. All I care about is to handle
  the special case of 10 (exactly). I have seen cases where the
  equality fails, but the number is really 10. It probably has
  something to do with how this is done. I never could get a 
  small case to fail.

> Back to the problem:
> 
> The calculation goes only into an endless loop if the exponent is
> "-huge(integer)-1".
> 
> Test case:
>       program test
>       implicit none
>       print *, 5.0**(-huge(i)-1)
>       end program
> 
> The problem is that in pow_r4_i4 the following operation is done:
> 
>       if (n < 0)
>         {
>           n = -n;
> 
> The problem is: For n == -huge(n) - 1 exists no positive number. Thus
>  n == -n == -2147483648  (for integer(4))
> 
> 
> -- 
> 
> burnus at gcc dot gnu dot org changed:
> 
>            What    |Removed                     |Added
> ----------------------------------------------------------------------------
>              Status|UNCONFIRMED                 |NEW
>      Ever Confirmed|0                           |1
>            Keywords|                            |wrong-code
>    Last reconfirmed|0000-00-00 00:00:00         |2007-02-27 19:35:04
>                date|                            |
>             Summary|Program "Hangs"             |a ** exp fails for integer
>                    |                            |exponents if exp is "-
>                    |                            |huge()-1" (endless loop)
> 
> 
> http://gcc.gnu.org/bugzilla/show_bug.cgi?id=30981
> 
> ------- You are receiving this mail because: -------
> You reported the bug, or are watching the reporter.
> 
> -------------------+  Spam Scoring Results   +------------------------
> Content analysis details:   (-2.6 hits, 5.0 required)
> -2.6 BAYES_00               BODY: Bayesian spam probability is 0 to 1%
>                             [score: 0.0033]
> 
> -------------------+ End Spam Scoring Results +-----------------------
> 


Ultramarine, Inc.
http://www.ultramarine.com
Phone: 713-975-8146
Fax:   713-975-8179


-- 


http://gcc.gnu.org/bugzilla/show_bug.cgi?id=30981