Is storage of COMPLEX in fortran a guarantee of two REALs?
Many FFT algorithms use complex numbers stored in alternating real and imaginary parts of the array. By creating a COMPLEX array and passing it to the FFT procedure, is it guaranteed that it can be moved to a REAL array (twice as large) with real and imaginary component variables?
subroutine fft (data, n, isign)
dimension data(2*n)
do 1 i=1,2*n,2
data(i) = ..
data(i+1) = ..
1 continue
return
end
...
complex s(n)
call fft (s, n, 1)
...
(and, btw, is measurement data (2 * n) just like saying it is REAL?)
source to share
I am only writing this answer because experience has taught me that as soon as I write such an answer, one of the real Fortran experts put me together to fix me.
I don't think the current standard or any of its predecessors explicitly states that it complex
should be implemented as two adjacent in memory reals
. However, I believe that this implementation is a necessary consequence of the standard definitions of equivalence
and common
. I don't think I have ever come across an implementation that has complex
n't been implemented as a pair reals
.
The standard ensures that it complex
can be converted to a pair reals
. So, given some definitions:
complex :: z
complex, dimension(4) :: zarr
real :: r1, r2
real, dimension(8) :: rarr
the following will do what you expect
r1 = real(z)
r2 = aimag(z)
Both of these features are elemental, and this is where the wrinkle comes in:
real(zarr)
returns a 4 element array of reals, like
aimag(zarr)
while
[real(zarr), aimag(zarr)]
is an 8 element array of reals with real parts zarr
followed by complex parts. Maybe,
rarr(1:8:2) = real(zarr)
rarr(2:8:2) = aimag(zarr)
will be good for you. I'm not sure if I can do anything neater.
Alexander is not the only one who can quote the standard! The part he quotes got me thinking about non-standard complex scalars. So I read on and I think that paragraph 6 of the section that he points to us is native
a scalar nonpointer object of any type not specified in items (1) - (5) occupies one unsettled store, which is different for each case and each set of type parameter values, and this is different from unspecified item storage units (4),
I don't think this has any impact on any of the answers here.
source to share
To add to the Mark answer , this is indeed stated in the standard: Section 16.5.3.2 "Storage Sequence":
2 In the context of a storage association
[...]
(2) a scalar nonpointer that is a double precision implementation or default complex occupies two contiguous numeric storage blocks ,
my accent.
As for storage association context: Section 16.5.3.1 "General" reads
1 Storage sequences are used to describe existing relationships among variables, shared blocks, and result variables. Storage location association is the union of two or more data objects that occurs when two or more storage sequences are shared or aligned with one or more storage units.
So this happens for blocks common
, explicit equivalence
and result variables. As Mark said , there is no explicit statement for the general case. I guess it is most convenient to always follow this to ensure compatibility.
Thanks to IanH for pointing this out!
source to share