Generic call with rank mismatch

Help
1

Dear all,

does anybody know an easy/elegant way to “work around” the fact, that the F2018 standard does not allow rank mismatch in generic calls? In the program below, I have type specific F77 (could be also C) routines, which copy elements from one array to an other. They are simple type specific memcopy routines, expecting a source and target address and the nr. of elements to copy. They are rank-agnostic though. I have defined explicit interfaces for those routines and would like now to find a way, that the caller can call all of them using the same generic name.

Note, I can not change the F77 routines themselves, as they are part of a library out of our control (ScaLAPACK), which we develop some modern wrappers for (in ScalapackFx).

subroutine copy_integer(source, dest, nelem)
  integer, intent(in) :: source(*)
  integer, intent(out) :: dest(*)
  integer, intent(in) :: nelem
  print *, "Integer copy"
  dest(1:nelem) = source(1:nelem)
end subroutine copy_integer

subroutine copy_real(source, dest, nelem)
  real, intent(in) :: source(*)
  real, intent(out) :: dest(*)
  integer, intent(in) :: nelem
  print *, "Real copy"
  dest(1:nelem) = source(1:nelem)
end subroutine copy_real

module testmod
  implicit none

  interface copy

    subroutine copy_integer(source, dest, nelem)
      integer, intent(in) :: source(*)
      integer, intent(out) :: dest(*)
      integer, intent(in) :: nelem
    end subroutine copy_integer

    subroutine copy_real(source, dest, nelem)
      real, intent(in) :: source(*)
      real, intent(out) :: dest(*)
      integer, intent(in) :: nelem
    end subroutine copy_real

  end interface

end module testmod


program testprog
  use testmod, only : copy, copy_integer, copy_real
  implicit none

  integer :: array1(2, 2), array2(4, 3)

  array1(:,:) = reshape([1, 2, 3, 4], [2, 2])
  array2(:,:) = 0

  ! Specific call with rank mismatch works nicely
  call copy_integer(array1(1, 2), array2(2, 2), 2)
  print "(A)", "Array2"
  print "(4I4)", array2

  ! Generic call with rank-mismatch seems not to be allowed (due to C.9.6/5, I guess)
  ! Is there any workaround, so that the caller does not have to use the type-specific call?
  ! call copy(array1(1, 2), array2(2, 2), 2)

end program testprog
2

@aradi, per current Fortran standard (2018), unfortunately there is not a clean way to do this exclusively using a Fortran processor. That’s why the effort to include better generic facilities in Fortran 202X (cc: @everythingfunctional here as to whether such use cases will be covered in the current iteration of the Generics effort).

Presuming a companion C processor with your Fortran processor (not a bad one in practice on most platforms), any number of options do open up using the C interoperability facilities.

One option that is close to your current line of thinking is this: apply your generic interface approach but to wrappers that implement type-specific dummy arguments of assumed-rank arrays introduced in Fortran starting with the 2018 revision. That is, instead of a generic interface to current procedures (in ScalaPack, etc.) that employ assumed-size arrays in dummy arguments.

See an example below toward the integer type; elided are wrappers for other intrinsic types.

module testmod

  use, intrinsic :: iso_c_binding, only : c_loc, c_f_pointer
  implicit none

  interface

    subroutine copy_integer(source, dest, nelem)
      integer, intent(in) :: source(*)
      integer, intent(out) :: dest(*)
      integer :: nelem
    end subroutine copy_integer

    subroutine copy_real(source, dest, nelem)
      real, intent(in) :: source(*)
      real, intent(out) :: dest(*)
      integer, intent(in) :: nelem
    end subroutine copy_real

  end interface

  interface copy !<-- Generic interface to new wrapper subprograms that use C-interoperability
   module procedure copy_int
  end interface 

contains
  
   subroutine copy_int( src, dest )
      integer, intent(in), target :: src(..)
      integer, intent(inout), target :: dest(..)
      integer, pointer :: rhs(:)
      integer, pointer :: lhs(:)
      call c_f_pointer( cptr=c_loc(src), fptr=rhs, shape=[ size(src) ] )
      call c_f_pointer( cptr=c_loc(dest), fptr=lhs, shape=[ size(src) ] )
      call copy_integer( source=rhs, dest=lhs, nelem=size(src) )
   end subroutine 

