Bug 18902 - Naive (default) complex division algorithm
Naive (default) complex division algorithm
Status: RESOLVED FIXED
Product: gcc
Classification: Unclassified
Component: middle-end
4.0.0
: P2 normal
: 4.0.0
Assigned To: Richard Henderson
:
: 19138 (view as bug list)
Depends on: 19486 19609 19953 19974
Blocks: 5900 19292 19693
  Show dependency treegraph
 
Reported: 2004-12-09 11:07 UTC by Paolo Carlini
Modified: 2005-02-24 14:03 UTC (History)
8 users (show)

See Also:
Host: Any
Target: Any
Build: Any
Known to work:
Known to fail:
Last reconfirmed: 2004-12-09 13:01:54


Attachments
A trivial testcase (279 bytes, text/plain)
2004-12-09 11:09 UTC, Paolo Carlini
Details
complex division with scaling (2.89 KB, text/plain)
2005-01-19 08:27 UTC, Thomas Koenig
Details

Note You need to log in before you can comment on or make changes to this bug.
Description Paolo Carlini 2004-12-09 11:07:56 UTC
Apparently, the default algorithm (the other one available is selectable via
flag_complex_divide_method = 1) is the naive one, which is not able to deal
correctly with "large" denominators. For some additional details, see:

  http://gcc.gnu.org/ml/gcc-bugs/2004-12/msg00820.html

I'm attaching a trivial, pure C, testcase, showing at least inconsistency in the
behaviors for real vs complex types. Swithing to flag_complex_divide_method = 1
fixes simple testcases, but maybe we really want the full C99-conforming algo,
as sketched in Annex G. In that case this would be a (high priority ;) request
for enhancement.
Comment 1 Paolo Carlini 2004-12-09 11:09:13 UTC
Created attachment 7712 [details]
A trivial testcase
Comment 2 joseph@codesourcery.com 2004-12-09 11:17:42 UTC
Subject: Re:  New: Naive (default) complex division
 algorithm

On Thu, 9 Dec 2004, pcarlini at suse dot de wrote:

> I'm attaching a trivial, pure C, testcase, showing at least inconsistency in the
> behaviors for real vs complex types. Swithing to flag_complex_divide_method = 1
> fixes simple testcases, but maybe we really want the full C99-conforming algo,
> as sketched in Annex G. In that case this would be a (high priority ;) request
> for enhancement.

The Annex G example still comes with a warning that it may yield undue 
overflow: it illustrates how to get the treatment of infinities expected 
in that informative Annex, not how to avoid excess overflow in general.

Comment 3 Paolo Carlini 2004-12-09 11:26:27 UTC
> The Annex G example still comes with a warning that it may yield undue 
> overflow: it illustrates how to get the treatment of infinities expected 
> in that informative Annex, not how to avoid excess overflow in general.

I see, probably I should have pointed out that I'm reading in detail Annex G
for the first time, and I'm not a floating point expert. Nvertheless, in the
light of some basic notions of numerical analysis, seems rather obvious to me
that "we have a problem" which, moreover, nobody apparently discusses in detail.
Which is your take on the issue, more generally?
Comment 4 Paolo Carlini 2004-12-09 11:33:21 UTC
A naive idea: would make sense swithing from flag_complex_divide_method == 0
to flag_complex_divide_method == 1 basing on -ffast-math or other, finer grained,
floating point, switch?!?
Comment 5 joseph@codesourcery.com 2004-12-09 11:46:38 UTC
Subject: Re:  Naive (default) complex division algorithm

On Thu, 9 Dec 2004, pcarlini at suse dot de wrote:

> A naive idea: would make sense swithing from flag_complex_divide_method == 0
> to flag_complex_divide_method == 1 basing on -ffast-math or other, finer grained,
> floating point, switch?!?

I'd suggest naming the switch -fcx-limited-range (including in 
-ffast-math).  With due warning in the documentation that even the better 
algorithm currently implemented doesn't avoid all excess overflow.  Given 
that C99 does expect something better than the naive algorithm it might 
still make sense to switch to the better algorithm we have as the default 
(with this flag to control it, and in future hopefully the 
CX_LIMITED_RANGE pragma) for C.

