@@ -35,20 +35,44 @@ at some point:
35
35
! ...
36
36
allocate( r1(n), r2(n), r3(n) )
37
37
! ... some computations on r1(:) and r2(:) as reals
38
- call rfft(r1)
39
- call rfft(r2)
40
- call fftconvol(r1,r2,r3,n/2) ! called without any interface
38
+ call rfft(r1) ; call rfft(r2)
39
+ call fftconvol(r1,r2,r3,n/2) ! called without any interface,
40
+ ! hence type mismatchs
41
41
call irfft(r3)
42
42
! ... some computations on r3(:) as real
43
43
end program
44
44
45
+ ! dangling routine (not contained and not in a module)
45
46
subroutine fftconvol(c1,c2,c3,n)
46
47
integer, intent(in) :: n
47
48
complex, intent(in) :: c1(n), c2(n)
48
49
complex, intent(out) :: c3(n)
49
50
c3(:) = c1(:) * c2(:)
50
51
end subroutine
51
52
53
+ With modern Fortran another trick has got quite popular, by using the
54
+ C interoperability features. Nonetheless the trick is non standard:
55
+
56
+ program foo
57
+ use iso_c_binding
58
+ real, allocatable :: r1(:), r2(:), r3(:)
59
+ integer, parameter :: n = 100000
60
+ complex, pointer :: c1(:), c2(:), c3(:)
61
+ ! ...
62
+ allocate( r1(n), r2(n), r3(n) )
63
+ ! ... some computations on r1(:) and r2(:) as reals
64
+ call rfft(r1) ; call rfft(r2)
65
+ ! non-standard trick
66
+ call c_f_pointer( c_loc(r1), c1, [n/2] )
67
+ call c_f_pointer( c_loc(r2), c2, [n/2] )
68
+ call c_f_pointer( c_loc(r3), c3, [n/2] )
69
+
70
+ c3(:) = c1(:) * c2(:)
71
+
72
+ call irfft(r3)
73
+ ! ... some computations on r3(:) as real
74
+ end program
75
+
52
76
53
77
3. Proposed solution
54
78
@@ -63,11 +87,8 @@ vice-versa. The above code would become:
63
87
! ...
64
88
allocate( r1(n), r2(n), r3(n) )
65
89
! ... some computations on r1(:) and r2(:) as reals
66
- call rfft(r1)
67
- call rfft(r2)
68
- c1 => r1
69
- c2 => r2
70
- c3 => r3
90
+ call rfft(r1) ; call rfft(r2)
91
+ c1 => r1 ; c2 => r2 ; c3 => r3
71
92
c3(:) = c1(:) * c2(:)
72
93
call irfft(r3)
73
94
! ... some computations on r3(:) as real
@@ -80,26 +101,26 @@ used, with additional rules and restrictions:
80
101
81
102
`c => r`
82
103
83
- - `r` is a *contiguous* real array, which has either the target or the
84
- pointer attribute
85
- - `c` is a complex array pointer of the same kind as `r`, and of the same
86
- rank by default (but pointer rank remapping can be used)
87
- - `c` can also be a complex scalar pointer, in the case where r is a
104
+ - `r` shall be a *contiguous* real array, which has either the target
105
+ or the pointer attribute
106
+ - `c` shall be a complex array pointer of the same kind as `r`, and of
107
+ the same rank by default (but pointer rank remapping can be used)
108
+ - `c` could also be a complex scalar pointer, in the case r is a
88
109
rank-1 array of size 2
89
110
- the size of the first dimension of `r` shall be even
90
- - `c%re` refers to the same storage as `r(1::2)` (rank-1), or
111
+ - `c%re` shall refer to the same storage as `r(1::2)` (rank-1), or
91
112
`r(1::2,:)` (rank-2), etc...
92
- - `c%im` refers to the same storage as `r(2::2)` (rank-1), or
113
+ - `c%im` shall refer to the same storage as `r(2::2)` (rank-1), or
93
114
`r(2::2,:)` (rank-2), etc...
94
115
95
116
`r => c`
96
117
97
118
- the exact opposite
98
- - `c` is a *contiguous* complex array or a complex scalar, which has
99
- either the target or the pointer attribute
100
- - `r` is a real array pointer of the same kind as `c`, and of the same
101
- rank by default (but pointer rank remapping can be used)
102
- - if `c` is a scalar, then `r` is a rank-1 pointer of size 2
119
+ - `c` shall be a *contiguous* complex array or a complex scalar, which
120
+ has either the target or the pointer attribute
121
+ - `r` shall be a real array pointer of the same kind as `c`, and of the
122
+ same rank by default (but pointer rank remapping can be used)
123
+ - if `c` is a scalar, then `r` shall be a rank-1 pointer of size 2
103
124
- same other rules as above
104
125
105
126
3.2 Alternative syntaxes
@@ -126,9 +147,16 @@ r => real :: c
126
147
```
127
148
128
149
Something generic (that is, not `complex_pointer()` or `real_pointer()`)
129
- may be desirable in case other inter-type pointer association would be
150
+ may be desirable in case other inter-type pointer associations would be
130
151
allowed in future versions of the standard.
131
152
153
+ 3.3 Prototyping, etc...
154
+
155
+ I think that this proposal doesn't need to have a preliminary prototype
156
+ implementation, as it essentially consists in standardizing an already
157
+ existing and common practice. A prototype implementation would do
158
+ nothing else than mimicing the `c_f_pointer()` trick.
159
+
132
160
133
161
4. Issues / Objections / Limitations
134
162
@@ -158,16 +186,18 @@ In this case, `c%re` would refer to the same storage as `r(1:n/2)`.
158
186
159
187
Allowing a real actual argument to be associated to a complex dummy
160
188
argument -and vice-versa- has also been considered, but it would raise
161
- backward compatibility issues. So this part has been dropped.
189
+ backward compatibility issues. So this part has been dropped from the
190
+ proposal.
162
191
163
192
4.3. Alignment
164
193
165
194
Considering for instance a default real type stored on 4 bytes, the
166
195
default complex type is stored on 8 bytes. Compilers may prefer/want to
167
- align complex arrays on 16 bytes boundaries, which cannot be guaranteed
196
+ align the complex arrays on 8 bytes boundaries, which cannot be guaranteed
168
197
if a complex pointer is associated to an arbitrary real array
169
198
(e.g. `c => r(2:)`). If this is a problem, the pointer association may be
170
- allowed only in the other way (real pointer associated to a complex array).
199
+ allowed only the other way (real pointer associated to a complex array),
200
+ where alignement should not be a problem.
171
201
172
202
173
203
5. References
0 commit comments