Летом 1998 года в ЛВТА появился 8-процессорный компьютер SPP-2000 производства фирмы Hewlett-Packard. Однако обычные наши фортранные программы не в состоянии использовать одновременно более одного его процессора. Настоящий опус и призван показать, как можно это сделать, используя инструментальный пакет MPI, входящий в состав программного обеспечения SPP-2000.

Начнем с так называемого закона Амдаля.

Пусть 0 <= S <= 1 - доля вычислительных операций Вашей программы, которые должны совершаться сугубо последовательно. Тогда при одновременном использовании P процессоров Вы можете ускорить свою программу максимум в

В частности, если S = 0.1, то K < 10 при любом P. Если же S = 0 (чего в природе, вообще говоря, не наблюдается), то K = P. Напомним, что у нас P = 8 и в ближайшее время вряд ли увеличится. Так что решайте сами - окупятся Ваши труды, или нет! Впрочем, пакет MPI позволяет распараллеливать программу для выполнения на любом количестве разнотипных компьютеров, соединенных сетью EtherNet.

Тех, кто хочет получить более подробную информацию о распараллеливании вычислений или подробное описание пакета MPI, мы отсылаем к специализированному серверу МГУ.

Наиболее доступным примером, пожалуй, является программа умножения матриц, где все элементы матрицы-произведения можно вычислять параллельно, т.е. S = 0.

Вот "нераспараллеленный" вариант:

      program mumu         ! matrix multiplication
      parameter (N=400)    ! matrix dimension
      real*8 A(N,N),B(N,N),C(N,N)
      real*8 t
C...     мы опустили формирование исходных матриц А и В
      do i=1,N
        do j=1,N
          t=0.0
          do k=1,N
            t=t+A(i,k)*B(k,j)
          enddo
          C(i,j)=t
        enddo
      enddo
      end

Теперь попытаемся проделать эту же работу коллективом из нескольких процессов. Процессы, пронумерованные от 0 до P - 1, исполняют один и тот же программный код, используя независимо работающие процессоры. Процесс 0 распределяет работу между всеми исполнителями, пересылая им обе исходные матрицы А и В. Каждый исполнитель (в том числе и сам процесс 0) вычисляет "свои" столбцы матрицы С, после чего пересылает результат своей работы обратно процессу 0.

      program mumu         ! matrix multiplication (parallel version)
      parameter (N=400)    ! matrix dimension
      include 'mpif.h'     ! здесь описаны нужные нам MPI-обьекты
      integer status(MPI_STATUS_SIZE) ! важно сказать,что это массив!
      integer comm,typ,tag,ierr,myProcess,P
      data tag/0/, typ/MPI_DOUBLE_PRECISION/, comm/MPI_COMM_WORLD/
      real*8 A(N,N),B(N,N),C(N,N)
      real*8 t

C...   инициализация MPI : запрос номера "своего" процесса

      call MPI_Init(ierr)
      call MPI_Comm_rank(comm,myProcess,ierr)  ! кто я ?
      call MPI_Comm_size(comm,P,ierr)          ! сколько всего нас ?

      nc=N/P           ! какие столбцы матрицы С должен я посчитать ?
      nrest=mod(N,P)   ! (а надо поделить их приблизительно поровну)
      if(myProcess.lt.nrest) then
        nc1=1+(nc+1)*myProcess     ! это номер первого столбца
        nc2=nc1+nc                 ! а это номер последнего
      else
        nc1=1+(nc+1)*nrest+nc*(myProcess-nrest)
        nc2=nc1+nc-1
      endif
      write(*,*) ' Process',myProcess,'  of',P,
     -           '  started for columns from ',nc1,' till ',nc2

C...    Процесс 0 рассылает остальным обе исходные матрицы
C...          (  мы опустили их формирование )
      if (myProcess.eq.0) then
        do i=1,P-1
          call MPI_Send(A,N*N,typ,i,tag,comm,ierr)
          call MPI_Send(B,N*N,typ,i,tag,comm,ierr)
        enddo
      else
        call MPI_Recv(A,N*N,typ,0,tag,comm,status,ierr)
        call MPI_Recv(B,N*N,typ,0,tag,comm,status,ierr)
      endif

C--- Все начали работать ...

        do i=1,N
          do j=nc1,nc2 ! каждый вычисляет только свою часть матрицы С
            t=0.0
            do k=1,N
              t=t+A(i,k)*B(k,j)
            enddo
            C(i,j)=t
          enddo
        enddo
      enddo

