如何构建包含9个小矩阵的大矩阵

Question

如何构建包含9个小矩阵的大矩阵

5

我有九个矩阵，它们的维度为（N × N）

A1（i，j），A2（i，j），A3（i，j），A4（i，j），A5（i，j），A6（i，j），A7（i，j），A8（i，j），A9（i，j）

然后我想构建一个更大的矩阵（3N × 3N），其中包括这九个矩阵：

A = [A1 A2 A3
     A4 A5 A6
     A7 A8 A9]

在Fortran中，我可以使用命令行作为

do i=1,FN
   do j=1,FML
      A(i,j) = [A1(i,j),A2(i,j),A3(i,j);A4(i,j),A5(i,j),A6(i,j);A7(i,j),A8(i,j),A9(i,j)]
   end do
end do

- Jeremy_Tamu

你想用什么语言来构建？最后的命令行并没有太多意义。使用那个命令行是什么意思？ - Vladimir F Героям слава

如果您需要，您可以始终将子数组用作 A(1:3,1:3) = A1。 - Vladimir F Героям слава

1

@VladimirF 谢谢，先生。 - Jeremy_Tamu

3个回答

3

尽管Fortran在数组操作方面很有帮助，但创建块状矩阵并不像您从其他语言中期望的那样优雅。

使用数组构造器来创建所需的矩阵是可能的，就像使用标量元素（链接1）一样。也就是说，RESHAPE([A1, A2, A3, ..., A9],[3*N,3*N])将给您一个3*Nx3*N矩阵。只是这不是您想要的那个。

与其他问题/答案一样，数组构造器[...]考虑数组元素顺序来创建长度为9*N**2的秩-1数组，然后将其reshape为一个方阵。在其他示例中使用标量元素，在此处您需要用于构造器中数据的数组。构造器的元素本身按数组元素顺序取出，这相当于

[A1(1,1), A1(2,1), ..., A1(1,2), A1(2,2), ..., A2(1,1), ... ]

这是不需要的。

因此，构造函数可能是类似于这样的东西

[A1(:,1), A4(:,1), A7(:,1), A1(:,2), ..., A6(:,3), A9(:,3)]

这个方法虽然可行，但是不太方便。

也许有其他技巧可以在构造函数中实现更加“优雅”的方式，但正如 Vladimir F 的评论所说，直接将各个块进行直接赋值可能会更加简洁明了:

A(1:N,1:N) = A1
A(1:N,N+1:2*N) = A2
A(1:N,2*N+1:3*N) = A3
A(N+1:2*N,1:N) = A4
....

- francescalus

如果您可以将第一个集合设置为9xnxn矩阵，则至少可以在循环中执行该赋值。 - agentp

0

program reshape_test
   implicit none
   integer, parameter :: N = 2
   integer, dimension(N,N) :: A1,A2,A3,A4,A5,A6,A7,A8,A9
   character(20) fmt
   integer A(3*N,3*N)
   A1 = reshape([11,21,12,22],[N,N])
   A2 = reshape([13,23,14,24],[N,N])
   A3 = reshape([15,25,16,26],[N,N])
   A4 = reshape([31,41,32,42],[N,N])
   A5 = reshape([33,43,34,44],[N,N])
   A6 = reshape([35,45,36,46],[N,N])
   A7 = reshape([51,61,52,62],[N,N])
   A8 = reshape([53,63,54,64],[N,N])
   A9 = reshape([55,65,56,66],[N,N])
   write(fmt,'(*(g0))') '(a/',N,'(i2:1x))'
   write(*,fmt) 'A1 = ',transpose(A1)
   write(*,fmt) 'A2 = ',transpose(A2)
   write(*,fmt) 'A3 = ',transpose(A3)
   write(*,fmt) 'A4 = ',transpose(A4)
   write(*,fmt) 'A5 = ',transpose(A5)
   write(*,fmt) 'A6 = ',transpose(A6)
   write(*,fmt) 'A7 = ',transpose(A7)
   write(*,fmt) 'A8 = ',transpose(A8)
   write(*,fmt) 'A9 = ',transpose(A9)
   A = reshape([reshape([A1,A4,A7,A2,A5,A8,A3,A6,A9],[N,3,N,3],order=[1,3,2,4])],[3*N,3*N])
   write(fmt,'(*(g0))') '(a/',3*N,'(i2:1x))'
   write(*,fmt) 'A = ',transpose(A)
end program reshape_test

超出 4774 个字符的限制

