The method as provided by chw21 worked. However I quickly ran into problems for deep nested loops. The workload became too heavy to handle. Luckily in my case only indices that were distinct were of interest. Moreover the order did not matter, i.e. the indices 5,3,2 would yield the same result as 2,3,5. Basically the problem degenerated to providing all the combinations on a lottery ticket. If your problem is of such a nature the code below might be of interest.
program indices
implicit none
integer, dimension(:), allocatable :: ns
integer :: i,j,k,ni,np,nt,ntmp
integer*8 :: nc
print *, "Number of towns to visit"
read(*, *) np
allocate(ns(np))
print *, "Total number of towns"
read(*, *) nt
if (nt<=0) then
print*,' Error: Please provide a positive value'
stop
endif
if(nt<np) then
print*,' Error: Number of towns to visit must be less'
print*,' than or eqaul to total number of towns.'
stop
endif
! Initialize .....
do i=1,np
ns(i)=i
enddo
ntmp=nt-np
nc=0
! ................
print*,' Combinations of towns to visit..:'
do
print*,ns(:)
! Do the appropriate work with ns here.
! .......
!
! Provide a new combination:
nc=nc+1
if (ns(np)<nt) then
ns(np)=ns(np)+1
elseif (ns(1)==(ntmp+1)) then
exit
else
do i=2,np
if(ns(i)==(ntmp+i)) then
ni=ns(i-1)
k=0
do j=i-1,np
k=k+1
ns(j)=ni+k
enddo
exit
endif
enddo
endif
enddo
print*,' Number of combinations..........:',nc
end program indices
elemental
function. For others a recursive approach is appropriate. So, why do you want to do this ?