end module testmod

Driver code:

program testprog
  use testmod, only : copy
  implicit none

  integer :: array1(2, 2), array2(4, 3)

  array1(:,:) = reshape([1, 2, 3, 4], [2, 2])
  array2 = 0
  call copy( array1, array2 )
  print "(A)", "Array1"
  print "(4I4)", array1
  print "(A)", "Array2"
  print "(4I4)", array2

end program testprog

Program response:

C:\temp>ifort /standard-semantics z.f90
Intel(R) Fortran Intel(R) 64 Compiler Classic for applications running on Intel(R) 64, Version 2021.5.0 Build 20211109_000000
Copyright (C) 1985-2021 Intel Corporation. All rights reserved.

Microsoft (R) Incremental Linker Version 14.30.30706.0
Copyright (C) Microsoft Corporation. All rights reserved.

-out:z.exe
-subsystem:console
z.obj

C:\temp>z.exe
Integer copy
Array1
1 2 3 4
Array2
1 2 3 4
0 0 0 0
0 0 0 0

C:\temp>

Edit: as pointed out correctly by @aradi, the dummy arguments in the wrapper require the TARGET attribute for standard conformance. The example has been updated.

3 Likes
3

Although I’m not sure if this is valid (in the standard), passing an array slice to dummy assumed-size arrays may also be another approach?

call copy(array1(1:1, 2), array2(2:2, 2), 2)
4

@septc, I think you want to do

call copy(array1(1:, 2), array(2:, 2), 2)

don’t you

5

@sept Yes, that seems to be valid Fortran, and this was our first approach indeed. However, I don’t see, how it can be warranted, that this passes the address of that selected element in the array, not the address of a temporary copy (of the wrong size), which the compiler in theory would be allowed to create.

6

@rwmsu Yes, I think, if we passed array subobjects, which contain at least as many elements as nelem, we may have a solution, which is warranted to work, even if the compiler makes temporary copies.

7

@rwmsu yes, thanks much!
call copy(array1(1:, 2), array(2:, 2), 2)

@aradi

how it can be warranted, that this passes the address of that selected element in the array,

This is also my concern… I guess if we remove all intent(in) and intent(out) from the dummy argument list (such that it is closer to F77 interfaces), compiler may be less likely to make a copy (as long as the actual argument is a contiguous array section), but not sure such a behavior is “guaranteed” by the standard.

8

@FortranFan Thank you, that is a very neat idea! I also experimented around a bit with wrappers using assumed rank dummy arguments, but I did not think about the trick to convert the dummy argument to a c_ptr and back with c_f_pointer before passing it to the F77 copy routine. Cool, thanks!

One remark only: I think, the dummy arguments src, and dest in copy_int need the target attribute, don’t they?

9

Yes indeed, per the standard the TARGET attribute is required though the onus is on the code author, a processor is not required to include a capability to detect this. Thanks much for pointing it out.

10

Before reading @FortranFan’s code I was about to suggest a very similar approach. You’re essentially doing explicitly what the assumed-size arguments of F77 did implicitly. Reinterpret the array as a 1D array in element order. The use of c_loc c_f_pointer is one way to do this. The (now obsolete) equivalence would be another way. I feel like you should be able to do this directly with Fortran pointers, but maybe there’s a constraint somewhere that prevents this.

At any rate, a rank-agnostic copy procedure is not currently one of our use cases for the generics facility, but I will bring up to subgroup the possibility of including it.

11

The standard places tight restrictions on where assumed-rank arrays, which have to be dummy arguments, can appear. It disallows them as a designator or in an expression.

And in the absence of assumed-rank arrays, the construction of the generic interface will either be prohibitively laborious or limiting.

Hence my suggestion above re: the use of the companion C processor.

12

Right. Is something like the following allowed?

integer, target :: rhs(a, b)
integer, pointer :: lhs(:)

lhs => rhs

If so, it should be possible for the standard to specify (perhaps for the limited case that the pointer is 1D) that the left hand side assumes size of product(size(rhs)), and then allow assumed rank variables on the right hand side. This would allow one to - completely on the Fortran side of things - emulate the behavior of assumed size arguments.

13

