That is, you don’t need to know the length of an array to get its imaginary part (nor the real). And the real and imaginary part should be spread in memory at fixed interval (in order to be associated by a pointer).
You can still have [I1,R1,I2,R2,…] but I don’t see any advantage in doing it compared to [R1,I1,…], or to have some alignment bit, but I don’t see why one have to have these alignment bits in the complex and not in the corresponding real array.
That’s the way equivalence works. Defining a variable that is equivalenced to a variable of a different type causes that other variable to become undefined. It is not standards conforming to make the assumption that the bit pattern for a particular value has any meaning when interpreted as a different type.
Seems reasonable, but I’d suggest that it be stated that it is then invalid to refer to z or x within the body of the associate construct.
It breaks code that followed the rules. I.e.
interface foo
subroutine foo_complex(c)
complex :: c(*)
end subroutine
subroutine foo_real(r)
real :: r(*)
end subroutine
end interface
real :: r(10)
complex :: c(10)
call foo(r)
call foo(c)
end
The fact that existing codes rely on non-standard behavior is not the fault of the language. It is the fault of the programmers. I don’t think we should poke holes in the language because some people didn’t realize their code was invalid (or at least relying on unspecified behavior).
exactly.
Let’s talk about some hypothetical (but potentially applicable) example. Let’s say I represent my floating point values in 6 bits, but the numeric storage unit is 8 bits (cause I use 8 bits for integers). Well, I could store my complex values in the first 12 bits of the 16 bits of storage they “occupy”. There are real-world examples of this kind of partial occupation in 80 bit reals.
We were talking about the move_alloc() type of generalization for real and complex entities. As I said previously, pointers remain defined and associated after move_alloc(), so it seems like the real and complex generalization to move_alloc() might retain that feature.
Regarding EQUIVALENCE, let’s discuss the following code:
program equiv
implicit none
real :: x, y, a(2)
complex :: z
equivalence (x, z, a(1)), (y,a(2))
x = 2.
y = 3.
write(*,*) a
write(*,*) z
a = [4.,5.]
write(*,*) x, y
write(*,*) z
z = (6.,7.)
write(*,*) x, y
write(*,*) a
end program equiv
$ gfortran equiv.f90 && a.out
2.00000000 3.00000000
(2.00000000,3.00000000)
4.00000000 5.00000000
(4.00000000,5.00000000)
6.00000000 7.00000000
6.00000000 7.00000000
Apart from the list directed output, I think the behavior of that whole program is specified by the fortran standard. Namely, the value of the complex variable z is defined when either the x and y pair of values are set, or when the array a(:) is set. Also, the values of the x and y pair and of the array a are defined when z is defined. And for completeness, the array a(:) is defined when either the x and y pair or the complex z is defined. I don’t think there is any ambiguity regarding whether x is associated with a(2) or with z%im, I think storage sequence association defines everything in that program.
As I stated previously, this is a different situation when equivalence of say, an integer and a real occur. In that case, the standard states the the integer becomes undefined when the real variable is defined, and the real variable becomes undefined when the integer is defined.
The situation with real and complex is exactly the opposite. The standard specifies that z becomes defined when its components become defined through equivalence (or also through the z%re and z%im syntax). Another situation, not demonstrated above, involves unformatted i/o. The standard requires that z becomes defined when two real values are read, and it requires that two real values become defined when a complex z is read. If common blocks are involved, the same things are required to occur there. If z is set in one place, then x and y (or a(1:2)) are defined through storage sequence association elsewhere.
These EQUIVALENCE rules are limited to only default real and complex entities. This is for historical reasons, before f90 there was only one complex kind, corresponding to the default real. With f90, the EQUIVALENCE rules were not updated to account for the KIND system, and now since f2003, the c_loc() and c_f_pointer() facility has similarly not been updated to account for the KIND system.
That means that your imaginary part cannot start in the middle of a byte.
It should start at a byte boundary, likely at the byte multiple corresponding to it’s size.
Your code compiles fine on Intel Fortran version 19. But so does mine here .
The [real_1, imag_1, real_2, imag_2, …] convention is also followed therefore points 1 and 2 are satisfied. Thus I ask again: what is broken?
FFTW is written in C and has no obligation to follow Fortran conventions. In order to satisfy the Fortran standard a ‘good’ programmer would have to take their complex array, copy it to a real matrix, perform the FFT and then copy it back to the original complex array. Nobody does this in practice because it would incur an unacceptable inefficiency.
If the equivalenced objects have differing type or type parameters, the EQUIVALENCE statement does not cause type conversion or imply mathematical equivalence.
An EQUIVALENCE statement, a COMMON statement, or an ENTRY statement can cause storage association of storage sequences.
Partial association shall exist only between
…
an object that is default complex, double precision real, or of numeric sequence type and an object that is default integer, default real, default logical, double precision real, default complex, or of numeric sequence type.
When a numeric storage unit becomes defined, all associated numeric storage units of the same type become defined. When an entity of double precision real type becomes defined, all totally associated entities of double precision real type become defined.
When a default complex entity becomes defined, all partially associated default real entities become
defined.
When both parts of a default complex entity become defined as a result of partially associated default real or default complex entities becoming defined, the default complex entity becomes defined.
When a scalar variable of intrinsic type becomes defined, all totally associated variables of different type become undefined. When a double precision scalar variable becomes defined, all partially associated scalar variables become undefined. When a scalar variable becomes undefined, all partially associated double precision scalar variables become undefined.
This is as much of the relevant bits of the standard as I could find. It seems that yes, the storage association you show is valid, but it still does not say what the actual output of that program should be. It says that defining x and y (or a) does define z, but not with what value. And it says that defining z does define x, y and a, but again, not with what values.
It doesn’t necessarily say that it must start at a byte boundary, but I will grant that I don’t think there’s any possible implementation with the bits mixed together. I.e. an 80 bit real usually takes up 128 bits (16 bytes), and usually a numeric storage unit is considered to be 8 bytes, so an 80 bit real occupies 2 storage units, and thus a complex with the same kind will occupy 4 storage units (32 bytes). I don’t think there’s anything in the standard that says the imaginary part couldn’t be stored in the first 10 bytes, and the real part stored in the next 10 bytes, (i.e. all the data packed into 20 bytes, with the unused 12 bytes at the end) and still be able to have a pointer to the real part.
But it would be invalid if we allowed passing real arrays to complex arguments and vice-versa. The currently valid generic interface would be invalid because it would not be possible to distinguish which specific procedure to call.
Your example is not standards conforming, whether ifort 19 compiles it or not.
$ gfortran -pedantic invalid.f90
invalid.f90:5:11:
5 | call foo(z)
| 1
Error: Type mismatch in argument ‘x’ at (1); passed COMPLEX(8) to REAL(8)
$ ftn -eDC invalid.f90
subroutine foo(x)
^
ftn-1616 ftn: WARNING FOO, File = invalid.f90, Line = 11, Column = 16
Procedure "FOO" is also referenced at line 5 (invalid.f90). Argument 1 is scalar. This argument is array-valued.
FWIW I’ve found that the Intel compilers are the least reliable when it comes to enforcing standards conformance.
Does it use float _Complex, which is promised to be interoperable with complex, or does it do something that is not promised to be interoperable with complex?
Ok, but that doesn’t mean it’s not relying on unspecified behavior.
Yes, but this is the point: there exists a widely used compiler which violates the standard and breaks nothing except the letter of the standard. This is evidence that the de facto standard could be made the official standard without undue harm.
This only means that for an FFTW function with prototype
void some_func(int n, fftw_complex vals[n]);
The equivalent Fortran interface would be
subroutine some_func(n, vals) bind(C, name="some_func")
use iso_c_binding, only: c_double, c_int
integer(c_int), intent(in), value :: n
real(c_double) :: vals(n*2)
end subroutine
And that FFTW expects you to interpret those vals as R1, I1, R2, I2, …
But I went down the rabbit hole a bit, and this is where things get interesting. That page also says
if you have a C compiler (such as gcc ) that supports the C99 revision of the ANSI C standard, you can use C’s new native complex type (which is binary-compatible with the typedef above).
But that parenthetical is not true. The C standard doesn’t say that double _Complex has that format. Maybe most implementations do, but again, it’s relying on unspecified behavior.
So if you used that interface to call into an FFTW library compiled using the new C _Complex type, you’re potentially not binary compatible. They say on that page
C++ has its own complex<T> template class, defined in the standard <complex> header file. Reportedly, the C++ standards committee has recently agreed to mandate that the storage format used for this type be binary-compatible with the C99 type
But the C Standard doesn’t say what the binary format of that type is. So what are they doing? I looked it up, and yes, the C++ standard does say their Complex type is basically just an array of two values, real part first, imaginary part second. But if you scroll towards the bottom of that page, they link to a discussion about their decision. The initial summary says
The absence of explicit description of std::complex<T> layout makes it imposible to reuse existing software developed in traditional languages like Fortran or C with unambigous and commonly accepted layout assumptions.
But again, the Fortran and C standards do not actually dictate what the layout is.
So if the Fortran standard dictates what the storage format of the complex type is, it will have to add a caveat that it is only interoperable with the equivalent C types if the companion processor is implemented using the format that the Fortran standard has specified.
And at any rate, changing the Fortran standard to say that it is valid to pass complex array arguments to real array arguments would make currently valid programs invalid. Complex and real dummy arguments are currently distinguishable for the purposes of generic resolution. That change would make them not distinguishable, and it would be become invalid to have specific procedures in a generic interface that are distinguishable only on real vs complex dummy arguments, where as currently it is valid. The committee takes great pains to avoid making changes that cause currently valid programs to become invalid.
Isn’t that like saying the statement x=1.0 defines the variable x, but not with what value, or that x=1.0 defines a(1), but not with what value?
How exactly should the standard be worded to require what seems to be the obvious intent? Namely, that the assignment x=1.0 defines the value of x, the value of a(1), and the value of z%re, to all have the value of 1.0.
The argument above was that the layout was specified by the standard, but not the actual values. Now it is the layout itself that is not specified?
I would say that the equivalence statement:
equivalence (x, z, a(1)), (y,a(2))
does fully specify the layout of z through storage sequence association. However, as stated previously, this only applies to default real and complex kinds, not to all possible real and complex kinds.
And, why the Fortran, and C++ committees don’t meet and agree once and for all about the format of a complex number?
It looks like it is always [Real, Imag] of the corresponding real type, in all the compilers and in all the codes, or I’m wrong?
Why nobody want to put a sign on it?
Is there a better and viable alternative?
Should the IEEE754 define also a complex datatype?
That’s not what it’s saying. x=1.0 defines x with the value of the expression on the rhs, and in your example, a(1) is unambiguously storage associated and of the same type, so defines it to have the same value.
At present, all that is stated by the standard is that z becomes “partially defined” when x=1.0 is executed. It doesn’t say which part of z. It could be the real part, or it could be the imaginary part, or if as I pointed out in an example above the parts aren’t exactly aligned with separate numeric storage units, something different.
I never said that the layout was specified. It is because the layout is unspecified that partially defining a complex object via storage association is undefined.
I’m not entirely sure, but my guess was that there was some thought that magnitude and angle might have been a desirable representation for complex numbers. When Fortran added the %re and %im components, they drastically reduced the chances of anybody actually using that representation. And then C++ also closed the door to a different representation than real part and imaginary part. A C compiler could still do it and be valid though.
A Fortran compiler could in theory use a different representation (I think), but it would be really hard to get pointer association of real variables to parts of a complex variable to work. I.e. it could insert conversion code for something like r = c%re, but it would have to remember that a real pointer is associated to a complex part and do the conversion on every access, I.e. r => c%re; x = r; r = 1.
So, if multiple programming languages all wanted to decide that “we’re all going to store complex numbers in this format, aligned this way, to help with interoperability” that would be fine. What the Fortran standard should not do though is violate type safety and break backwards compatibility.
And then?
Copying in copying out?
To all C function as well?
This already happens with arrays. They can be internally implemented with a descriptor or direct pointer to contiguous data, and indeed they have to be copied if a non-contiguous array is passed to a function (such as a C function) that only accepts a pointer.
I know of course
Let assume that some optimized compiler (like LFortran, just to take one ) transform on the fly an Interleaved-Complex to a split one (may be a block split) before any complex instruction…
Well I’m just joking, cheers.
I think it is possible for a compiler to support these alternative representations, and they obviously can’t be used in all cases, but I think it is possible to robustly know for the compiler when they can be used and when not. One can also consider a language extension (or pragmas) to allow the user to select it.
Many other alternative representations are possible for arrays: row order, column order, as well as various device accelerator arrays, the polar representation of complex numbers (didn’t even think of this one until @everythingfunctional mentioned it above), etc. The trick is to allow this without breaking any existing code, but I think it’s possible.