C---  об'единяем результаты в памяти процесса 0
C     (поскольку матрицы хранятся в памяти по столбцам,
C     подряд идущие столбцы можно пересылать за одно обращение!)

      if(myProcess.eq.0) then
        do i=1,P-1         ! 0-й процесс обращается ко всем остальным
          if(i.lt.nrest) then
            nc1=1+(nc+1)*i ! надо вспомнить, кто какие столбцы считал
            k=nc+1
          else
            nc1=1+(nc+1)*nrest+nc*(i-nrest)
            k=nc
          endif
          call MPI_Recv(С(1,nc1),N*k,typ,i,tag,comm,status,ierr)
        enddo
      else                  ! а ненулевые процессы это знают
          call MPI_Send(C(1,nc1),N*(nc2-nc1+1),typ,0,tag,comm,ierr)
      endif

      write(*,*) ' Process',myProcess,'  finished.'
      call MPI_Finalize(ierr)   ! Ну, все...
      end

Как заставить работать это произведение?

Во-первых, Вы должны получить доступ к пакету MPI. Для этого надо дополнить свои переменные окружения PATH и MANPATH:

         setenv PATH /opt/mpi/bin:$PATH
         setenv MANPATH /opt/mpi/share/man:$MANPATH 

Добавьте это заклинание к своему .login -файлу.

Во-вторых, вместо транслятора f77 вызывайте mpif77:

         mpif77 example.f -o primer

В-третьих, при запуске программы указывайте, на сколько процессов Вы хотите ее распараллелить:

         primer -np 3                 - в данном случае на троих

Если вызовете без параметров - будет работать в одиночку. Кстати - ничто не мешает Вам указывать число процессов большее, чем количество имеющихся в наличии процессоров! Просто Ваши процессы будут простаивать в очереди к процессорам, напрасно расходуя ресурсы системы.

Мои эксперименты с этой программой при N=400 и разным числом процессов P показали, что при всех 1< P <6 работа выполняется в P раз быстрее, чем в однопроцессном варианте. При P > 5 уменьшения времени уже нет, наоборот - процессы начинают "толкаться" в памяти компьютера.

Действительно, в этой программе не требуется межпроцессных коммуникаций во время счета, обмен информацией требуется только в начале и в конце работы. Поэтому, пока системе хватает ресурсов памяти, ускорение и должно линейно зависеть от числа процессов (три землекопа выкопают ту же самую яму втрое быстрее, чем один). Если же все землекопы сразу в одну яму не поместятся - они будут только мешать друг другу!

Итак, мы видим, что даже в простейших случаях распараллеливание программы требует изрядных усилий. Более того, в ходе вычислений как правило необходимы межпроцессные коммуникации, которые могут вообще "съесть" весь эффект от распараллеливания. Здесь все зависит от соотношения цены программы и стоимости Вашего труда:

• Если программа легко распараллеливается - почему бы не сделать это?
• Если программа не очень нужная - стоит ли мучиться?
• Если программа просто незаменима - может быть стоит потрудиться?

Еще одно замечание. Интуиция подсказывает, что достаточно легко могут быть распараллелены так называемые Монте-Карловские программы, где вычислительной обработке подвергаются независимые события, сгенерированные с помощью датчика случайных чисел. Здесь важно обеспечить, чтобы каждый из параллельно работающих процессов получал свою, независимую от остальных процессов, серию случайных чисел. Для этого каждый из процессов, начиная свою работу, должен как-то по-своему инициализировать датчик.

На наш взгляд, идеальным датчиком для использования в распараллеленной программе является датчик, предложенный G.Marsaglia, способный выдавать до 32000 независимых серий равномерно распределенной на [0,1] случайной величины. Лучше всего инициализировать серию номером своего процессора:

        ...
      call MPI_Comm_rank(comm,myProcess,ierr)  ! кто я ?
      call RandomInitiate(myProcess,myProcess) ! начинаем свою серию
        ...

Датчик входит в состав нашей библиотеки JINRLIB. Немаловажным обстоятельством является то, что это самый быстрый из известных нам датчиков: всего 5 сложений и ни одного умножения с плавающей запятой!

И последнее: использование пакета MPI вовсе не уменьшит мобильность Вашей программы. На машинах, где нет MPI, Вы можете использовать заглушку:

  файл mpif.h :
         parameter(MPI_COMM_WORLD=0)
         parameter(MPI_DOUBLE_PRECISION=0)
         parameter(MPI_STATUS_SIZE=10)
  файл mpi.for :
         subroutine MPI_Comm_rank(comm,myProcess,ierr)
         myProcess=0   ! наш процесс имеет номер 0
         return
         subroutine MPI_Comm_size(comm,nProcs,ierr)
         nProcs=1      !    и он единственный
         return
           ...         ! все остальные "MPI-программы" - пустые !!!

Применение этой заглушки позволит запускать Вашу программу в однопроцессном режиме, не меняя ее текста.

Здесь можно посмотреть полный текст программ на языках Фортран и С.