How to concisely declare a function result as a multi-dimensional automatic array?

I want to declare the function result to have the same shape as the input multi-dimensional array dummy argument. For rank 2, I do this:

function f(x)
  real, intent(in) :: x(:,:)
  real :: f(size(x,1),size(x,2)) ! <- more concise? how?

For higher ranks, it gets worse, as you can imagine.

Is this the best I can do?

elemental is not an option as I want this function to be pointed to by a procedure pointer.

Effectively yes, with the current standard - depending on the needs, the other intrinsics LBOUND, UBOUND may come in handy.

If I recall correctly though, with the next standard revision Fortran 2023 (? 202X?), the following likely becomes an option:

function f(x)
  real, intent(in) :: x(:,:)
  real :: f(shape(x))
1 Like

I asked myself a similar question a while ago. Ultimately, I decided to change my code to have a flattened array (i.e., a vector) as input and output. Of course, this approach has drawbacks, and it is not always applicable. Optionally, we can also internally reshape (and finally flatten again), to facilitate some operations.
Of course, this is not a solution to OP’s question, rather a workaround.

1 Like

Fortran needs to be extended for this. I was thinking something like this:

function f(x)
  dim :: n, m
  real, intent(in) :: x(n,m)
  real :: f(n,m)

Which has many advantages, consider a matmul, currently:

function matmul(A, B) result(C)
  real, intent(in) :: A(:,:), B(:,:)
  real :: C(size(A,1),size(B,2))

Note that you can’t use shape easily for this and note that dimensions are not tied, so the compiler cannot check that A and B are conforming. Now after the extension:

function matmul(A, B) result(C)
  dim :: m, n, k
  real, intent(in) :: A(m,n), B(n,k)
  real :: C(m, k)

The dim behaves just like A(:), in other words, the array still gets passed in with a descriptor. But you can now tie the different dimensions together and the compiler can check (at compile time and/or runtime) that the arrays A, and B are conforming with their dimensions, so you do not need to do this check by hand inside the subroutine. And it is consistent if you want to use the old fashioned arrays:

function matmul(m, n, k, A, B) result(C)
  integer, intent(in) :: m, n, k
  real, intent(in) :: A(m,n), B(n,k)
  real :: C(m, k)

Note 1: instead of dim :: m, n, k, we could use integer, infer :: m, n, k, to be consistent with integer, intent(in) :: m, n, k.

Note 2: real, intent(in) :: A(:,:), B(:,:) would be equivalent to integer, infer :: m, n, k, l; real, intent(in) :: A(m, n), B(k,l)

2 Likes

A slight different syntax could be:

function matmul(A, B) result(C)
  dim :: m, n, k
  real, intent(in) :: A(m=>:, n=>:), B(n=>:,k=>:)
  real :: C(m, k)

That is stressing the fact that A and B are assumed shape arrays. Even though is more verbose.
The first dim declaration could be even unnecessary.

There is a problem in situation like this one:

subroutine strange(a, b, c)
real :: a(n=>:), b(n=>:), c(n=>:) 

As the values can be inconsistent at run time.
One possibility is to introduce a new intrinsic function, let call it consistent that has as argument this dim values and return .TRUE. or .FALSE. if they are consistent or not.
Then then the programmer can error stop or do whatever he/she likes. With the assumption that, if one of the dim variable is inconsistent it cannot be used anywhere (almost), but the processor is not asked to report or detect it.
Basically the following:

subroutine strange(a, b, c)
real :: a(n=>:), b(n=>:), c(n=>:)

if (.not. consistent(n)) error stop 

will be equivalent

subroutine strange(a, b, c)
real :: a(:), b(:), c(:)

associate(n=> size(a))
if (.not. (n==size(b) .and. n==size(c)) error stop 

end associate

If you like it it could be transformed in a proposal for Fortran 202z.

1 Like

@egio excellent ideas! Thank you. I really feel there is something here. Btw, once we have something that more people (than just me :slight_smile: ) like, this would be easy to implement in LFortran as an optional extension (for now) and we can start using it and see if we actually like it. So that the actual proposal to the standard committee will be very robust.

Couple questions.

Does it make sense to require to define “n” to be an integer?

subroutine strange(a, b, c)
integer, dim :: n                   ! Possibly without the "dim"
real :: a(n=>:), b(n=>:), c(n=>:) 

Regarding:

subroutine strange(a, b, c)
real :: a(n=>:), b(n=>:), c(n=>:) 

That seems analogous to:

subroutine strange(n, a, b, c)
integer :: n
real :: a(n), b(n), c(n) 

The only difference is that in the first case “n” is inferred from “a”, while in the second case it is passed in.

In both cases, the compiler should check this using bounds checks at runtime automatically in Debug mode. In Release mode no checks will be done for performance.

So there is no issue from my point of view.

I don’t know.
Let’s say that if not declared it will be of the kind of a standard integer, while if declared it could be of a different kind, that may be useful in some situation.

I will keep the dim.

Your suggestion about including it in LFortran is wonderful!
:slight_smile:

1 Like

With assumed shape declarations (:), the dummy argument need not be contiguous. With explicit shape arguments, the dummy argument must be contiguous. That means in some cases, copy-in/copy-out is required by the compiler. That semantic difference has performance consequences.

As for this proposal in general, as described above it is basically just syntactic sugar that achieves exactly the same thing as the currentl syntax using size(). If you want to make the proposal appeal to a wider range of programmers, you need to add some useful new capabilities. The ability to ensure/assert consistency of dimensions over several arrays, or even within a single array (e.g. square matrices), might do this. That would be most useful if it could be done at compile time. However, this assertion most likely could be done only at run time for many cases, e.g. involving separate compilation of caller and callee, in which case we are back to just syntactic sugar for what the programmer can already do in a straightforward way.

Most of the time when I use assumed shape array I add a comment after the array with the name of the dimension just to remind me which are:

subroutine alpha( a, b)
  real :: a(:,:) ! m, n
  real :: b(:,:,:) ! k, n, m
  integer :: m,n,k

  m = size(a,1)
  n = size(a,2)
  k = size(b,1)
  if (m/= size(b,3) .or. n /= size(b, 2) ) error stop

Well I think it is definitely error prone compared to my proposal

subroutine alpha( a, b)
  real :: a(m=>:,n=>:) 
  real :: b(k=>:,n=>:,m=>:)

  if (.not. consistent(m,n) ) error stop
1 Like

You could make the function result allocatable and then write

allocate(f,mold=x)

3 Likes

I was thinking the same thing. Is there any practical, observable difference between an automatic and allocatable function result, other than perhaps where it is allocated (stack or heap)?

3 Likes

Just to make it clear, the proposal above is not just syntactic sugar, the added value is that the compiler will check (at runtime) that the array sizes are compatible, which is not possible to do today, e.g. it is not possible today for the compiler to check automatically for you at runtime that the following arrays A and B are compatible with the dimension “n”:

function matmul(A, B) result(C)
  dim :: m, n, k
  real, intent(in) :: A(m,n), B(n,k)
  real :: C(m, k)

The other added value is that the code is self-documenting in terms of which dimensions must be compatible.

Regarding non-contiguous arrays, indeed, it would allow it.

(The code is also more readable and concise, which is just “syntactic sugar”, but that’s a good thing in my opinion.)

1 Like

Ah, I certainly hope so! That would be so convenient.