2

Is there an intrinsic function in Fortran to compute the square modulus |z|^2 of a complex number z?

If not, is there a simpler or better way to compute it than the following?

REAL(z, precision_specifier)**2 + AIMAG(z)**2
phuclv
  • 37,963
  • 15
  • 156
  • 475
norio
  • 3,652
  • 3
  • 25
  • 33
  • 1
    You could use `abs(z)**2` (see [here](https://gcc.gnu.org/onlinedocs/gcc-3.4.4/g77/Abs-Intrinsic.html)) – d_1999 Mar 15 '17 at 13:37
  • 1
    I think the point is to avoid taking an expensive sqrt. – Ross Mar 15 '17 at 13:53
  • That sounds inefficient because it would compute a square root of a real number and then square it again. – norio Mar 15 '17 at 13:53
  • 2
    Yes, that's true -- it might be worth mentioning in the question if this is a concern for your application (it might help make some of the answers clearer as comments can be transient). – d_1999 Mar 15 '17 at 17:07
  • Not only is sqrt to be avoided, but, in compile modes other than complex-limited-range, the compiler will expand code or call a library function to avoid overflow in multiplication (which you appear not to be concerned about). – tim18 Mar 15 '17 at 21:51

4 Answers4

5

No, the Fortran standard defines no such intrinsic routine. Nor, as far as I am aware, does any of the current crop of widely-used compilers provide such a routine.

If OP really wants to avoid an expensive sqrt, and doesn't like her existing solution, OP could try:

real :: rslt
real, dimension(2) :: parts
complex :: z
...
parts = transfer(z, parts)
rslt = dot_product(parts, parts)

or, given the same declarations, this might be preferred

rslt = dot_product(transfer(z, parts), transfer(z, parts))

As always, if performance matters to you, measure it.

Simpler or better ? You decide.

As always, transfer is not to be operated by persons under the age of 18, or persons under the influence of alcohol or other performance-impairing drugs.

High Performance Mark
  • 77,191
  • 7
  • 105
  • 161
2

Just for fun, here is some attempt to use a (good old??) statement function for defining abs2() with the hope of being inlined... (not recommending its use, just experiment!) 

program main
    implicit none
    integer, parameter :: dp = kind(0.0d0)
    complex(dp) z

    real(dp) abs2
    abs2( z ) = real(z)**2 + aimag(z)**2    ! (could be included or macro-ed??)

    ! abs2( z ) = conjg( z ) * z          ! maybe slower
    ! abs2( z ) = z % re**2 + z % im**2   ! near future
    ! intrinsic :: abs2                   ! best solution   

    z = ( 1.0_dp, 2.0_dp )

    print *, abs( z )**2
    print *, abs2( z )

    print *, abs( z + 1.0_dp )**2
    print *, abs2( z + 1.0_dp )
end program

Result:
   5.0000000000000009     
   5.0000000000000000     
   8.0000000000000018     
   8.0000000000000000

I hope z% re and z% im will be coming soon... (not yet with gfortran-6 and ifort-16)

roygvib
  • 7,218
  • 2
  • 19
  • 36
  • If you want to make this prettier (and isn't that the goal with statement functions?): then `real(z,dp)` can be replaced by `real(z)` with the same effect. – francescalus Mar 15 '17 at 15:27
  • real(z) also worked... :o I didn't know this returns the same kind as z... – roygvib Mar 15 '17 at 15:48
2

Just exploit well-known formula (a+i*b)*(a-i*b) = a^2 + b^2

 SquaredModulus = real part of (Z * CONJG(Z))

CONJG is complex conjugate

MBo
  • 77,366
  • 5
  • 53
  • 86
0

https://gcc.gnu.org/onlinedocs/gfortran/ABS.html

You are probably looking for CABS() or CDABS() and then squaring it. I just checked and it worked with f95 compiler.

RLLL
  • 21
  • 2
  • 1
    Welcome, I suggest you to take the [tour]. There is no point using CABS or CDABS. The generic ABS provided since Fortran 77 is just fine. The point of the question was probably to avoid the square root done by ABS and the squaring it again. See the comments under the question. – Vladimir F Героям слава May 08 '19 at 15:59