The following function should do what you want, I think. It is as simple as possible, taking as input a rank 1 array containing the cardinalities of the data sets and returning a rank 2 array, one column for each set, containing the indices for that set. The expression 1 + mod((i-1)/rep, N)
represents the i
th element of a sequence of integers 1,2,...,N
, with each element repeated rep
times.
! (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
Suppose you have defined a dynamic length vector as follows, you could directly obtain the matrix of combinations:
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
A short test program:
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
Result:
1 3
1 4
1 5
2 3
2 4
2 5