Is it appropriate to think of rank as a compile time concept, rather than a run time concept? And that is why you can’t create a function that returns an array with a rank that depends on a runtime parameter. That is, you can’t write something like the reshape intrinsic in Fortran itself.
Mostly yes, but Fortran procedures can have assumed-rank array arguments whose rank can be queried at runtime using the select rank construct.
2 Likes
Is the question specific to Fortran, or more general? In a compiled (like Fortran) AND dynamically typed (unlike Fortran) langage, having runtimes varying rank would be possible
This is specific to Fortran, since I am trying to learn the language.
Is it possible to pass an array to h whose rank depends on a run time value?
program assumed_rank
implicit none
integer :: i, j
integer, allocatable :: shp(:), a(:)
read *, i
allocate(shp(i))
shp = [(j, j=2,i+1)]
allocate(a(product(shp)))
a = 1
call h(reshape(a, shp))
contains
subroutine h(r)
integer, intent(in) :: r(..)
select rank(r)
rank(1)
print *, "vector"
rank(2)
print *, "matrix"
end select
end subroutine h
end program assumed_rank
The program above complains that:
assumed_rank.f90(16): error #7959: The ‘SHAPE’ argument of intrinsic function ‘RESHAPE’ shall be of type integer, rank one, and constant positive size less than 8. [SHP]
call h(reshape(a, shp))
----------------------^
…which I’m thinking means that the size of SHAPE argument need to be known at compile-time (which implies that the rank is fixed at compile time). Maybe there is another intrinsic other than reshape that would return an array whose rank isn’t known until run time?
When it comes to Fortran and its standard document which influences a lot of the language, the concept is whether something is declared or dynamic. That’s: compile-time or run-time are not part of the “official” jargon.
So, for example, you will know with an object of CLASS(foo-t), the declared type is “foo_t” (what you might view as compile-time type) but in a given execution statement, the dynamic type may be any extended type of “foo_t” (during run-time as you will see it).
There is no such distinction with RANK and KIND, the other two that are part of TKR triumvirate and established distinguishable characteristics of objects in Fortran, particularly when it comes to authoring generic methods and containers.
This is a crucial limitation in the language standard the work toward enhanced Generics starting with Fortran 202Y must address but thus far I am yet to notice a sufficiently compact and usable solution for practitioners that is put forth by the subgroup working on this feature.
Until such a time, many “generic” intrinsics (FINDLOC is another one such case) can be implemented via “compiler magic” (as one may see them) but not by users themselves which is most unfortunate when it comes to the advancement of Fortran as the much-needed lingua franca of scientific and technical computing.
1 Like
Thanks. I came across this quote:
Data objects may be dynamic in size, shape, type, or length type parameters, but not rank or kind type parameters.
…which looks like it comes from the book The Fortran 2003 Handbook: The Complete Syntax, Features and Procedures. With a quick search, I could not find something like that in the 2023 draft standard, but there is probably something like that in there(?). And that is likely a historical feature for performance reasons? It seems like (with my extremely limited Fortran knowledge), you could create a custom data type that for all appearances, looked like a first-class “dynamic”-rank array, which would really be wrapper functions around a rank-1 vector.
At the implementation level, it seems like it might almost be more work for compiler writers to keep rank separated from being “dynamic”, since I’m assuming that there is a header of sorts that has the rank and bounds for each array(?).
In your code, the line
allocate(a(product(shp))
has mismatched parentheses, and the type in
real, intent(in) :: r(..)
does not match the caller, where an integer array is passed. Changing the code to
program assumed_rank
implicit none
integer :: i, j
integer, allocatable :: shp(:), a(:)
read *, i
allocate(shp(i))
shp = [(j, j=2,i+1)]
allocate(a(product(shp)))
a = 1
call h(reshape(a, shp))
contains
subroutine h(r)
integer, intent(in) :: r(..)
select rank(r)
rank(1)
print *, "vector"
rank(2)
print *, "matrix"
end select
end subroutine h
end program assumed_rank
gfortran says
xxassumed_rank.f90:11:22:
11 | call h(reshape(a, shp))
| 1
Error: 'shape' argument of 'reshape' intrinsic at (1) must be an array of constant size
and ifort says
xxassumed_rank.f90(11): error #7959: The 'SHAPE' argument of intrinsic function 'RESHAPE' shall be of type integer, rank one, and constant positive size less than 8. [SHP]
call h(reshape(a, shp))
----------------------^
compilation aborted for xxassumed_rank.f90 (code 1)
which agrees with your conjecture that
the size of SHAPE argument need to be known at compile-time (which implies that the rank is fixed at compile time).
1 Like
Whoops, thanks for catching those errors. Too much copy/paste in creating the minimal example. Thanks!
In current Fortran, rank is a compile-time concept. As others have said, the assumed-rank feature lets you declare a dummy argument that takes its rank from the actual argument, but once inside the procedure you can’t do much with it unless you use a SELECT RANK construct, where within each select block the array has a fixed rank.
In Fortran 2023, you can declare an array whose rank is taken from that of another array using the RANK clause, where the rank is a constant expression, thus it is still known at compile-time.
3 Likes
With a quick glance at the gfortran implementation, I see the reshape intrinsic:
typedef GFC_FULL_ARRAY_DESCRIPTOR(GFC_MAX_DIMENSIONS, char) parray;
static void
reshape_internal (parray *ret, parray *source, shape_type *shape,
parray *pad, shape_type *order, index_type size)
…where it looks like arrays are stored as an instance of (from libgfortran.h):
#define GFC_FULL_ARRAY_DESCRIPTOR(r, type) \
struct {\
type *base_addr;\
size_t offset;\
dtype_type dtype;\
index_type span;\
descriptor_dimension dim[r];\
…with the dimension/rank being stored as an array of upper bound, lower bound, and stride for each dimension:
typedef struct descriptor_dimension
{
index_type _stride;
index_type lower_bound;
index_type _ubound;
}
descriptor_dimension;
…so it might be the case that rank information is available to the run-time system.
Any time you move information from compile-time to run-time, you take a significant performance hit. I don’t think many of us want Fortran to turn into an interpreted language.
2 Likes
RE performance, I also remember this comment by Brad Chamberlain (the lead developer of Chapel in HPE) about compile-time vs run-time rank, though the comment is 6-years ago and so not very sure about the latest situation… (I guess probably the same).
I also wonder how Julia manages the rank of an array for efficiency (but there may be no problem as long as the JIT compiler can understand the rank for a given piece of code).
RE the rank of Fortran arrays, I use array pointers of different ranks (1D vs 2D) so that I can access the same data (in a derived type) in either the linear or Cartesian way, using the rank remapping feature. This is convenient for me because the same data can be accessed linearly (e.g. for simple model systems) or Cartesian (xyz + n-particles) for atomic systems.
So, not only is Fortran strongly-typed, it is strongly-kound and strongly-ranked as well.
What other property of a named Fortran entity is never allowed to vary during execution?
By the way, “arity” is a better name than “rank” (which has a very different meaning for matrices, e.g. What Is a Rank-Revealing Factorization? – Nick Higham), but I guess it’s too late for that.