[Bug libfortran/24313] New: complex sqrt function does not return principal value

[Bug libfortran/24313] New: complex sqrt function does not return principal value

[mick at nag dot co dot uk]

This program demonstrates a bug in the complex sqrt function
in gfortran. When complex number x = (0.0,-1.0), the result comes
back as the negative of the expected result. Of course, the negative
of the expected result is still a square root of the original
number, but the f95 standard says that the principal square root
should be returned, i.e. the real part of the result should be positive.

      complex x, y
      x = cmplx(0.0,-1.0)
      y = sqrt(x)
      write (*,*) x, y
      x = cmplx(0.0+1.0e-38,-1.0)
      y = sqrt(x)
      write (*,*) x, y
      x = cmplx(0.0-1.0e-38,-1.0)
      y = sqrt(x)
      write (*,*) x, y
      end

The results I get are

% gfortran -v
Using built-in specs.
Target: x86_64-unknown-linux-gnu
Configured with: ../gcc/configure --prefix=/var/tmp/gfortran-20051007/irun
--enable-languages=c,f95
Thread model: posix
gcc version 4.1.0 20051007 (experimental)

% gfortran sqrt.f
% ./a.out
 (  0.000000    , -1.000000    ) (-0.7071068    , 0.7071068    )
 ( 9.9999994E-39, -1.000000    ) ( 0.7071068    ,-0.7071068    )
 (-9.9999994E-39, -1.000000    ) ( 0.7071068    ,-0.7071068    )

Using a tiny number in place of the zero real part shows the
correct results.

This also affects the double complex version of sqrt.


-- 
           Summary: complex sqrt function does not return principal value
           Product: gcc
           Version: 4.1.0
            Status: UNCONFIRMED
          Severity: normal
          Priority: P2
         Component: libfortran
        AssignedTo: unassigned at gcc dot gnu dot org
        ReportedBy: mick at nag dot co dot uk
GCC target triplet: x86_64-unknown-linux-gnu


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[pinskia at gcc dot gnu dot org]

------- Comment #1 from pinskia at gcc dot gnu dot org  2005-10-11 17:31 -------
Hmm, I think this is a bug in glibc's csqrt. as I get the following on
powerpc-darwin:
 (  0.000000    , -1.000000    ) ( 0.7071068    ,-0.7071068    )
 ( 9.9999994E-39, -1.000000    ) ( 0.7071068    ,-0.7071068    )
 (-9.9999994E-39, -1.000000    ) ( 0.7071068    ,-0.7071068    )


-- 

pinskia at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
 GCC target triplet|x86_64-unknown-linux-gnu    |x86_64-*-linux-gnu


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[pinskia at gcc dot gnu dot org]

------- Comment #2 from pinskia at gcc dot gnu dot org  2005-10-11 17:33 -------
This is a bug in glibc's csqrt as gfortran just calls it.  Please report it to
them.


-- 

pinskia at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|UNCONFIRMED                 |RESOLVED
         Resolution|                            |WORKSFORME


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[kargl at gcc dot gnu dot org]

------- Comment #3 from kargl at gcc dot gnu dot org  2005-10-11 18:48 -------
I don't think this is a glibc problem.  gfortran will use a GCC
builtin function.  See gcc/fortran/mathbuiltins.def  To confirm 
this, I get

troutmask:kargl[210] ./z
 (  0.000000    , -1.000000    ) (-0.7071068    , 0.7071068    )
 ( 1.1754944E-38, -1.000000    ) ( 0.7071068    ,-0.7071068    )
 (-1.1754944E-38, -1.000000    ) ( 0.7071068    ,-0.7071068    )

on amd64-*-freebsd, which does not have a csqrtf in libm.

Note, I've replaced Mick's subnormal number with TINY(1.).


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[kargl at gcc dot gnu dot org]

------- Comment #4 from kargl at gcc dot gnu dot org  2005-10-11 19:56 -------
This is definitely a gcc bug.

gfortran -fdump-tree-original csqrt.f 

yields (with removal of IO crap and shorting of long constants)

MAIN__ ()
{
  complex4 y;
  complex4 x;

  x = __complex__ (0.0, -1.0e+0);
  y = __builtin_csqrtf (x);

 x = __complex__ (1.17549435082228750796873653722224e-38, -1.0e+0);
  y = __builtin_csqrtf (x);

  x = __complex__ (-1.175494350822287507968736537222e-38, -1.0e+0);
  y = __builtin_csqrtf (x);

}


Note the use of __builtin_csqrtf().


-- 

kargl at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|RESOLVED                    |UNCONFIRMED
         Resolution|WORKSFORME                  |


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[pinskia at gcc dot gnu dot org]

