Fortran + open mp

[DaPhil]

Hi,
I wrote the following little example:

program main
use omp_lib
implicit none

integer :: n, l, m, j, nDEQ, nk
character*13 :: problemType
real*8, allocatable, dimension) :: x0, xIn
real*8, allocatable, dimension, :: chiProbeArray, chiDiffArray

n = 3
m = 3

allocate(xIn, x0, chiProbeArray(m, n), chiDiffArray(m, n))

!$omp parallel do default(private), shared(chiProbeArray, chiDiffArray)
do j = 1, m
print *, omp_get_thread_num()
xIn = x0
do l = 1, n
chiProbeArray(j, l) = l * m
chiDiffArray(j, l) = l * m +1
end do
end do
!$omp end parallel do

print *, sum(chiProbeArray)
print *, sum(chiDiffArray)

deallocate(xIn, x0, chiProbeArray, chiDiffArray)

end program main

It compiles fine. But it never enters the loop, "print *, omp_get_thread_num()" is never invoked. What am I doing wrong? I noticed that it runs when I change default(private) to default(shared), but it does not give the right result (for sum(chiProbeArray) it should be 36). Anybody an idea?

Cheers
Phil

 

[gummibaer]

How do you know ?
Add a printout of the loopcounter as first line in the loop to see if it gets executed.
And add a printout of the contents of chiProbeArray to see if you have the right values there as you expectthem to be.

Norbert


The optimist believes we live in the best of all possible worlds - the pessimist fears this might be true.

 

[melmacianalf]

Shouldn't one keep the indices as the shared data as follows?

!$omp parallel do default(shared), shared(l,j)
do j = 1, m
print *, omp_get_thread_num()                                                                                                                                           
xIn = x0
do l = 1, n
chiProbeArray(j, l) = l * m
chiDiffArray(j, l) = l * m +1
end do
end do
!$omp end parallel do

 

[FJacq]

as mentioned by melmacianalf, declaring the arrays chiProbeArray and chiDiffArray private is simply wrong (a private allocatable array must always be allocated within the most external parallel loop). About the loop indexes, they are automatically private.

On the contrary, the array xln should be declared 'private' (or the instruction xln=x0 could be deleted because is never used and x0 is never initialized).

So I propose the following simple version :

!$omp parallel do
do j = 1, m 
  print *, j,omp_get_thread_num()
  do l = 1, n 
    chiProbeArray(j, l) = l * m 
    chiDiffArray(j, l) = l * m +1 
  end do 
end do 
!$omp end parallel do

François Jacq

 

[DaPhil]

Hi,
first, thanks for the reply!

@gummibaer: I did that, I didn't print anything.

@FJacq: All loop indices are private automatically? Even in nested loops? I did not declare chiDiffArray and chiProbeArray as private (or at least I think I didn't). The statement

!$omp parallel do default(private), shared(chiProbeArray, chiDiffArray)

means that everything is declared private except chiProbeArray and chiDiffArray, doesn't it?

My original code is much more complicated and this code only reflects the general structure. But in principle I do this: Within the outer loop I calculate Runge-Kutta steps (therefore another loop). I store the result in chiDiffArray and after the outer loop I calculate the mean value (I know I can do that within the parallel region via REDUCTION, but I wanted to try a simple case first!). So I thought making every variable (and I have plenty of them) private except the ones that I need for the mean value would do the job! Is that right or am I wrong?

 

[melmacianalf]

Why not start from something simple. Following is the code snippet which I use to calculate the dot product of two double precision 1D arrays ('X' and 'Y') of length 'N'

 dotproduct = 0.0d0
!$OMP PARALLEL DO &
!$OMP DEFAULT(SHARED) PRIVATE(i) &
!$OMP REDUCTION(+:dotproduct)
  do i = 1, N
    dotproduct = dotproduct + dconjg(X(i)) * Y(i)
  enddo
!$OMP END PARALLEL DO

Another simple task is to Assign a constant value 'c' to a 1D array

!$OMP PARALLEL DO &
!$OMP DEFAULT(SHARED) &
!$OMP PRIVATE(i)
  do i = 1, N
    X(i) = c
  enddo
!$OMP END PARALLEL DO

I hope these examples give you some idea

 

[FJacq]

With OpenMP, I advise simple but genaral constructs like the following one :

!$OMP PARALLEL DO
Do i=1,n
  CALL stuff-routine(i,...)
ENDDO

In stuff-routine, you just have to avoid saved local variables.

The variable are "shared" by default and this is the best choice because, with OpenMP, most variables must be "shared". A private variable is always expensive because you need to allocate it for each thread.

About the loop indexes, these are always "private" and this is true for nested loops.

François Jacq