[Bug fortran/107753] gfortran returns NaN in complex divisions (x+x*I)/(x+x*I) and (x+x*I)/(x-x*I)

["weslley.pereira at ucdenver dot edu" <gcc-bugzilla@gcc.gnu.org>]

https://gcc.gnu.org/bugzilla/show_bug.cgi?id=107753

--- Comment #7 from Weslley da Silva Pereira <weslley.pereira at ucdenver dot edu> ---
(In reply to anlauf from comment #3)
> I guess the reporter assumes that gcc uses a clever algorithm like Smith's
> to handle such extreme cases of complex division.  Not sure if that one is
> available by some compilation flag, and I think it would impact performance.
> 
> In any case, if the reporter wants to get robust results and in a portable
> way, I would advise him to change/fix his algorithm accordingly.  It appears
> that a few other compilers behave here like gfortran.

Thanks for the suggestion of changing the algorithm that needs such a division.
What are the ranges for nominator and denominator where one should rely on the
intrinsic complex division? Maybe this is the good question to ask. Then we can
build algorithms that attend such requisites. LAPACK has cladiv and zladiv,
which are routines for complex division that avoids unnecessary under and
overflow. They are used in many parts of the code. This implementation is for
sure less efficient than the intrinsic complex division, but we rely on it
because of robustness.

More data for the discussion:
1. In a Ubuntu 18.04.5 LTS, using GNU Fortran 7.5.0, I tested optimization
flags `-O` but still reproduce the wrong result for complex divisions with huge
numbers. See
https://github.com/Reference-LAPACK/lapack/issues/575#issuecomment-910616816
that used the code from
https://github.com/Reference-LAPACK/lapack/blob/master/INSTALL/test_zcomplexdiv.f.
This is the test currently in LAPACK 3.11.0.
2. I have just reproduced what was reported in
https://github.com/Reference-LAPACK/lapack/issues/575#issuecomment-910616816 in
my Ubuntu 20.04.5 LTS, using GNU Fortran 9.4.0.
3. I noticed that the optimization flag is unable to target divisions like
`x/x` depending on where they are inside a program.
4. My Ubuntu 20.04.5 LTS with compiler ifort 2021.7.1 computes the complex
division `x/x` accurately even for the case of huge numbers. Scenarios tested:
   - I tested the program in
https://github.com/Reference-LAPACK/lapack/blob/master/INSTALL/test_zcomplexdiv.f
and the one in https://godbolt.org/z/b3WKWodvn.
   - I tested ifort with flags -fp-model precise and -fp-model fast. The latter
enables more aggressive optimizations on floating-point data.
   - I tested compilation with optimization flags -O0, -O, -O1, -O2, -O3. 

Here is the implementation of the complex division in LAPACK if it somehow
helps the discussion:
https://netlib.org/lapack/explore-html/d8/d9b/group__double_o_t_h_e_rauxiliary_gad1c0279ec29e8ac222f1e319f4144fcb.html