I am working with some very complex, and rather old fortran code (most of it was written in fortran77 but some was earlier). I need to have the memory dynamically allocated and am therefore converting "COMMON" block statements to modules with use statements.
I have come to a problem where the rank of a given variable within the common block changes between routines (the overall memory footprint is the same however).
A highly simplified example of this would be
Subroutine sub1
COMMON /cblock/array(L,N,M)
!code using the array in a 3D manner
End
Subroutine sub2
COMMON /cblock/array(L*N*M)
!code using the array in a 1D manner
End
As you can see the array is of rank 3 in the first routine and of rank 1 in the second routine.
I have searched high and low for a method of making this function with allocatable arrays, but have not been successful.
It seems like with all the additional functionality built into the modern fortran standards there should be a way to duplicate what the common blocks did in an instance like this.
Thank you.
P.S. L,N,M would have been compile time parameters defined using include statements but will now be run time variable integers that are defined and used to allocate the arrays before sub1 or sub 2 would be called. I have sucsessfully converted a number of other common blocks to the modules and understand that portion of the process. It is only where these variables are treated as different rank that are giving me problems.
链接已复制
14 回复数
You could do something like this:
[fortran]
subroutine sub2
use mymod
real, pointer :: array1d(:)
array1d(L*M*N) => array3d
[/fortran]
You'll need to give array3d the TARGET attribute. Also you'll need the 2013 compiler for this.
Tim,
Using 1D array addressing was used in 70's and 80's to improve performance, as manageing the array subscript increment was typically faster than what the compiler could achieve. With modern compilers, optimisation makes this unnecessary.
You need to be careful with the integer arithmetic, as the actual index for the 1D array might be more complex than it first appears.
Converting to a 3D array would be more robust for further development.
.
I am surprised by the definition "COMMON /cblock/array(L*N*M) " as this was not available in most old compilers.
.
John
I sympathize with Tim because I too often used this old "trick." Even more often I would change the rank of arrays via argument passing. Apparently this is now illegal also, or at least obsolent and will be in the future. Current advice I can find regarding workarounds don't really suffice for me, at least given my limited understanding of things.
There were two reasons for changing the rank of matching arrays: (1) to intentionally convert 2-D susbscript management into 1-D ("inline") management, which is faster (used to be important, maybe not now?), simpler, and in places more intuitive. (2) to enable a sort of internal dynamic memory allocation, whereby a temporary array, regardless of rank, could be temporarily "equivalenced" to some point in the middle of a very large 1-D array. (I personally feel there were some real advantages to doing it this way, as opposed to true dynamic allocation, but I'll leave that for another time.)
Sometimes it was convenient to use EQUIVALENCE to accomplish the same thing, but of course that is now frowned upon too.
In short, I feel that all these modern rules and "enhancements" have actually eliminated some functionality, and at the very least, cause some real headaches in upgrading many old perfectly valid programs to fit the modern mantra.
I never thought of using pointers to accomplish this. But, if this is the best recommendation for fixing things up, and if 2013 is needed, where does that leave those of us that are stuck at say version 2011?
Steve,
This method looks promising. Is this type of pointer specific to the intel compiler or will it work with others as well (I am trying to make the code as portable as possible and avoid using compiler specific code)?
Also I have tested this out on a simple case and it seems to work for me. But here is the real question.
My common block is not as simple as a single array. It is instead a mix of anywhere from 6 to 20 arrays of various ranks and sizes (different declarations of the common block have differing numbers of arrays and different ranks of arrays, and there are 12 unique uses of the same common block, with dozens of overall uses). To make matters worse each array does not necessarily have the same memory footprint. An extremely simplified example would be:
/cblock/array1(n,j+1,k,L),array2(k),array3(n,j+1),array4(k,l,n,m,o,j),array5(i,i,n,j),rjnkm(1)
in one routine
/cblock/array6(n,k,L),array7(k,n,i),array8(l),array9(k,l,o,j),array10(i,j),rjnkm(n*k*l*j+...+i*n*j)
in another
The array rjnkm is a "junk" memory placeholder added at the end of each common block version to ensure that each one had the same overall memory footprint (the length of the rjnkm array is determined by a very mathematically simple yet drawn out calculation which I have only alluded to in the above example).
From my understanding a common block acts as a pointer to the first memory location of a block of contiguous memory which is accessed in each routine based on the variables declared in the common block declaration of that routine. Thus in the above example array7 in the second usage would really overlap a lot of the memory of array 1 from the first usage. Luckily in this case they at least used all of the memory as reals and not some real and some integers and some logicals (I hate to think what is really going on with the memory in that case).
From what I can tell the code was written in this manner to reuse the memory for machines that had limited resources. I understand this is horrible coding practice but it would be nye impossible to correct the code at this point. I want to mimic the old common block behavior so as to avoid any unintended changes to the code behavior.
Based on your previous response it seems like an option might be to define a 1D array in the module, lets call it Rmem(:)
and I allocate it to be of the length of the longest of any of the common blocks without the rjnkm array included.
in this example it would be of length (n*(j+1)*k*L+k+n*(j+1)+(k*l*n*m*o*j)+(i*i*n*j)) (you should see what it really looks like with the non-simplified example)
then in each routine I carve up the memory according to how the usage is in that common block
i.e. in the first routine I have pointers
real, pointer:: array1(:,:,:,:),array2(:),array3(:,:),array4(:,:,:,:,:,:),array5(:,:,:,:)
array1(1:n,1:j+1,1:k,1:L)=>Rmem(1:n*(j+1)*k*L)
array2(1:k)=>Rmem(n*(j+1)*k*L+1:n*(j+1)*k*L+k)
array3(1:n,1:j+1)=>Rmem(n*(j+1)*k*L+k+1:n*(j+1)*k*L+k+n*(j+1))
array4(1:k,1:l,1:n,1:m,1:o,1:j)=>Rmem(n*(j+1)*k*L+k+n*(j+1)+1:n*(j+1)*k*L+k+n*(j+1)+(k*l*n*m*o*j))
array5(1:i,1:i,1:n,1:j)=>Rmem(n*(j+1)*k*L+k+n*(j+1)+(k*l*n*m*o*j)+1:n*(j+1)*k*L+k+n*(j+1)+(k*l*n*m*o*j)+(i*i*n*j))
and the other routine would have different pointers and limits.
1 Would this even work in theory (it will take me days to implement if so)?
2 Is this really the best way to mimic the old behavior? It seems like there should be some workaround that allows each routine to work off of a single memory pointer like with the common blocks and let the routine determine the starting index and shape/size of each array.
Thanks again,
~Tim
dboggs,
To expand on Steve's suggestion:
module fooData
real, allocatable, target :: array3d(:,:,:)
...
end module fooData
-------------------- new file-----------
module as_3D
use fooData
real, pointer :: array(:,:,:)
end module as_3D
-------------------- new file-----------
module as_1D
use fooData
real, pointer :: array(:)
end module as_1D
-------------------- new file-----------
At program start you will need to initialize the two pointers (or more pointers)
Then in source files include the module in the representation you desire.
The above uses the same name for the array.
I think it would be clearer to use different names for the pointers. e.g.
arrayAs1D and arrayAs2D
or
type As
real, pointer :: As1D(:)
real, pointer :: As2D(:,:)
real, pointer :: As3D(:,:,:)
end type As
...
type(As) :: array
...
array%As1D(i) = ...
I will leave it to you as to fill in the allocation and pointer assignments (hint, contains in type or module).
Jim Dempsey
The method I describe is standard Fortran 2003. Whether all compilers you want to use implemented this feature, I don't know.
The Fortran standard also supports "Sequence Association" where you pass a single element of an array to a routine where the corresponding dummy argument is declared with a different rank. All compilers you use should support this.
Steve,
I take it you are referring to something like
call sub(x(1,4,7))
subroutine sub(dummyarray)
real dummyarray(i,j,k)
so now the sub is using memory starting at the location of x(1,4,7) and assuming that it has a shape of (i,j,k)
if not please expound on "Sequence Association"
thanks again,
~tim
I'm very glad to hear a recommendation for this technique, for it means that a good share of my old practice will still work. But whether it is actually advised, or will continue to work, is still questionable. Here are the rules for argument association that I find in Intel's manual:
A scalar dummy argument can be associated with only a scalar actual argument.
The type and kind parameters, and rank of the actual argument must match those of its associated dummy argument.
I find the same rules (of course) in two other textbooks. To me this means that a single element of an actual array argument cannot be associated with a array dummy argument. Is there some fine point of interpretation that I'm missing?
The text you cite is correct, but is the opposite direction of what we're talking about here, which is an array dummy argument associated with a scalar array element. Here is the text from the standard:
29 12.5.2.11 Sequence association
30 1 An actual argument represents an element sequence if it is an array expression, an array element designator, a
31 default character scalar, or a scalar of type character with the C character kind (15.2.2). If the actual argument is
32 an array expression, the element sequence consists of the elements in array element order. If the actual argument
33 is an array element designator, the element sequence consists of that array element and each element that follows
34 it in array element order.
35 2 If the actual argument is default character or of type character with the C character kind, and is an array
36 expression, array element, or array element substring designator, the element sequence consists of the storage
37 units beginning with the first storage unit of the actual argument and continuing to the end of the array. The
38 storage units of an array element substring designator are viewed as array elements consisting of consecutive
39 groups of storage units having the character length of the dummy array.
40 3 If the actual argument is default character or of type character with the C character kind, and is a scalar that is
41 not an array element or array element substring designator, the element sequence consists of the storage units of
42 the actual argument.
NOTE 12.32
Some of the elements in the element sequence may consist of storage units from different elements of the
original array.
1 4 An actual argument that represents an element sequence and corresponds to a dummy argument that is an array
2 is sequence associated with the dummy argument if the dummy argument is an explicit-shape or assumed-size
3 array. The rank and shape of the actual argument need not agree with the rank and shape of the dummy
4 argument, but the number of elements in the dummy argument shall not exceed the number of elements in the
5 element sequence of the actual argument. If the dummy argument is assumed-size, the number of elements in the
6 dummy argument is exactly the number of elements in the element sequence.
Stated elsewhere in the standard is the restriction that sequence association does not apply if the actual argument is a pointer. See also Doctor Fortran in "I've Come Here For An Argument"
Thanks very much for this explanation, Steve. It makes perfect sense and is good news to me. The important concept is that an actual argument that is an element of an array (which will most often be just A(1)) is called a "sequence" argument and is not at all the same thing as a "scalar" argument. What seems strange to me--given the vast importance of this concept--is that a sequence is not mentioned in the book and manual reading I do (and I do a lot of it). I never heard of a sequence argument until now!
In Fortran 77, subroutine arguments were called by reference, in that they transferred a memory address. You were allowed (able) to manage the memory address how you liked. This often involved dividing a storage array for multiple uses, eg
.
! main program
real*8 a(1000000)
...
n3 = n2 + length ! a calculated position in A
call sub (A(n3), n,m)
...
Subroutine Sub ( array, n,m)
real*4 array(n,m)
.
Most of this functionality has been replaced by ALLOCATE, except for the resize option, where you initially provided a very large array (all that was available) then resized it to what was required after reading what was needed.
Mixing array sizes and types is now non-standard, but still essential for these old codes.
The standard has moved away from those approaches to Fortran use, and probably why us formula translaters have moded away from modern Fortran.
.
John
Ian,
.
There were many significant codes written in the 70's, which a lot of programmers, like me, used and copied at that time.
Their legacy is their contribution to development of computational methods in many fields, not their being described as non-standard.
I don't see some of the recent changes to the Fortran standard contributing further to computational methods in these areas where I work.
.
John
