Examples for the new @ rank agnostic array notation

At one of the last J3 meetings, the following feature has passed:

Rank-agnostic array element and section denotation · Issue #157 · j3-fortran/fortran_proposals · GitHub

The issue has some examples. What other examples would you use this in? I would like to see several diverse examples so that we can write a tutorial on this feature, as well as see if we need to submit some proposals to make this feature more useful.

For example, right now it would allow to write code like this (see the issue for context):

subroutine max_at(x, q1, q2, q3, vec)
real, intent(in) :: x(..)
real, intent(in) :: q1(..)
real, intent(in) :: q2(..)
real, intent(in) :: q3(..)
real, intent(out) :: vec(3)
integer, allocatable :: idx(:)
select rank(x)
rank(2)
    idx = maxloc(x)
    vec = [q1(@idx), q2(@idx), q3(@idx)]                             
rank(3)
    idx = maxloc(x)
    vec = [q1(@idx), q2(@idx), q3(@idx)]
rank(4)
    idx = maxloc(x)
    vec = [q1(@idx), q2(@idx), q3(@idx)]
end select
end subroutine

Wouldn’t it make sense to submit a proposal to allow writing it like this?

subroutine max_at(x, q1, q2, q3, vec)
real, intent(in) :: x(..)
real, intent(in) :: q1(..)
real, intent(in) :: q2(..)
real, intent(in) :: q3(..)
real, intent(out) :: vec(3)
integer :: idx(rank(x))
idx = maxloc(x)
vec = [q1(@idx), q2(@idx), q3(@idx)]                             
end subroutine

(No allocatable, thus faster, and no select rank construct, thus shorter.)

1 Like

Note that the syntax could be just q1(idx) (no @).

This feature needs a lot more work, IMO.

2 Likes

I would also prefer if it was just q1(idx), but I think this syntax is already taken, I believe currently
q1(idx) is equivalent to [q1(idx(1)), q1(idx(2)), ..., q1(idx(size(idx)))].

I also feel this feature needs a lot more work.

q1(idx) is a vector-valued subscript applied to an assumed-rank array, which is not valid today, so it can be given a meaning without invalidating any extant code.

2 Likes

I also thought of suggesting this, but then x(idx) would have a very different meaning depending on the rank of x. If x is 1D, x(idx) returns the number of elements in idx, but if x is >1D, x(idx) would return a single value. Maybe this inconsistency is worth it for a simpler syntax.

1 Like

Unless the committee passes another paper that changes it, yes @ is already determined and it will be part of 202X. However, if we don’t like it, we can propose a paper that can change the behavior.

(For the record, I voted no on this paper, I would like to see it more developed. But the majority disagreed with me.)

1 Like

This feature will be more useful if rank is allowed to be allocatable. For example, it can be used for generalizing ‘reduce/fold’ functions as shown below. (I have left some details for brevity.)

pure function general_reduce(binary_op,x,dim) result (r)
  real   ,intent(in) :: x(..)
  integer,intent(in) :: dim
  real,allocatable  :: r(..)
  integer :: i, x_shape(RANK(x)), idx(RANK(x)), idr(RANK(x)-1)
  interface
    real pure function binary_op(x,y)
      real,intent(in) :: x, y
    end function binary_op
  end interface

  x_shape = SHAPE(x)
  r = ALLOCATE(r([x_shape(:dim-1),x_shape(dim+1:)]))
  do i = 1,PRODUCT(x_shape)
    idx = counter_to_tuple(i)
    idr = [idx(:dim-1),idx(dim+1:)]
    if(idx(dim)==1) then
      r(@idr) = x(@idx)
    else
      r(@idr) = binary_op(r(@idr),x(@idx))
    end if
  end do
contains
  pure function counter_to_tuple(n) result (t)
  ! find index tuple from 1D counter. e.g. 1->(1,1,1), 2->(2,1,1) etc
    integer,intent(in) :: n
    integer :: t(RANK(x))
    ...
  end function counter_to_tuple
end function general_reduce

One can still achieve this Fortran by:
a) Using PACK and RESHAPE. But, this will require memory duplication.
b) Defining an array type. This can be as close to the proposed syntax as x.AT.[1,2,1] etc.

1 Like