Comment 6 Paolo Carlini 2004-12-09 13:01:54 UTC
Thanks Joseph. I will try to come up with a patch as soon as possible, but
please be gentle while reviewing it, would be my first one for the compiler
proper ;)
Comment 7 joseph@codesourcery.com 2004-12-09 13:25:49 UTC
Subject: Re:  Naive (default) complex division algorithm

On Thu, 9 Dec 2004, pcarlini at suse dot de wrote:

> Thanks Joseph. I will try to come up with a patch as soon as possible, but

There is no regression here, so I recommend holding off until 4.0 has 
branched.

Comment 8 Paolo Carlini 2004-12-09 13:31:15 UTC
> There is no regression here, so I recommend holding off until 4.0 has 
> branched.

Definitely. Thanks again.
Comment 9 Gabriel Dos Reis 2004-12-09 16:03:55 UTC
Subject: Re:  Naive (default) complex division algorithm

"joseph at codesourcery dot com" <gcc-bugzilla@gcc.gnu.org> writes:

| ------- Additional Comments From joseph at codesourcery dot com  2004-12-09 11:46 -------
| Subject: Re:  Naive (default) complex division algorithm
| 
| On Thu, 9 Dec 2004, pcarlini at suse dot de wrote:
| 
| > A naive idea: would make sense swithing from flag_complex_divide_method == 0
| > to flag_complex_divide_method == 1 basing on -ffast-math or other, finer grained,
| > floating point, switch?!?
| 
| I'd suggest naming the switch -fcx-limited-range (including in 
| -ffast-math).  With due warning in the documentation that even the better 
| algorithm currently implemented doesn't avoid all excess overflow.  Given 

I suppose you don't mean a warning each time that flag is used, but
rather a note in our documentation.

| that C99 does expect something better than the naive algorithm it might 
| still make sense to switch to the better algorithm we have as the default 
| (with this flag to control it, and in future hopefully the 
| CX_LIMITED_RANGE pragma) for C.

I strongly support the idea that the non-grammar school algorithm
should be the default.

-- Gaby
Comment 10 Wolfgang Bangerth 2004-12-23 23:05:00 UTC
A special case of what is covered in this PR is PR 19138. 
 
W. 
Comment 11 Paolo Carlini 2004-12-23 23:22:50 UTC
*** Bug 19138 has been marked as a duplicate of this bug. ***
Comment 12 Thomas Koenig 2005-01-19 08:27:48 UTC
Created attachment 7989 [details]
complex division with scaling

Here's an example of an algorithm with scaling, as C source code.

Because it uses frexp() to scale floting point numbers, it will
introduce nasty loss of precision on machines where FLT_RADIX != 2.
(Is gcc targetting any systems like that?  What about S/390?)
It also doesn't handle NaN and infinity.  I expect it to be slow
unless somebody writes a clever inline implementation of frexp() and
ldexp().  This is not meant to be production code.

Apart from that, it should be pretty good at preventing over- and
underflow, and it does not introduce any additional precision
problems due to additional divisions or multiplications.

	Thomas
Comment 13 Thomas Koenig 2005-01-20 10:45:44 UTC
For gfortran, this is a regression against g77.
See this test program:

$ cat complex-scale.f
      program main
      implicit none
      complex ca,cb,cc
      data ca /(+2.3955909e+19,-1.2258349e-38)/
      data cb /(+0.0000000e+00,-1.7649423e+19)/
      cc = ca/cb
      print *,cc
      print *,-real(ca)/aimag(cb)
      end
$ gfortran complex-scale.f
$ ./a.out
 (  0.000000    ,     +Infinity)
   1.357320
$ g77 complex-scale.f
$ ./a.out
 (0.,1.35731959)
  1.35731959
