fortran - 对不同形状的数组使用单个虚拟参数

Question

假设我们有几个不同的连续数组

real :: a1(2),a2(2),a3(2)
real :: b1(2,2), b2(2,2),b3(2,2)
real :: c1(3,3,3), c2(3,3,3),c3(3,3,3)

和一个子程序

subroutine mtpy(x,y,z)
  real,contiguous, intent(in) :: x(:), y(:)
  real, intent(out) :: z(size(x))

  z=x*y
end subroutine mtpy

如何在以下一系列调用中使用 mtpy:

call mtpy(a1,a2,a3)
call mtpy(b1,b2,b3)
call mtpy(c1,c2,c3)

显然，这会导致编译器错误，因为实际参数和虚拟参数的形状不匹配。在这样的情况下，我曾经声明几个特定的​​过程，每个过程处理一个特定的形状，然后使用一个接口(interface)包装所有这些过程。然而，这是相当乏味的(想象一下，有大量的简单基本函数和纯过程，将多维数组(最多三个维度)视为单维数组，然后为每个数组提供诸如 sub_1d、sub_2d、sub_3d 等实例。尽管他们实际上都在做同样的工作)。

我认为，部分解决方案是使用 RESHAPE

call mtpy(reshape(b1,[4]),reshape(b2,[4]),bb)

但是，我可以确定编译器(我主要对 gfortran 和 ifort 感兴趣)不会开始创建 1d 临时变量来保存 reshape 的 b1 和 b2 数组吗？

现在，我也知道可以声明一个数组指针，例如

real, pointer, contiguous :: p1(:),p2(:),p3(:)

并进行以下指针赋值，例如

p1(1:size(c1))=>c1

但是，这种方法有一个缺点，我需要将原始数组声明为目标。这不会影响编译器能够执行的优化吗？

我认为，另一个解决方案是使用假定大小的数组，但我注意到 Metcalf 等人称其用法“已弃用”，并且我再次不确定对优化的影响。

那么，是否有一种简单的方法可以将多维 Fortran 数组视为单维数组(在子例程或函数中)，并且不会对该数组强加不必要的假设(例如 TARGET)？如果我可以使用 RESHAPE 而不必担心创建临时数组(我只处理连续的数组)，我会这样做。有什么建议吗？

Answer 1

future 的 Fortran 2018 标准将提供假定等级数组(允许接收任何等级的数组)和一个选择等级构造，该构造将轻松地允许使用一个假定来解决类似的情况 -排名数组(例如，参见 The new features of Fortran 2018 从第 16 页开始)，并且使用多个假定排名数组会更加困难。

假定大小的数组虽然不流行或不推荐，但它们是当前标准(Fortran 2008)以及下一个标准草案(Fortran 2018)的有效且未过时的功能，因此可以在需要时使用它们。由于许多 Fortran 77 代码都依赖于此，并且其中大部分已有数十年历史，因此我希望它在大多数编译器中得到显着优化。

但是，您不需要使用假定大小的数组，您可以使用显式形状数组(具有显式维度的数组)，并且只要数组实际参数有足够的元素，代码就有效，因为，根据2008年标准第12.5.2.11节第4段，

An actual argument that represents an element sequence and corresponds to a dummy argument that is an array is sequence associated with the dummy argument if the dummy argument is an explicit-shape or assumed-size array. The rank and shape of the actual argument need not agree with the rank and shape of the dummy argument, but the number of elements in the dummy argument shall not exceed the number of elements in the element sequence of the actual argument. If the dummy argument is assumed-size, the number of elements in the dummy argument is exactly the number of elements in the element sequence.

所以你可以

call mtpy(a1,a2,a3,size(a3))
call mtpy(b1,b2,b3,size(b3))
call mtpy(c1,c2,c3,size(c3))
...
subroutine mtpy(x,y,z,n)
  integer, intent(in) :: n
  real, intent(in) :: x(n), y(n)
  real, intent(out) :: z(n)
  z=x*y
end subroutine mtpy