! [Sigh]
program reshape_test
!
! PROBLEM: Insert J*K M X N matrices AI(I)%A into a J X K block matrix A.
! The blocks AI(I)%A may be input in column-major or row-major order, and
! may be transposed or direct.
!
! SOLUTION: The RESHAPE intrinsic may perform the required pivot
! operation, if the command is designed carefully.
! Step 1: Work out the strides the inputs in array element order will
! be seen in output matrix.
! Step 2: The first 3 extents of the intermediate array may then be
! obtained by sorting the strides ascending and taking the ratios
! between adjacent strides. The 4th extent is such as to match the
! size of the intermediate array with the inputs.
! Step 3: The order of the indices is set so that each index steps
! through the inputs at its proper stride.
!
! EXAMPLE: In the column-major, direct case the second element of the
! inputs will be AI(1)%A(2,1) which should map to A(2,1) in the output,
! at offset of 1 from A(1,1), thus the first stride is 1.
! The second level block will be AI(1)%A(:,2) which should map to
! A(1:M,2) in the output, at offset of M*J from A(1:M,1), thus the
! second stride is M*J.
! The third level block will be AI(K+1)%A which should map to
! A(M+1:2*M,1:N), at offset M from A(1:M,1:N), so the third stride
! is M.
! The fourth level block is AI(2:K*(J-1)+2:K)%A which should map to
! A(1:M*J,N+1:2*N) at offset of M*N*J from A(1:M*J,1:N), so the fourth
! stride is M*N*J.
! Now strides = [1,M*J,M,M*N*J]
! Sort ascending: [1,M,M*J,M*N*J]
! Take ratios: [M,J,N]
! (M*N)*(J*K)/(M*J*N) = K, so SHAPE = [M,J,N,K].
! The first stride is 1, and that is stride take by the first index.
! The second stride is M*J, taken by the third index.
! The third stride is M, taken by the second index.
! The fourth stride is M*N*J, taken by the fourth index of the
! intermediate array. Hence ORDER = [1,3,2,4]
! Finally, to specify the input blocks in column major order such
! that they look read as row ajor in the output as specified,
! SOURCE = [((AI(IJ*(K-1)+IK)%A,IJ=1,J),IK=1,K)]
! Now that the elements of A are in the right order, a second RESHAPE
! casts the collection into the desired shape.
!
! EXERCISE: Permuting the elements of an array in bit-reversed order is
! another class of pivot operation. Write a RESHAPE invocation that
! performs this pivot operation. What are its limitations?
!
   implicit none
   integer M,N ! Each A_I is an M X N matrix
   integer J,K ! A is a J X K block matrix
   type A_type
      character(:), allocatable :: A(:,:) ! A block!
   end type A_type
   type(A_type), allocatable :: AI(:) ! Input blocks
   character(:), allocatable :: A(:,:) ! Output block matrix
   character(20) test_string ! To find max length of index
   integer LI ! Max length of block index
   integer LM ! Max length of row index
   integer LN ! Max length of column index
   integer I ! Block index
   integer IM ! Row of A(I)%A
   integer IN ! Column of A(I)%A
   character(40) fmt ! Variable FORMAT
   integer IJ ! Row of A
   integer IK ! Column of A

! Define problem dimensions
   M = 2
   N = 2
   J = 3
   K = 3
! Get max index lengths
   write(test_string,'(i0)') J*K
   LI = len_trim(test_string)
   write(test_string,'(i0)') M
   LM = len_trim(test_string)
   write(test_string,'(i0)') N
   LN = len_trim(test_string)
! Create FORMAT for array element label
   write(fmt,'(4(a,i0))') "('A',i0.",LI,",'(',i",LM,",',',i",LN,",')')"
! Create blocks
   allocate(AI(J*K))
   do I = 1, J*K
      allocate(character(4+LI+LM+LN)::AI(I)%A(M,N))
      do IM = 1,M
         do IN = 1,N
            write(AI(I)%A(IM,IN),fmt) I,IM,IN
         end do
      end do
   end do