------- Comment #5 from pinskia at gcc dot gnu dot org  2005-10-11 19:58 -------
(In reply to comment #4)
> Note the use of __builtin_csqrtf().

But that just calls into csqrt which is part of libm.  Look at the asm output.


-- 

pinskia at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
             Status|UNCONFIRMED                 |RESOLVED
         Resolution|                            |INVALID


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[kargl at gcc dot gnu dot org]

------- Comment #6 from kargl at gcc dot gnu dot org  2005-10-11 20:11 -------
FreeBSD *does not* have a csqrt or csqrtf in libm!  The problem is
in intrinsics/c99_functions.c.  I'm testing a patch.


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[pinskia at gcc dot gnu dot org]

------- Comment #7 from pinskia at gcc dot gnu dot org  2005-10-11 20:18 -------
(In reply to comment #6)
> FreeBSD *does not* have a csqrt or csqrtf in libm!  The problem is
> in intrinsics/c99_functions.c.  I'm testing a patch.

But if csqrtf in c99_functions.c was pulled from glibc, then there is still a
bug in glibc as this was reported with the target being x86_64-linux-gnu.


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[kargl at gcc dot gnu dot org]

------- Comment #8 from kargl at gcc dot gnu dot org  2005-10-11 20:21 -------
Yes, I agree glibc has the bug.  But, libgfortran also has the bug
and will be present on every OS that does not have a libm with csqrt
and csqrtf.  I have a patch for libgfortran.


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[kargl at gcc dot gnu dot org]

------- Comment #9 from kargl at gcc dot gnu dot org  2005-10-11 20:46 -------
Patch for libgfortran 

http://gcc.gnu.org/ml/gcc-patches/2005-10/msg00626.html


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[cvs-commit at gcc dot gnu dot org]

------- Comment #10 from cvs-commit at gcc dot gnu dot org  2005-10-11 23:35 -------
Subject: Bug 24313

CVSROOT:        /cvs/gcc
Module name:    gcc
Changes by:     kargl@gcc.gnu.org       2005-10-11 23:35:27

Modified files:
        libgfortran    : ChangeLog 
        libgfortran/intrinsics: c99_functions.c 
        gcc/testsuite  : ChangeLog 
Added files:
        gcc/testsuite/gfortran.dg: csqrt_2.f 

Log message:
        PR libgfortran/24313
        * c99_functions.c (csqrtf, csqrt): Fix choice of branch cut.  Note
        csqrt{f} were imported from glibc, and this bug is still present
        there. glibc PR is 1146.

Patches:
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/libgfortran/ChangeLog.diff?cvsroot=gcc&r1=1.320&r2=1.321
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/libgfortran/intrinsics/c99_functions.c.diff?cvsroot=gcc&r1=1.16&r2=1.17
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/testsuite/gfortran.dg/csqrt_2.f.diff?cvsroot=gcc&r1=NONE&r2=1.1
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/testsuite/ChangeLog.diff?cvsroot=gcc&r1=1.6170&r2=1.6171


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[pinskia at gcc dot gnu dot org]

------- Comment #11 from pinskia at gcc dot gnu dot org  2005-11-24 16:41 -------
*** Bug 25017 has been marked as a duplicate of this bug. ***


-- 

pinskia at gcc dot gnu dot org changed:

           What    |Removed                     |Added
----------------------------------------------------------------------------
                 CC|                            |harald dot vogt at desy dot
                   |                            |de


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[harald dot vogt at desy dot de]

------- Comment #12 from harald dot vogt at desy dot de  2006-01-25 21:07 -------
(In reply to comment #9)
> Patch for libgfortran 
> 
> http://gcc.gnu.org/ml/gcc-patches/2005-10/msg00626.html
> 
It is a useful strategy to use libm but...

It seems that glibc has still problems with csqrt. See also
http://sources.redhat.com/bugzilla/show_bug.cgi?id=2181
http://sources.redhat.com/bugzilla/show_bug.cgi?id=2182 .
Therefore it may be useful to check in intrinsics/configure if libm is broken
in the architecture used, and if so use the code in intrinsics/c99_functions.c
.
A possible solution for csqrt is shown in
http://www-zeuthen.desy.de/~hvogt/gfortran/acinclude.m4.diff .
It is better to be independend from glibc coming with the OS. This enables to
install gcc4(gfortran) without patching libm if it has bugs.


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[burnus at gcc dot gnu dot org]

------- Comment #13 from burnus at gcc dot gnu dot org  2006-10-29 11:57 -------
Subject: Bug 24313

Author: burnus
Date: Sun Oct 29 11:56:56 2006
New Revision: 118142

URL: http://gcc.gnu.org/viewcvs?root=gcc&view=rev&rev=118142
Log:
2006-10-28  Tobias Burnus  <burnus@net-b.de>

       PR libgfortran/24313
       * gfortran.dg/csqrt_2.f: Remove xfail *-*-linux-gnu.


Modified:
    trunk/gcc/testsuite/ChangeLog
    trunk/gcc/testsuite/gfortran.dg/csqrt_2.f


-- 


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

[Bug libfortran/24313] complex sqrt function does not return principal value

[burnus at gcc dot gnu dot org]

------- Comment #14 from burnus at gcc dot gnu dot org  2006-10-29 11:58 -------
For completeness, the two problems are in the meanwhile fixed in glibc:

- csqrt bug, fixed 2005-10-13,
http://sources.redhat.com/bugzilla/show_bug.cgi?id=1466

- cacosh bug, fixed 2006-08-03,
http://sources.redhat.com/bugzilla/show_bug.cgi?id=2182


-- 


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