I wanted to play with assumed-rank limitations, and
so slightly modified the example; and got three different
results with three compilers.

Should gfortran enforce passing assumed rank to only
assumed rank for an external procedure?

module M_copy
implicit none

interface copy
   module procedure :: copy_integer
end interface

contains

subroutine copy_integer(lhs,rhs,count)
integer :: lhs(..),rhs(..)
integer,intent(in) :: count
   call copy_integer_flat(lhs,rhs,count)
end subroutine copy_integer

end module M_copy

subroutine copy_integer_flat(source, dest, nelem)
integer, intent(in) :: source(*)
integer             :: dest(*)
integer, intent(in) :: nelem
  dest(1:nelem) = source(1:nelem)
end subroutine copy_integer_flat

program testprog
use M_copy, only: copy
implicit none
integer :: array1(2, 2), array2(4, 3)
  array1(:,:) = reshape([1, 2, 3, 4], [2, 2])
  array2(:,:) = 0
  call copy(array1(1:, 2:), array2(2:, 2:), 2)
  print "(A)", "Array2"
  print "(4I4)", array2
end program testprog

$ gfortran testit.f90

   13 |    call copy_integer_flat(lhs,rhs,count)
      |                           1
Error: Assumed-rank argument requires an explicit interface at (1)

$ nvfortran testit.f90
$ ./a.out
Array2
   0   3   4   0
   0   0   0   0
   0   0   0   0

$ ifort testit.f90
$ ./a.out
Array2
   0   0   0   0
   0   3   4   0
   0   0   0   0
14

@urbanjost

gfortran 11 doesn’t give a compile error if you move copy_integer_flat in front of copy integer in your module. However, ifort 2021.5 gives the following errors.

test.f90(20): error #8779: If the actual argument is assumed rank, the corresponding dummy argument must also be assumed rank. [LHS]
call copy_integer_flat(lhs,rhs,count)
--------------------------^
test.f90(20): error #8779: If the actual argument is assumed rank, the corresponding dummy argument must also be assumed rank. [RHS]
call copy_integer_flat(lhs,rhs,count)

15

For practically every (if not every, I’ll have to check if there are any exceptions) modern Fortran facility starting with Fortran 90 revision, the standard stipulates explicit interfaces.

Additionally with assumed-rank arrays, the standard has a numbered constraint whereby, “An assumed-rank variable name shall … appear … as an actual argument that corresponds to a dummy argument that is assumed-rank …”

The code you have posted does not conform and gfortran is right to detect and report the error.

16

Do you think it would be sensible to add a contiguous attribute to src and dest?
Although it may not compile in general cases.

17

@egio That’s a very good catch! Without the contiguous attribute (or without an explicit is_contiguous() check), the example is very fragile. If the caller passes a non-contiguous region, it won’t work as expected.

For completeness, below the accordingly corrected example of @FortranFan (which also corrects a missing intent(in) in the interface):

subroutine copy_integer(source, dest, nelem)
  integer, intent(in) :: source(*)
  integer, intent(out) :: dest(*)
  integer, intent(in) :: nelem
  print *, "Integer copy"
  dest(1:nelem) = source(1:nelem)
end subroutine copy_integer

module testmod
  use, intrinsic :: iso_c_binding, only : c_loc, c_f_pointer
  implicit none
  
  interface
      subroutine copy_integer(source, dest, nelem)
      integer, intent(in) :: source(*)
      integer, intent(out) :: dest(*)
      integer, intent(in) :: nelem
    end subroutine copy_integer
  end interface
  
  interface copy
    module procedure copy_int
  end interface copy
  
contains
  
  subroutine copy_int(src, dest)
    integer, intent(in), contiguous, target :: src(..)
    integer, intent(inout), contiguous, target :: dest(..)
    integer, pointer :: rhs(:)
    integer, pointer :: lhs(:)
    call c_f_pointer(cptr=c_loc(src), fptr=rhs, shape=[size(src)])
    call c_f_pointer(cptr=c_loc(dest), fptr=lhs, shape=[size(src)])
    call copy_integer(source=rhs, dest=lhs, nelem=size(src))
  end subroutine copy_int
  
end module testmod


