I think you are doing great job in popularising and explaining the necessity of your proposal. 
1 Like
If it is contiguous you can always use pointer remapping to see a (2,N) as a (2N) array. I don’t see it as a problem.
The problem is how is going to be stored a complex array: is it [real_1, imag_1, real_2, imag_2, …] or [real_1, real_2, …, imag_1, imag_2, …].
It is generally stored as in the first options (as far as I know), but is the Standard prescribing it?
That, I think, is the real problem.
Your proposal want to allow to point a real pointer to a complex target variable and vice versa. Someone objected that it will invalidate type safety. Well, why not protect it inside a sort of select statement? That would solve the type safety issue. For your case it will become:
complex, target :: wz(N)
real, pointer :: wr(2*N)
select complex as real (pr => wz)
wr(1:2*N) => pr
end select
And you have got it without worrying about type safety. The converse can be also allowed in case of simple contiguity as in your proposal.
But, there is a but, it assumes the Standard fix a preferred storage for the complex numbers like [R,I,R,I,R,I…].
If one wants to keep open the possibility of different way to store complex number…
Well, I don’t know.
Actually, it could be even more interesting (if you allow me to do some brainstorming).
Let add a new statement: view. To be used like the select statement:
view <type one> as <type two> (a => b)
...
end view
Where <type one> and <type two> should be in some way compatible. Where, for example, one is real and the other complex, or between two kind of integers when the size of one is a multiple of the size of the other…
But just some brainstorming…
I think I’m starting to understand the use case a bit better. I still think it violates strong typing, but it may be sufficiently specified by the standard already to expect this to work. However, I’d want to be a bit more strict about it and do something like this.
complex, pointer, contiguous :: c(:)
real, pointer, contiguous :: r(:)
allocate(c(100))
! can work with c
call reinterpret(from = c, to = r)
! now can work with r
call reinterpret(from = r, to = c)
! now can work with c again
Basically something like move_alloc for pointers, but the actual arguments must be contiguous, rank 1, and either complex or real with the same kind. There might still be some loopholes I haven’t thought of that make this still not a good idea, but I think this at least breaks things the least.
It’s not a problem per se, it’s just that I can’t see which use case such a shape would enable. So at the end it would be just an additional step towards a flatter pointer.
One of the purpose of the proposal is to standardize for good this de facto standard layout.
I’d more favor an extension of the current pointer behavior to a kind of “checked storage association by a pointer whenever it is safely possible”. We would probably only need one or two intrinsic procedures and no changes in the current strong type checking mechanism anywhere (which is definitely worth of being kept in my opinion).
real, allocatable, target :: r(:)
complex, pointer :: c(:)
integer :: stat
allocate(r(100))
! Built in function to check, whether two types can be safely storage associated.
if (.not. is_storage_compatible(r, c)) error stop "Can not safely storage associate"
! Let's c point to the storage provided by r, like a normal pointer
! The target must be allocated at this point, have the right type and right size
! The compiler can make further checks, whether the association is safe.
call associate_with_storage(r, c, stat=stat)
if (stat /= 0) then
error stop "Storage association failed"
end if
! From here, you can use c as a normal complex pointer wherever you want to access
! the data of the real array as complex.
The above mechanism could also allow for other storage associations like char to integers (to manipulate something byte-wise), etc.
The compiler had the chance to refuse the storage association, whenever it is not possible to do it safely. (For example, if the storage of an item of one type does not correspond to the storage of multiple items of the other kind, or if the size of the allocated storage is not commensurable with an integer number of items of the other type – e.g. you have an odd number of reals and try to associate a complex pointer with it).
EDIT: Here is a Fortran prototype implementation of the associate_with_storage() subroutine. But this would have to be an intrinsic one, so that the compiler cancheck the safety of the storage association (part of it already at compile time, part of it at runtime).
Fortran prototype implementation pointer based storage association
module storage_assoc
use, intrinsic :: iso_c_binding
implicit none
interface associate_with_storage
module procedure associate_with_storage_real1_cmplx1
module procedure associate_with_storage_cmplx1_real1
end interface associate_with_storage
contains
subroutine associate_with_storage_real1_cmplx1(trg, ptr, stat)
real, pointer, intent(in) :: trg(:)
complex, pointer, intent(out) :: ptr(:)
integer, intent(out) :: stat
if (mod(size(trg), 2) /= 0) then
stat = 1
return
end if
call c_f_pointer(c_loc(trg), ptr, [size(trg) / 2])
stat = 0
end subroutine associate_with_storage_real1_cmplx1
subroutine associate_with_storage_cmplx1_real1(trg, ptr, stat)
complex, pointer, intent(in) :: trg(:)
real, pointer, intent(out) :: ptr(:)
integer, intent(out) :: stat
call c_f_pointer(c_loc(trg), ptr, [2 * size(trg)])
stat = 0
end subroutine associate_with_storage_cmplx1_real1
end module storage_assoc
program test
use storage_assoc
implicit none
complex, pointer :: cc(:)
real, allocatable, target :: rr(:)
integer :: ii, stat
! If you try this with an odd nr. of elements, program should stop with error
rr = [(ii, ii = 1, 10)]
call associate_with_storage(rr, cc, stat)
if (stat /= 0) error stop "Could not assocaite rr with cc"
write(*, "(2F8.2)") cc
end program test
1 Like
This is offered as an alternative syntax in the proposal; the name and precise arguments of the routine are of course still subject of debate.
Alternative syntax 2 (similar to the
c_f_pointer()subroutine):call complex_pointer(r,c[,shape]) call real_pointer(c,r[,shape])
1 Like
I have just realized that the standard don’t allow a storage like:
[R1,R2,R3,…,I1,I2,I3,…]
Because where the imaginary part starts for an assumed size dummy argument like this:
complex :: wz(*)
So we are left with not many options:
[R1,<alignment bits>, I1, <alignment bits>, R2, <alignment bits>, I2, <alignment bits>,…]
or the converse with imaginary part first:
[I1,<ab>, R1,<ab>,…]
But it will be likely that also in the corresponding real there will be alignment bits like the infamous real(10) :: .
I don’t see other options that will be compatible with assumed size array.
So I don’t think will be a big issue to specify that the real part come first and then the imaginary part.
You may need to check for contiguity. But the same routines can be easily written in C using ISO_Fortran_binding.h, even though in a not conforming way as you may need to directly manipulate the CFI_cdesc_t structure and that is forbidden.
Then, even if your target is not contiguous, everything will work if it can work (like a real pointer to a complex target, or in other situation).
By the way I still prefer to have an additional dimension of size two for the real array for clarity even though it may mean an additional step.
One may also need an additional concept, that of contiguity up to a given dimension.
That may mean for example that an array:
A(N,M,K,H)
is contiguous up to 2 dimension if all A(:,:,k,h) for all k and h (as index) are contiguous, I think one may just use the first index and check that A(:,:,1,1) is contiguous (if 1 are the lower bounds in that dimensions).
Then, following your idea without the additional dimension, one can request that a real array (the target),
WR(N,M, …), contiguous up to the first dimension and where N is even have a complex pointer array WZ(N/2,M,…) (the pointer) associated to it that will still be contiguous up to the first dimension.
Why I like this additional dimension? Because I can always associate (now I’m using the correct phraseology) a real pointer to a complex target even if the complex target is not contiguous if I add an additional dimension. Without the additional dimension it will not be always possible.
If I have a complex target WZ(M,N,K) I can always associate a real pointer to it like WR(2, M, N, K) which is contiguous in the first dimension. While I may not be able to associate to it a real pointer like WR(2*M,N,K) (unless WZ was contiguous up to the first dimension).
As before some brainstorming
I would prefer something more modest which doesn’t require any new functions or statements. Namely:
- The de facto convention of storing complex numbers as [real_1, imag_1, real_2, imag_2, …] is standardised.
- The following code compiles and runs
program test
implicit none
complex(8) z
call foo(z)
print *,z
end program
subroutine foo(x)
implicit none
real(8), intent(out) :: x(2)
x(1)=1.d0
x(2)=2.d0
end subroutine
… without error. A warning or information statement would be fine (and even welcome). Right now to make this work we have to include -fallow-argument-mismatch with GFortran or hide foo by putting it in a separate file.
This is an interesting idea. With EQUIVALENCE or with the the c_loc()+c_f_pointer() approach, the storage can be referenced with two variable names. With a move_alloc() kind of approach, only one variable name would be allowed at any time (ignoring extraneous pointer assignments for the moment). A possible problem is that some kind of attribute would be required, similar to the allocatable attribute with move_alloc(), to designate which entities could appear in the new statement. If that attribute could be attached to dummy, automatic, module, and allocatable arrays, I think that might be acceptable to most programmers.
The final result would be that at any point in the program, the storage could be accessed through its real variable name, or through its complex variable name, but not both simultaneously.
Then there is the issue of extraneous pointer assignments. if rp=>r is a real pointer assignment, and then the real_complex_move_alloc(from=r,to=c) statement is executed, what should happen to the previous rp pointer? It would then be illegal to reference r, but should rp still be associated to the same storage as it was before? The same issue then arises the other way. After a cp=>c pointer assignment, and another real_complex_move_alloc(from=c,to=r) statement, should cp still be associated with its prior storage?
This is not modest. It breaks type safety and backwards compatibility (as outlined in my example above).
I do kind of like this name better. As for the implementation I would rather see
subroutine associate_with_storage_real1_cmplx1(trg, ptr, stat)
real, pointer, intent(inout) :: trg(:)
complex, pointer, intent(out) :: ptr(:)
integer, intent(out) :: stat
if (mod(size(trg), 2) /= 0) then
stat = 1
return
end if
call c_f_pointer(c_loc(trg), ptr, [size(trg) / 2])
stat = 0
nullify(trg)
end subroutine associate_with_storage_real1_cmplx1
And then, given that this seems possible to implement as a library, why not just add this functionality to stdlib and not worry about adding to the standard? Ignoring the minor bit about the c_f_pointer line being technically not standards conforming today, since I don’t know of any implementation where this would break.
This is exactly the kind of loop-hole I was alluding to. IMO it should be invalid to use other pointers associated with the original variable/storage after executing such a statement. I.e. the association status of other pointers referring to the target becomes undefined. Of course most implementations probably would just ignore that and not detect any problems.
@everythingfunctional Changing the attribute of the target variable trg to from intent(in) to intent(inout) as you suggest, would prohibit us to pass a non-pointer actual argument with the target attribute, which would be important use case in my opinion.
For sure, one could implement this in stdlib. The problem is, how do we check whether the functionality works reliably on a given architecture (whether the memory pattern of the data in real, dimension(2) is equivalent to that of complex, dimension(1)? There is no warranty by the standard, is there? If the standard warranted at least this, one could safely use the Fortran-wrapper, and won’t necessarily need an intrinsic routine for it.
Yes, but it also prohibits (at least easily) interpreting the same bit pattern two (or more) different ways at the same time. That’s what type safety is. You could of course still get around it with
real, target :: a(...)
real, pointer :: b(:)
complex, pointer :: c(:)
b => a
call associate_storage_with(from = b, to = c)
but one of the design principles for languages is to make dangerous things harder.
Here’s what I think are the relevant pieces to connect the dots:
The values of a complex type are ordered pairs of real values. The first real value is called the real part, and the second real value is called the imaginary part.
a nonpointer scalar object that is double precision real or default complex occupies two contiguous numeric storage units,
The values of C_FLOAT_COMPLEX, C_DOUBLE_COMPLEX, and C_LONG_DOUBLE_COMPLEX shall be the same as those of C_FLOAT, C_DOUBLE, and C_LONG_DOUBLE, respectively.
A Fortran intrinsic type with particular type parameter values is interoperable with a C type if the type and kind type parameter value are listed in the table on the same row as that C type.
…
C_FLOAT_COMPLEX → float _Complex
etc.
I don’t think either the C or Fortran standards quite dictate exactly what the bit storage of the real and imaginary parts are (real part first or imaginary part first, or heck even interleaved), but they better agree, and the size of a scalar complex object better be exactly the size of a size 2 real array.
I don’t think the standard should dictate exactly the bit patterns for complex storage. So really, what you want to do if you want to store real, then complex data in the same memory, is have a temporary copy, something like the following psuedo code.
do i = 1, size(r), 2
tmp_c%re = r(i)
tmp_c%im = r(i+1)
c => r(i:i+1)
c = tmp_c
end do
So I guess really we should call it convert_to_complex_pointer or convert_to_real_pointer, so it can do that transformation for you (although for sane implementations, and every one that I know of, there probably is no “transformation” needed).
Can you explain why you think this should be the semantics?
That’s not the way move_alloc() works. Pointers remain defined/associated to their storage after move_alloc().
I am a bit uneasy about providing this in the language initially as an intrinsic. Can we “jail” it instead by extending ASSOCIATE construct’s syntax
associate-name => selector
selector is expr
or variable
or complex-real-selector
ASSOCIATE ( zre => REAL z, xco => COMPLEX x)
ifx version 2024.0.2 20231213
Gives an internal compiler error on:
program test_inc_pointers
implicit none
integer, parameter :: N = 20
complex, target :: wz(N)
real, pointer :: wr(:)
integer :: i
wr => wz%re
print *, wr
end program
Moreover this program, that I think is standard conforming is not working as expected in gfortran 11.4.0, I haven’t tested in in a more recent version.
program wzerror
implicit none
integer, parameter :: N = 20
complex, target :: wz(N)
real, pointer :: wr(:)
integer :: i
wr => wz%re
wz = 0
wr = [(i,i=1,N)]
wr = wr + wz(N:1:-1)%re
print *, wr
end program
The error is in the expression wr = wr + wz(N:1:-1)%re where the compiler wouldn’t recognize that there is aliasing that should be permitted in an expression with pointers.
ifx doesn’t compile it at all.
What code does it break? In fact, it’s quite the opposite: without the de facto standard of [real_1, imag_1, real_2, imag_2, …] FFTW wouldn’t work, which in turn would break most, if not all, solid-state electronic structure codes. No compiler writer would ever do that. So even if the convention is not officially standardised, it makes no practical difference because the de facto standard is inviolable.