CODEWITHSUNDEEP
fortrando-loops
TheWhitestOfFangs
TheWhitestOfFangs

Reputation: 465

Fortran DO loop over generic set of indices?

I have an N-dimensional dataset (say real numbers) which is stored as a 1D array with an additional dimension array which specifies the original dimensions.

In addition, the functions to deduce the 1-D index from an N-D index and vice-versa are given.

I'm trying to figure out how do I make a do-loop (or equivalent) for the generic N-dimensional index (which will be transformed to 1D index of course) from some set of limiting lower indices to a set of upper indices. So I need an "N-dimensional" loop which does not go over all of the values - only a portion of the array, therefore doing the linear index of an equivalent 1D array is not relevant (at least without modifications).

This is a schematic of my problem:

subroutine Test(Array,Dims,MinIndex,MaxIndex)

implicit none
real   , dimension(1:), intent(inout) :: Array
integer, dimension(1:), intent(in)    :: Dims,MinIndex,MaxIndex

integer, dimension(size(Dims)) :: CurrInd
integer :: 1Dindex

! size(Dims) can be 1, 2, 3 ,..., N
! size(MinIndex)==size(MaxIndex)==size(Dims)
! size(Array)==Product(Dims)
! 1Dindex=Get1dInd(NDindex,Dims)
! NDindex=GetNdInd(1Dindex,Dims)

! How do I actually preform this?
do CurrInd=MinIndex,MaxIndex
  1Dindex=Get1dInd(CurrInd,Dims)
  <Some operation>
enddo


end subroutine

I figured it is possible to loop over the Dims array and use an inner loop but I don't manage to write the procedure down properly.

Another option that didn't work for me (maybe because I use it incorrectly?) is FORALL, as that requires each index to be specified separately.

Upvotes: 1

Views: 1449

Answers (2)

TheWhitestOfFangs
TheWhitestOfFangs

Reputation: 465

So this is the actual procedure I ended up with. I have verified it and hope it will be useful to you as it was for me. (I know this is not well commented, but the question has all the details and my time is extremely short at the moment)

Also, thanks ripero for helping out! I decided not to use the CYCLE approach as I assume the code will work more efficiently when only the actual amount of loops are performed.

  !-----------------------------------------------------------
  subroutine Test(Array,Rank,Dims,InitInd,FinInd)
    implicit none
    real,    dimension(1:), intent(inout) :: Array
    integer,                intent(in)    :: Rank
    integer, dimension(1:), intent(in)    :: Dims
    integer, dimension(1:), intent(in)    :: InitInd,FinInd
    !-----------------------------------------------------------    
    integer :: nOuter
    integer :: i,j, OneDInd
    integer, dimension(Rank) :: Curr
    !-----------------------------------------------------------

    ! Check how many repetition for the outer-loop
    Curr=FinInd-InitInd
    nOuter=1
    do i=2,Rank
       nOuter=nOuter*(Curr(i)+1)
    enddo

    !-----------------------------------------------------------
    ! Actual looping:
    Curr=InitInd
    do j=1,nOuter
       ! Update minor indices (>1):
       do i=1,Rank
          if (Curr(i).GT.FinInd(i)) then
             ! Update next index:
             Curr(i)=InitInd(i)
             Curr(i+1)=Curr(i+1)+1
          endif
       enddo

       ! Loop over major index:
       do i=InitInd(1),FinInd(1)
          !OneDInd=Get1dInd(Curr,Dims)
          !<operation>


          ! Advance major index:
          Curr(1)=Curr(1)+1
       enddo

    enddo

  end subroutine Test

Upvotes: 2

jme52
jme52

Reputation: 1123

If you knew the dimensions of your arrays at compilation time you could do a series of nested DO loops, each of them running between pairs of components of MinIndex and MaxIndex. Since you don't know the dimensions, that's not possible.

The easiest strategy I can think of is to loop with a single DO loop over all the 1D indices. For each of them, compute the N-dimensional index, and check if it is within the bounds provided by MinIndex and MaxIndex: if it is, continue doing what you need; if it is not, discard that 1D index and go to the next one. If the indices are sequential you may be able to do something smarter that skips blocks of indices that you know you will not be interested in.

do OneDindex = 1, size(Array)
   CurrInd = GetNDInd(OneDindex, Dims)
   if ((any(CurrInd<MinIndex)) .or. (any(CurrInd>MaxIndex))) cycle
   ! <Some operation>
end do

Note that, as far as the index manipulation is concerned, this strategy is compatible with parallelising the loop.

Also note that Fortran variables must start with a letter: 1Dindex is not a valid Fortran variable name.

Upvotes: 2

Related Questions