创建所有可能组合的算法

Question

我正在编写备用网格代码，需要将 N 个一维网格点（以向量形式编写）组合成所有可能点的数组。例如，可以将两个向量 (a,b) 与 (c,d,e) 混合，给出以下几点：

(a,c) (a,d) (a,e) (b,c) (b,d) (b,e)

Matlab 有一个叫做 combvec 的函数：

http://www.mathworks.co.uk/help/nnet/ref/combvec.html

我正在用 FORTRAN 编写这段代码，但是我找不到底层算法。代码需要接收 N (N>1) 个向量（即 2,3...N），每个向量可以是不同的长度。有人知道算法吗？

Answer 1

我不知道 Fortran，但既然你说你找不到底层算法，我假设你一旦知道算法就可以自己编写这个。其实很容易。伪代码是这样的（假设没有重复）：

index = 0   ! or 1
for each element in first vector
    for each element in second vector
        matrix(index,1) = current element of first vector
        matrix(index,2) = current element of second vector
        index = index + 1
    end for
end for

这应该会给你一个类似于使用combvec 得到的矩阵。 可能有更有效的方法可以做到这一点，但由于我不知道 Fortran 的详细信息，很遗憾我无法为您提供帮助。在 Matlab 中，您当然会对其进行矢量化处理。

祝你好运=）

Answer 2

我认为下面的函数应该做你想做的事。它尽可能简单，将包含数据集的基数的 rank 1 数组作为输入，并返回 rank 2 数组，每个集合一列，包含该集合的索引。表达式1 + mod((i-1)/rep, N) 表示整数序列1,2,...,N 的第i 个元素，每个元素重复rep 次。

! (requires explicit interface)
pure function cartprod(v) result(cp)
  integer, intent(in) :: v(1:)
  integer :: cp(product(v), size(v,1))

  integer :: i, j, p, rep

  p = product(v)

  do j = 1, size(v,1) 
    rep = p / product(v(1:j))
    do i = 1, p
      cp(i,j) = 1 + mod((i-1)/rep, v(j))
    enddo
  enddo
end function

假设你定义了一个动态长度向量如下，你可以直接得到组合矩阵：

module dynamic_vector
  implicit none

  type :: d_vector
    integer, allocatable :: val(:)
  end type

contains

  pure function combvec(v) result(cv)
    type(d_vector), intent(in) :: v(1:)
    integer, allocatable :: cv(:,:)

    integer :: i, j, prod, rep, len, sizes(size(v,1))

    len = size(v,1)

    ! Determine sizes of the vectors, loop is necessary because we may not
    ! reference v%val if v is an array and val is allocatable.
    do i = 1, len
      sizes(i) = size(v(i)%val,1)
    enddo

    prod = product(sizes)

    ! Allocate and fill the output matrix
    allocate(cv(prod, len))

    do j = 1, len
      rep = prod / product(sizes(1:j))
      do i = 1, prod
        cv(i,j) = v(j)%val(1 + mod((i-1)/rep, sizes(j)))
      enddo
    enddo
  end function
end module

一个简短的测试程序：

program test
  use dynamic_vector
  implicit none

  type(d_vector) :: foo(2)
  integer :: i, bar(:,:)
  allocatable :: bar

  allocate(foo(1)%val, source = [1,2])
  allocate(foo(2)%val, source = [3,4,5])

  bar = combvec(foo)
  write(*,'(2(I0,X))') (bar(i,:), i = 1, 6)
end program

结果：

1 3
1 4
1 5
2 3
2 4
2 5