program testprog
  use testmod, only : copy
  implicit none

  integer :: array1(2, 2), array2(4, 3)

  array1(:,:) = reshape([1, 2, 3, 4], [2, 2])
  array2(:,:) = 0
  ! Copies elements [1, 3] at positions [[1, 1], [1, 2]] in array1 into
  ! positions [[1, 2], [3, 2]] in array2.
  call copy(array1(1:1, 1:2), array2(1:3:2, 2))
  print "(A)", "Array1"
  print "(4I4)", array1
  print "(A)", "Array2"
  print "(4I4)", array2

end program testprog

If you execute the example, you get:

Array2
   0   0   0   0
   1   0   3   0
   0   0   0   0

while without the contiguous attributes, I obtain:

Array2
   0   0   0   0
   1   2   0   0
   0   0   0   0

which is clearly wrong.

18

Thinking more about it, one could maybe avoid the necessary detour through c_loc(), if the standard allowed contiguous assumed rank arrays being target for 1D pointers with specified bounds:

  subroutine copy_int(src, dest)
    integer, intent(in), contiguous, target :: src(..)
    integer, intent(inout), contiguous, target :: dest(..)
    integer, pointer :: rhs(:)
    integer, pointer :: lhs(:)
    rhs(1:size(src)) => src
    lhs(1:size(dest)) => dest
    call copy_integer(source=rhs, dest=lhs, nelem=size(src))
  end subroutine copy_int

Apparently, that is not allowed currently, (at least, no compiler I’ve tried was willing to compile it). Is there a fundamental reason, why such assignments should be forbidden, or was it, that nobody thought about it so far? I’d find this much more fortranic than the solution with c_loc()

2 Likes
19

Just as a try I have written a C function that will copy an array into another, not yet optimized for contiguous array.

#include <stdio.h>
#include "ISO_Fortran_binding.h"

int linear_to_c_index(CFI_dim_t dim[], int rank, CFI_index_t index[], size_t i_seq)
{
    for (int i = 0; i < rank; i++)
    {
        index[i] = i_seq % dim[i].extent + dim[i].lower_bound;
        i_seq /= dim[i].extent;
    }
    return i_seq;
}
size_t fortran_size(CFI_cdesc_t *v)
{
    size_t n_element = 1;
    for (int i = 0; i < v->rank; i++)
    {
        n_element *= v->dim[i].extent;
    }
    return n_element;
}

void display_info(CFI_cdesc_t *src)
{
    for (int i = 0; i < src->rank; i++)
    {
        printf(" %d", i);
        printf(" lower_bound %d", src->dim[i].lower_bound);
        printf(" extent %d", src->dim[i].extent);
        printf("stride %d\n", src->dim[i].sm);
    }
    printf("\nelem_len %d\n", src->elem_len);
}
void copy_out(CFI_cdesc_t *src, CFI_cdesc_t *dest, size_t n, size_t *copied)
{
    CFI_index_t ind_src[CFI_MAX_RANK] = {0};
    CFI_index_t ind_dest[CFI_MAX_RANK] = {0};
    if (src->type != dest->type)
    {
        *copied = 0;
        return;
    }
    size_t n_element_src = fortran_size(src);
    size_t n_element_dest = fortran_size(dest);

    size_t n_actual = n_element_src <= n_element_dest ? n_element_src : n_element_dest;

    n_actual = n_actual <= n ? n_actual : n;
    *copied = n_actual;

    for (size_t i = 0; i < n_actual; i++)
    {
        linear_to_c_index(src->dim, src->rank, ind_src, i);
        linear_to_c_index(dest->dim, dest->rank, ind_dest, i);

        void *ptr_src = CFI_address(src, ind_src);
        void *ptr_dest = CFI_address(dest, ind_dest);
        memcpy(ptr_dest, ptr_src, src->elem_len);
    }
}

Can be called with the following interface:

interface
     subroutine copy_out(src, dest, n, copied) bind(C)
        use iso_c_binding
        integer(c_int), intent(in) :: src(..)
        integer(c_int), intent(inout) :: dest(..)
        integer(c_size_t), value :: n
        integer(c_size_t),intent(out) :: copied
    end subroutine 
end interface

One have also to add checks that elem_len is the same in both src and dest, and correctly manage the rank=0.

1 Like