Цикл Fortran DO по общему набору индексов?

У меня есть N-мерный набор данных (например, реальные числа), который хранится в виде одномерного массива с дополнительным массивом измерений, который указывает исходные размеры.

Кроме того, даны функции для вывода одномерного индекса из N-D индекса и наоборот.

Я пытаюсь понять, как мне сделать цикл выполнения (или эквивалент) для общего N-мерного индекса (который, конечно, будет преобразован в 1D-индекс) из некоторого набора ограничивающих нижних индексов в набор верхних индексов . Поэтому мне нужен «N-мерный» цикл, который не проходит по всем значениям — только часть массива, поэтому выполнение линейного индекса эквивалентного одномерного массива не имеет значения (по крайней мере без модификаций).

Это схема моей проблемы:

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

Я подумал, что можно зациклиться на массиве Dims и использовать внутренний цикл, но мне не удается правильно записать процедуру.

Другой вариант, который у меня не сработал (может быть, потому, что я его неправильно использую?), — это FORALL, так как он требует указания каждого индекса отдельно.


fortran do-loops
person TheWhitestOfFangs    schedule 19.06.2018    source источник

Ответы (2)


arrow_upward
2
arrow_downward

Если бы вы знали размеры ваших массивов во время компиляции, вы могли бы выполнить ряд вложенных циклов DO, каждый из которых выполняется между парами компонентов MinIndex и MaxIndex. Поскольку вы не знаете размеров, это невозможно.

Самая простая стратегия, которую я могу придумать, - это цикл с одним циклом DO по всем одномерным индексам. Для каждого из них вычислите N-мерный индекс и проверьте, находится ли он в границах, заданных MinIndex и MaxIndex: если это так, продолжайте делать то, что вам нужно; если это не так, отбросьте этот 1D-индекс и перейдите к следующему. Если индексы последовательные, вы можете сделать что-то умнее, пропуская блоки индексов, которые, как вы знаете, вас не интересуют.

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

Обратите внимание, что в том, что касается манипулирования индексами, эта стратегия совместима с распараллеливанием цикла.

Также обратите внимание, что переменные Фортрана должны начинаться с буквы: 1Dindex не является допустимым именем переменной Фортрана.

person jme52    schedule 19.06.2018
comment
Спасибо. Это очень похоже на то, о чем я думал в конце концов. - person TheWhitestOfFangs; 20.06.2018

arrow_upward
2
arrow_downward

Так что это фактическая процедура, с которой я закончил. Я проверил это и надеюсь, что это будет полезно для вас, как это было для меня. (Я знаю, что это плохо прокомментировано, но в вопросе есть все детали, и на данный момент у меня очень мало времени)

Также спасибо ripero за помощь! Я решил не использовать подход CYCLE, так как предполагаю, что код будет работать более эффективно, когда выполняется только фактическое количество циклов.

  !-----------------------------------------------------------
  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
person TheWhitestOfFangs    schedule 19.06.2018
comment
Спасибо, что приняли мой ответ. Обратите внимание, что в вашей подпрограмме вам нужно явно передать ранг только в том случае, если вы хотите иметь возможность работать с измерениями первого ранга, в противном случае вы можете легко определить его, например, по размеру Dims. Также имейте в виду, что хотя вы выиграете в скорости, избегая вызовов GetNDInd и работая с меньшим количеством циклов, в вашем основном цикле каждая итерация зависит от предыдущих, поэтому вы не сможете распараллелить ее. - person jme52; 20.06.2018
comment
В моем предыдущем комментарии, пожалуйста, s/меньшее количество циклов/меньшее количество итераций/, извините. - person jme52; 20.06.2018
comment
Спасибо за комментарии. Я передаю ранг только для того, чтобы быстро определить Curr, так как я использую этот вызов много раз и могу вычислить ранг только один раз для многих вызовов. Кроме того, параллельная работа с этим не была проблемой с самого начала, но определенно имеет смысл для тех, кто будет использовать это в будущем. - person TheWhitestOfFangs; 21.06.2018