Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Fortran sqrt of complex number -1 gives different results

This code

print *, sqrt(cmplx(-1))
print *, sqrt(cmplx(-1,0))
print *, sqrt((-1,0))
print *, sqrt(-(1,0))

gives me this output

(0.00000000,1.00000000)
(0.00000000,1.00000000)
(0.00000000,1.00000000)
(0.00000000,-1.00000000)

I believe that the correct algebra is sqrt(-1)=i. Why the result of the last line?

The compiler version is GCC 7.3.0, running on Linux openSUSE 42.2 (x86_64).

EDIT

Following @francescalus answer I have tried more cases:

print *, sqrt((-1,-0))
print *, sqrt((-1,-0.))
print *, (-1,-0)
print *, (-1,-0.)

and I get

(0.00000000,1.00000000)
(0.00000000,-1.00000000)
(-1.00000000,0.00000000)
(-1.00000000,-0.00000000)

So, it seems that my compiler support negative zeros for real numbers. So, I guess it is important to care when working with variables like this:

complex             :: asd 
asd=(1.,0.)
print *, sqrt(-asd)

Here I get again the wrong result, but the zero negative thing is more difficult to predict. I have so many questions! Do you know some other exmple that can induce a mistake? Do you have an advice to avoid this mistakes? Do you now some compiler flag to turn off the negative cero support for the GCC compiler?

like image 487
alexis Avatar asked Oct 31 '18 15:10

alexis


1 Answers

Fortran 2008 (13.7.159) defines the result of the sqrt function, for argument X, as (my emphasis):

The result has a value equal to a processor-dependent approximation to the square root of X. A result of type complex is the principal value with the real part greater than or equal to zero. When the real part of the result is zero, the imaginary part has the same sign as the imaginary part of X.

Your square roots do indeed have zero real part, so let's look at the sign of the imaginary part of your argument. What is the sign of the imaginary component of -(1,0)? If your processor supports signed zero, then it could well be negative. In which case, the imaginary part of the result should be negative according to the standard's requirement.

In all other cases, there'd be no reason to expect the imaginary component of the argument to be negative, rather than positive, zero.

like image 179
francescalus Avatar answered Sep 23 '22 02:09

francescalus