Comment 14 Paolo Carlini 2005-01-20 11:50:22 UTC
Hi Thomas, can you possibly figure out from the ChangeLog the author of the
algorithm with scaling, which currently seems broken, and add him in CC? We
want to fix it, in the first place.
Comment 15 Paolo Carlini 2005-01-20 12:10:26 UTC
A first implementation of the algorithm was already in 3_4, under the name
expand_cmplxdiv_wide (in optabs.cc), then Rth rewrote it in the tree-ssa branch
as part of the new tree-complex.cc (It would be mildly interesting to see if
the version in 3_4 works) In the meanwhile, I'm adding Rth in CC of 19486.
Comment 16 Toon Moene 2005-01-21 20:23:44 UTC
Subject: Re:  Naive (default) complex division algorithm

pcarlini at suse dot de wrote:

> ------- Additional Comments From pcarlini at suse dot de  2005-01-20 12:10 -------
> A first implementation of the algorithm was already in 3_4, under the name
> expand_cmplxdiv_wide (in optabs.cc), then Rth rewrote it in the tree-ssa branch
> as part of the new tree-complex.cc (It would be mildly interesting to see if
> the version in 3_4 works) In the meanwhile, I'm adding Rth in CC of 19486.

Indeed (ChangeLog.1):

Tue May 18 03:53:37 1999  Craig Burley  <craig@jcb-sc.com>

         Improve open-coding of complex divide:
         * flags.h: Declare new front-end-malleable flag.
         * toplev.c: Define new flag.
         * optabs.c (expand_cmplxdiv_straight): New function to do original
         open-coding.
         (expand_cmplxdiv_wide): New function to do new open-coding,
         from Toon Moene, with changes (call to emit_barrier, dropping
         of spurious `ok = 1;', plus the obvious `break;' -> `return 0;').
         (expand_binop): A bit of spacing fixing, while at it.
         Use new functions instead of inlining the open-coding code.

Comment 17 Thomas Koenig 2005-01-26 09:32:06 UTC
The recent fixes in complex handling, and the
scaled division algorithm, have eliminated the
Lapack regressions with -O0 at least on ia64-unknown-linux-gnu
(see http://gcc.gnu.org/bugzilla/show_bug.cgi?id=5900#c23 ).

This is a huge step.  I think PR 18902 should now be considered
critical, at least for gfortran.

        Thomas
Comment 18 Paolo Carlini 2005-01-26 10:54:39 UTC
> This is a huge step.  I think PR 18902 should now be considered
> critical, at least for gfortran.

Ok, I'll start on it immediately (after all, it's rather straightforward, now).
Then the maintainers will decide...

Comment 19 Paolo Carlini 2005-02-02 12:54:17 UTC
Rth is actively working on these issues.
Comment 20 CVS Commits 2005-02-24 09:24:39 UTC
Subject: Bug 18902

CVSROOT:	/cvs/gcc
Module name:	gcc
Changes by:	rth@gcc.gnu.org	2005-02-24 09:24:18

Modified files:
	gcc            : ChangeLog c-opts.c common.opt opts.c toplev.c 
	gcc/doc        : invoke.texi 

Log message:
	PR middle-end/18902
	* c-opts.c (c_common_post_options): Set flag_complex_method to 2
	for c99.
	* common.opt (fcx-limited-range): New.
	* opts.c (set_fast_math_flags): Set flag_cx_limited_range.
	* toplev.c (flag_complex_method): Initialize to 1.
	(process_options): Set flag_complex_method to 0 if
	flag_cx_limited_range.
	* doc/invoke.texi (-fcx-limited-range): New.

Patches:
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/ChangeLog.diff?cvsroot=gcc&r1=2.7576&r2=2.7577
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/c-opts.c.diff?cvsroot=gcc&r1=1.135&r2=1.136
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/common.opt.diff?cvsroot=gcc&r1=1.61&r2=1.62
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/opts.c.diff?cvsroot=gcc&r1=1.93&r2=1.94
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/toplev.c.diff?cvsroot=gcc&r1=1.943&r2=1.944
http://gcc.gnu.org/cgi-bin/cvsweb.cgi/gcc/gcc/doc/invoke.texi.diff?cvsroot=gcc&r1=1.582&r2=1.583

Comment 21 Richard Henderson 2005-02-24 09:25:21 UTC
Fixed.