多维数组#

多维数组以列主序存储。这意味着最左边的(最内层的)索引连续寻址元素。从实践的角度来看,这意味着数组切片 V(:, 1) 是连续的,而切片 V(1, :) 中元素之间的步长是列的维度。当将数组切片传递给期望处理连续数据的过程时,这一点很重要。

内存的局部性对于根据您的应用程序进行考虑非常重要,通常在对多维数组执行操作时,顺序访问应始终以单位步长递增。

在下面的示例中,评估了两组点之间的反距离。请注意,这些点在数组 xyz1/xyz2 中连续存储,而最内层循环正在递增矩阵 a 的最左侧索引。

subroutine coulomb_matrix(xyz1, xyz2, a)
  real(dp), intent(in) :: xyz1(:, :)
  real(dp), intent(in) :: xyz2(:, :)
  real(dp), intent(out) :: a(:, :)
  integer :: i, j
  do i = 1, size(a, 2)
    do j = 1, size(a, 1)
      a(j, i) = 1.0_dp/norm2(xyz1(:, j) - xyz2(:, i))
    end do
  end do
end subroutine coulomb_matrix

另一个示例是三阶数组的第三维的收缩

do i = 1, size(amat, 3)
  do j = 1, size(amat, 2)
    do k = 1, size(amat, 1)
      cmat(k, j) = cmat(k, j) + amat(k, j, i) * bvec(i)
    end do
  end do
end do

连续数组切片可用于数组边界重新映射,以允许将更高阶数组用作较低阶数组,而无需重新整形并可能创建数组的临时副本。

例如,这可用于使用矩阵向量运算收缩三阶数组的第三维

subroutine matmul312(amat, bvec, cmat)
  real(dp), contiguous, intent(in), target :: amat(:, :, :)
  real(dp), intent(in) :: bvec(:)
  real(dp), contiguous, intent(out), target :: cmat(:, :)
  real(dp), pointer :: aptr(:, :)
  real(dp), pointer :: cptr(:)

  aptr(1:size(amat, 1)*size(amat, 2), 1:size(amat, 3)) => amat
  cptr(1:size(cmat)) => cmat

  cptr = matmul(aptr, bvec)
end subroutine matmul312