! Solve the 4 versions
   write(*,'(a)') 'Column-major, direct'
   write(*,'(a)') 'Input order'
   allocate(character(4+LI+LM+LN)::A(J*M,K*N))
   write(fmt,'(3(a,i0))') "(",J*K,"('A',i0.",LI,":','))"
   write(*,fmt) ((K*(IJ-1)+IK,IJ=1,J),IK=1,K)
   A = reshape([reshape([((AI(K*(IJ-1)+IK)%A,IJ=1,J),IK=1,K)],[M,J,N,K],order=[1,3,2,4])],[M*J,N*K])
   write(fmt,'(3(a,i0))') '(',N*K,'(a',len(A),':1x))'
   write(*,fmt) transpose(A)
   deallocate(A)
   write(*,'()')
   write(*,'(a)') 'Row-major, direct'
   write(*,'(a)') 'Input order'
   allocate(character(4+LI+LM+LN)::A(J*M,K*N))
   write(fmt,'(3(a,i0))') "(",J*K,"('A',i0.",LI,":','))"
   write(*,fmt) (I,I=1,J*K)
   A = reshape([reshape([(AI(I)%A,I=1,J*K)],[M,J,N,K],order=[1,3,4,2])],[M*J,N*K])
   write(fmt,'(3(a,i0))') '(',N*K,'(a',len(A),':1x))'
   write(*,fmt) transpose(A)
   deallocate(A)
   write(*,'()')
   write(*,'(a)') 'Column-major, transposed'
   write(*,'(a)') 'Input order'
   allocate(character(4+LI+LM+LN)::A(J*N,K*M))
   write(fmt,'(3(a,i0))') "(",J*K,"('A',i0.",LI,":','))"
   write(*,fmt) ((K*(IJ-1)+IK,IJ=1,J),IK=1,K)
   A = reshape([reshape([((AI(K*(IJ-1)+IK)%A,IJ=1,J),IK=1,K)],[N,J,M,K],order=[3,1,2,4])],[N*J,M*K])
   write(fmt,'(3(a,i0))') '(',M*K,'(a',len(A),':1x))'
   write(*,fmt) transpose(A)
   deallocate(A)
   write(*,'()')
   write(*,'(a)') 'Row-major, transposed'
   write(*,'(a)') 'Input order'
   allocate(character(4+LI+LM+LN)::A(J*M,K*N))
   write(fmt,'(3(a,i0))') "(",J*K,"('A',i0.",LI,":','))"
   write(*,fmt) (I,I=1,J*K)
   A = reshape([reshape([(AI(I)%A,I=1,J*K)],[N,J,M,K],order=[3,1,4,2])],[N*J,M*K])
   write(fmt,'(3(a,i0))') '(',M*K,'(a',len(A),':1x))'
   write(*,fmt) transpose(A)
   deallocate(A)
end program reshape_test

- user5713492

2

请为您的解决方案提供一些描述 - shivam

不要感到惊讶，如果你得到了负分。http://meta.stackoverflow.com/questions/276774/how-to-deal-with-post-containing-only-code - Vladimir F Героям слава

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- roygvib · Accepted Answer

仅供娱乐，您也可以使用do循环来制作大型A矩阵，如下：

do i = 1, N
    A( i,       : ) = [ A1( i,: ), A2( i,: ), A3( i,: ) ]
    A( i + N,   : ) = [ A4( i,: ), A5( i,: ), A6( i,: ) ]
    A( i + N*2, : ) = [ A7( i,: ), A8( i,: ), A9( i,: ) ]
enddo

这段代码按行主序填充A矩阵，因此小矩阵也以相同方式出现。如果真的非常必要，这也可以写成一行代码。

A = transpose( reshape(  &
        [ ( [ A1( i,: ), A2( i,: ), A3( i,: ) ], i=1,N ), &
          ( [ A4( i,: ), A5( i,: ), A6( i,: ) ], i=1,N ), &
          ( [ A7( i,: ), A8( i,: ), A9( i,: ) ], i=1,N ) ], [N*3, N*3] ))

这实际上是@francescalus答案中第二个数组构造函数的转置（以一行代码的形式）。

A = reshape(  &
        [ ( [ A1( :,i ), A4( :,i ), A7( :,i ) ], i=1,N ), &
          ( [ A2( :,i ), A5( :,i ), A8( :,i ) ], i=1,N ), &
          ( [ A3( :,i ), A6( :,i ), A9( :,i ) ], i=1,N ) ], [N*3, N*3] )

为了更进一步，我们可以像其他语言一样定义hcat和vcat程序（需要注意的是这里需要显式接口）：

function hcat( A, B, C ) result( X )
    integer, dimension(:,:) :: A, B, C
    integer :: X( size(A,1), size(A,2)+size(B,2)+size(C,2) )

    X = reshape( [ A, B, C ], shape( X ) )
endfunction

function vcat( A, B, C ) result( X )
    integer, dimension(:,:) :: A, B, C
    integer :: X( size(A,1)+size(B,1)+size(C,1), size(A,2) )

    X = transpose( reshape( &
            [ transpose(A), transpose(B), transpose(C) ], &
            [ size(X,2), size(X,1) ] ) )
endfunction

然后我们可以编写代码。

A = vcat( hcat( A1, A2, A3 ), hcat( A4, A5, A6 ), hcat( A7, A8, A9 ) )

这更接近于问题中所期望的形式：

A = [ A1 A2 A3 ; A4 A5 A6 ; A7 A8 A9 ]