Каким будет правильное назначение переменных (личных и общих) в распараллеленном цикле выполнения данной подпрограммы GAUSSLEG?

Я новичок в openmp. Я пытаюсь распараллелить цикл do в подпрограмме GAUSSLEG. Переменные Xg, Wg и Ng берутся из матрицы модуля. Я получаю неожиданные результаты. Я запутался в правильном назначении переменных (частных и общих). Кто-нибудь может мне помочь?

SUBROUTINE GAUSSLEG(f,a,b,s)
     USE OMP_LIB
     USE MATRIC , ONLY : XG ,WG , NG
     IMPLICIT DOUBLE PRECISION(A-H,O-Z)
     external f
     xm = 0.5d0*(b+a)
     xl = 0.5d0*(b-a)
     s = 0.d0
    !$omp parallel do reduction ( + : s) default(none)
    !$omp private(j) shared(xm,xl,wg,xg,ng,dx)

     do  j=1,ng
         dx = xl*xg(j)
         s = s + wg(j)*(func(xm+dx)+func(xm-dx))
     end do
   !$omp end parallel do

    s = xl*s/2.0
    return
  END

Привет, я использовал подпрограмму gaussleg для вычисления интегрирования sin(x) от 0 до pi, я получаю тот же результат (2,5464790894), независимо от того, делаю ли я dx частным или общим, но точный результат равен 2,0. Я также пытался поставить xl*xg(j) напрямую и удалить dx, все равно получая тот же результат, что и выше. Без опции -openmp в компиляции я получаю точный результат 2.0. Это целая программа.

  MODULE MATRIC
   IMPLICIT NONE
    INTEGER , PARAMETER :: NG = 40
    DOUBLE PRECISION , PARAMETER :: PI=2.0D0*ACOS(0.0D0)
    DOUBLE PRECISION ::  XG(60) , WG(60)
  END MODULE MATRIC
  program gauss
   use matric, only : xg,wg,pi
   implicit none
   double precision :: x1,x2,a,b,ans
   external :: f
   x1 = -1.0d0 ; x2 = 1.0d0
   a  = 0.0    ; b  = PI
   call gauleg(x1,x2)
   call gaussleg(f,a,b,ans)
   write(*,*)ans
  end program gauss
  !function to be integrated
  double precision function f(x)
        implicit none
        double precision, intent(in) :: x
        f = sin(x)
  end function f
  SUBROUTINE GAUSSLEG(func,a,b,ss)
     USE OMP_LIB
     USE MATRIC , ONLY : XG ,WG , NG
     double precision,intent(in) :: a , b
     double precision,intent(out)::ss
     double precision :: xm , xl , dx
     integer :: j
     double precision,external::func
     xm = 0.5d0*(b+a)
     xl = 0.5d0*(b-a)
     ss = 0.d0
    !$OMP PARALLEL DO REDUCTION( + : ss) default(none) &
    !$OMP PRIVATE(j,dx) SHARED(xm,xl,xg,wg)
       do  j=1,ng
           dx = xl*xg(j)
           ss = ss + wg(j)*(func(xm+dx)+func(xm-dx))
       end do
    !$OMP END PARALLEL DO
      ss = xl*ss/2.0
      return
      END

shared openmp private reduction
person rabeet singh    schedule 27.05.2015    source источник
comment
Можете ли вы также привести пример ожидаемого и фактического результата?   -  person Geradlus_RU    schedule 27.05.2015
comment
Привет, сэр, извините, я сделал ошибку, ng нет в общем списке. С ng я получал ошибку, поэтому мой общий список (xm, xl, wg, xg, dx). Я попытался проверить подпрограмму, рассчитав интегрирование sin(x) с ограничением x от 0 до pi. Если я скомпилирую программу без опции -openmp, я получу 2.00000 с точным ответом   -  person rabeet singh    schedule 28.05.2015

Ответы (2)


arrow_upward
0
arrow_downward

Ваш код включает каноническую гонку данных. Вы объявили dx общим, а затем написали

         dx = xl*xg(j)

чтобы все потоки могли обновлять одну и ту же общую переменную без какой-либо координации. Я думаю, но вы обязаны проверить это, что вы можете сделать dx закрытым и каждый поток будет следить за своим собственным значением переменной.

Кстати. НЕ ИСПОЛЬЗУЙТЕ implicit набор текста, вы просто напрашиваетесь на неприятности. Напрашиваться на неприятности, когда вы пытаетесь научиться использовать OpenMP, — это просто напрашиваться на дополнительные проблемы. ИСПОЛЬЗУЙТЕ implicit none. И не отвечайте О, я просто обновляю существующую кодовую базу, в которой используется неявная типизация. Если это то, что вы делаете, делайте это правильно.

person High Performance Mark    schedule 27.05.2015

arrow_upward
0
arrow_downward

Получил точные результаты следующим образом.

ПОДПРОГРАММА QGAUSSP(func,a,b,ss)

     USE OMP_LIB
     USE MATRIC , ONLY : XG ,WG , NG
     implicit none
     double precision, intent(in) :: a , b
     double precision, intent(out):: ss
     double precision :: xm , xl , dx , xgd , wgd
     double precision :: s(NG)
     integer :: j,tid
     double precision,external::func
     xm = 0.5d0*(b+a)
     xl = 0.5d0*(b-a)
     ss = 0.d0
  !$omp parallel do private(j,xgd,wgd,dx) shared(xm,xl,xg,wg,s) num_threads(15) 
   do  j=1,ng
       xgd=xg(j)
       wgd=wg(j)
       dx = xl*xgd
       s(j)=wgd*(func(xm+dx)+func(xm-dx))
   end do
  !$omp end parallel do
   ss=sum(s) *xl/2.0
   return
  END
person rabeet singh    schedule 04.06.2015