如何避免包含`a（i）= b（i，c（i））`的显式循环？

Posted 2021-04-16

您可以使用符号：索引或根本没有索引（在F95中）将您的普通矢量加法编码为do循环。我想这对编译器没有任何影响。不幸的是，我有一个嵌套的声明，如：

do i=1,n
  a(i)=b(i,c(i))
end do

当然a=b(c)是非法的。任何摆脱显式do循环的想法？

答案

更新：经过一些思考，你实际上可以只使用一个隐含的数组构造函数（这是不明确的）。这只是

a = (/ (b(i,c(i)), i=1,SIZE(c)) /)

要么

a = [ (b(i,c(i)), i=1,SIZE(c)) ]

参见标准的section 4.8。

下面使用重塑和指针的原始答案。

它可以在没有循环的情况下完成，但是如果这是一个好的做法，它是非常值得怀疑的。

首先，您需要知道数组存储在column-major order中。

这意味着在(/ n,m /)大小的数组上，元素(i,j)位于记忆位置i + n*(j-1)。所以你需要做的就是重塑和提取，你可以通过以下方式做到：

  d = RESHAPE(b,(/ SIZE(b) /))                    ! reshape
  a = d((/ (i,i=1,SIZE(c)) /) + (c(:)-1)*SIZE(c)) ! extract

它不需要循环，但你必须分配更多的内存。除非你破解你的方式进入记忆。这可以用指针完成

  TYPEOFB, POINTER, DIMENSION(:) :: p
  p(1:SIZE(b)) => b(:,1)             ! just close your eyes
  a = p((/ (i,i=1,SIZE(c)) /) + (c(:)-1)*SIZE(c))

示例程序：

PROGRAM reshape
  IMPLICIT NONE
  !-- declaration of variables -----------------------------
  INTEGER, DIMENSION(4,3),TARGET :: b  !< target for pointer
  INTEGER, DIMENSION(4) :: c,a,a2
  INTEGER, DIMENSION(12) :: d
  INTEGER, DIMENSION(:), POINTER :: p  !< ugly pointer
  INTEGER :: i,j
  !-- setup of the test case ------------------------------
  c(1)=3; c(2)=2; c(3)=1; c(4)=2
  DO i=1,4
     DO j=1,3
        b(i,j) = i + 4*(j-1)
     END DO
  END DO
  !-- reshape and extract ---------------------------------
  d = RESHAPE(b,(/ SIZE(b) /))                    ! reshape
  a = d((/ (i,i=1,SIZE(c)) /) + (c(:)-1)*SIZE(c)) ! extract
  !-- ugly pointer hack -----------------------------------
  p(1:SIZE(b)) => b(:,1)             ! just close your eyes
  a2 = p((/ (i,i=1,SIZE(c)) /) + (c(:)-1)*SIZE(c))
  !-- print output ----------------------------------------
  DO i=1,4
     PRINT *, i, b(i,c(i)), a(i), a2(i)
  END DO
END PROGRAM reshape

输出：

 % ./a.out 
 1 9 9 9
 2 6 6 6
 3 3 3 3
 4 8 8 8

另一答案

当然有可能。您可以使用EOSHIFT或CSHIFT将适当的元素放在第一列中，然后使用RESHAPE提取第一列。但它效率低下，因为编译器可能会转移数组的所有元素而不仅仅是大小（b，1）。这是我的例子，这对我来说很难阅读。

program permute
   implicit none
   integer, parameter :: N = 5
   dimension b(N,N)
   integer i
! Initialize matrix to arbitrary value
   character :: b = reshape([(achar(i+96),i=1,size(b))],shape(b),order=[2,1])
   character(30) fmt
   integer c(size(b,1))
! Output matrices
   character a(size(b,1)),z(size(b,1))
   write(fmt,'(*(g0))') '(a/(',size(b,2),'(a)))'
! Print out input matrix
   write(*,fmt) 'b =',transpose(b)
! A sample permutation
   c = [1,3,4,5,2]
   write(fmt,'(*(g0))') '(a/(',size(b,1),'(a)))'
! Matrix a will have the output we want to emulate
   do i = 1, size(b,1)
      a(i) = b(i,c(i))
   end do
   write(*,fmt) 'a =',a
! Matrix z gets its result as specified in text. Is it correct?
   z = reshape(eoshift(b,c-1,dim=2),shape(z))
   write(*,fmt) 'z =',z
end program permute