Silverfrost Forums

Welcome to our forums

Complex variables

17 Feb 2017 10:53 #18839

In ClearWin+ we have support for real numbers i.e. %rf

It would be really advantageous, and simplify coding if there was a similar control that supported complex variables i.e. A + jB (or A + iB to keep the mathematicians happy!), as it gets a bit messy working separately with real and imaginary parts - especially if they require updating when the code is running.

I've tied myself in knots a few times over this and never came up with a slick solution. Perhaps somebody can point the way for me?

Ken

17 Feb 2017 12:31 #18840

My understanding is that COMPLEX variables are stored as 2 consecutive REALs, so surely this is a case for EQUIVALENCE? Say EQUIVALENCE complex C to real D(2), and if you want A & B, equivalence A to D(1) and B to D(2)?

I would be happy to be corrected, because COMPLEX numbers are a closed book to me, but I'm always ready to learn.

Eddie

17 Feb 2017 1:15 #18841

There are two issues here. The first one is support for I/O of complex variables in ClearWin+. I have nothing to say about that.

The other is the equivalencing of REAL and COMPLEX variables. Here, care is necessary, because such equivalencing can introduce a bug that can remain dormant for a long time. The Fortran standard states (I am taking the following quote from Fortran 2008, but I think that the rule applies to Fortran 95 as well):

16.6.6 Events that cause variables to become undefined 7 1 Variables become undefined by the following events. 8 (1) With the exceptions noted immediately below, when a variable of a given type becomes defined, all 9 associated variables of different type become undefined. 10 (a) When a default real variable is partially associated with a default complex variable, the complex 11 variable does not become undefined when the real variable becomes defined and the real variable 12 does not become undefined when the complex variable becomes defined.

According to this rule, the output of the following program is 'undefined'.

program tst
   implicit none
   real r(2)
   double precision d(1)
   equivalence (r,d)
!
   r(1)=1.0; r(2)=2.0
   print *,d
end program
17 Feb 2017 1:47 #18842

Mecej4,

Isn't that what the first 4 lines of your extract say - which is very believable, and useful to point out.

But don't the last 3 lines say something else, i.e. that you CAN deal with a COMPLEX as 2 REALs without that particular issue?

Even I can see that if you equivalence a double precision to two single precisions and mess around with one of the latter, you will generate absolute nonsense, but if a complex genuinely is two consecutive reals then each can surely be manipulated independently.

As in:

       program tst 
       real r(2) 
       complex d(1) 
       equivalence (r,d) 
       r(1)=1.0
       r(2)=2.0 
       print *,d 
       r(1)=3.0
       print *,d
       r(2)=1.0 
       print *,d 
       end

Standard conforming or not, it works.

Eddie

17 Feb 2017 3:52 #18843

Thanks for the suggestions. I've noted the word of caution about equivalence but as the following shows it does appear to work. (Actually I've never used equivalence before in 30 years dabbling in Fortran).

Just need to fix the sign of j when the imaginary part is negative +j-0.5 looks a bit clumsy (and is misleading).

module test
implicit none
complex(kind = 2) a_c, b_c, c_c
real(kind = 2) a(2), b(2), c(2)
equivalence (a_c, a)
equivalence (b_c, b)
equivalence (c_c, c)

contains
  integer function init()
  complex(kind = 2) j
  j = cmplx(0.d0,1.d0)
  a_c = 0.5d0 - j*0.5d0
  b_c = 0.5d0 - j*0.5d0
  c_c = a_c*b_c
  init = 1
  end function init

  integer function c_m()
  include<windows.ins>
  c_c = a_c * b_c
  c_m = 1
  end function c_m

  integer function window()
  include <windows.ins>
  integer i
  i = winio@('%cn%ws&','Multiply two complex numbers')
  i = winio@('%2nl&')
  i = winio@('%cn%8.1ob&')
  i = winio@('%bg&',rgb@(250,250,250))
  i = winio@('%^rf&',a(1),c_m)
  i = winio@('%cb&')
  i = winio@('%ws&','+j')
  i = winio@('%^rf&',a(2),c_m)
  i = winio@('%cb&')
  i = winio@('%ws&','x')
  i = winio@('%cb&')
  i = winio@('%^rf&',b(1),c_m)
  i = winio@('%cb&')
  i = winio@('%ws&','+j')
  i = winio@('%^rf&',b(2),c_m)
  i = winio@('%cb&')
  i = winio@('%ws&','=')
  i = winio@('%cb&')
  i = winio@('%`rf&',c(1))
  i = winio@('%cb&')
  i = winio@('%ws&','+j')
  i = winio@('%`rf&',c(2))
  i = winio@('%cb')
  window = 1
  end function window
  
end module test

program main
use test
implicit none
integer i
 i = init()
 i = window()
end program main
17 Feb 2017 4:40 #18845

Quoted from LitusSaxonicum

But don't the last 3 lines say something else, i.e. that you CAN deal with a COMPLEX as 2 REALs without that particular issue?

Only in the exceptional case where the association is partial. In my example, the association is complete.

This property (of the inactive member of an EQUIVALENCE pair becoming undefined as a side effect of the active member being given a new value) is something that was a big surprise for me when I first saw it a few years ago. The context was a program that worked correctly with the NAG compiler without checking options, but gave an 'undefined variable' with -C=undefined. I complained, and the NAG people straightened me out by quoting references to the Fortran standard.

17 Feb 2017 4:57 #18846

I have a feeling that the COMPLEX case ought to always work because a COMPLEX variable is axiomatically made up of two REALs. A DOUBLE PRECISION is never made up of two REALs. An Integer might be made up of a signed integer and an unsigned one (both of smaller length), but it would matter if you were big-endian or little-endian which was which.

The standards extract probably means it will work with COMPLEX, may work with INTEGER and definitely won't work with REAL!

Of course, what matters is whether it works with FTN95 and any other compiler that exploits Clearwin+, because when you add 40,000 lines* to a 400 line Fortran program you cease to care whether it will work with a compiler that doesn't 'do' CW!

Eddie

*Not to forget the 300 or more icons you drew ...

21 Feb 2017 12:13 #18865

There are various ways to avoid EQUIVALENCE in this program. Here is one that defines a user operator for complex multiplication. (A simple function call is a less frightening alternative to a user operator.) For complex addition and subtraction one can just use the default + and - on the whole arrays.

module test 
implicit none 
real(kind = 2) a(2), b(2), c(2)
interface operator (.x.)
 module procedure cmult
end interface  

contains
  function cmult(a,b)
  real(kind = 2),intent(in)::a(2),b(2)
  real cmult(2)
  cmult(1) = a(1)*b(1) - a(2)*b(2)
  cmult(2) = a(2)*b(1) + a(1)*b(2)
  end function cmult
   
  integer function init() 
  a = (/0.5d0,-0.5d0/)
  b = (/0.5d0,-0.5d0/)
  c = a.x.b 
  init = 1 
  end function init 

  integer function c_m() 
  c = a.x.b 
  c_m = 1 
  end function c_m 

  integer function window() 
Please login to reply.