0

I'm trying to use OpenMP in Fortran 90 to parallelize a do loop with function call inside. The code listed first runs fine. The code listed next does not. I receive a segmentation fault.

First program: $ gfortran -O3 -o output -fopenmp OMP10.f90

program OMP10  

        !$ use omp_lib  
        IMPLICIT NONE  
        integer, parameter :: n = 100000  
        integer :: i  
        real(kind = 8) :: sum,h,x(0:n),f(0:n),ZBQLU01  
        !$ call OMP_set_num_threads(4)  
        h = 2.d0/dble(n)  
        !$OMP PARALLEL DO PRIVATE(i)  
        do i = 0,n  
                x(i) = -1.d0+dble(i)*h  
                f(i) = 2.d0*x(i)  
        end do  
        !$OMP END PARALLEL DO  
        sum = 0.d0  
        !$OMP PARALLEL DO PRIVATE(i) REDUCTION(+:SUM)  
        do i = 0,n-1  
                sum = sum + h*f(i)  
        end do  
        !$OMP END PARALLEL DO  
        write(*,*) "The integral is  ", sum  
end program OMP10  

Second program: $ gfortran -O3 -o output -fopenmp randgen.f OMP10.f90

program OMP10  

    !$ use omp_lib  
    IMPLICIT NONE  
    integer, parameter :: n = 100000  
    integer :: i  
    real(kind = 8) :: sum,h,x(0:n),f(0:n),ZBQLU01  
    !$ call OMP_set_num_threads(4)  
    h = 2.d0/dble(n)  
    !$OMP PARALLEL DO PRIVATE(i)  
    do i = 0,n  
            x(i) = ZBQLU01(0.d0)  
    end do  
    !$OMP END PARALLEL DO  
    sum = 0.d0  
    !$OMP PARALLEL DO PRIVATE(i) REDUCTION(+:SUM)  
    do i = 0,n-1  
            sum = sum + h*f(i)  
    end do  
    !$OMP END PARALLEL DO  
    write(*,*) "The integral is  ", sum  
end program OMP10  

In the above command, randgen.f is a library that contains the function ZBQLU01.

D.B.
  • 137
  • 6
  • 1
    Well given the only major difference is that the second calls ZBQLU01 while the first does not it is quite important that you provide the source for that routine - random number generators are very often not thread safe so there may be an issue there. Oh, and for efficiency, avoid creating more parallel regions that you need to. In fact if you are learning OMP I would recommend avoiding !$omp parallel do altogether, rather break it down into its constituent steps, it will help you understand what is going on. – Ian Bush Mar 31 '18 at 08:14
  • Oh, and real( kind = 8) is not portable, your first code will not compile e.g. with the NAG compiler because of it. See https://stackoverflow.com/questions/838310/fortran-90-kind-parameter – Ian Bush Mar 31 '18 at 08:16
  • Welcome. Please read [ask] and also read how to format your code in questions. See my edit. Please use tag [tag:openmp] for all OpenMP questions. Tags are very important. – Vladimir F Героям слава Mar 31 '18 at 12:15
  • 1
    The 2nd program doesn't define values of f(:) – tim18 Mar 31 '18 at 12:26
  • Thanks. I will try to investigate the random number generator routine. Don't worry about f(:) by the way. This is just some unrelated code for which I modified the first do loop simply to illustrate my problem. – D.B. Mar 31 '18 at 16:05

1 Answers1

0

You cannot just call any function from a parallel region. The function must be thread safe. See What is meant by "thread-safe" code? and https://en.wikipedia.org/wiki/Thread_safety .

Your function is quite the opposite of thread safe as is quite typical for random number generators. Just notice the SAVE statements in the source code for many local variables and for a common block.

The solution is to use a good parallel random number generator. The site is not for software recommendation, but as a pointer just search the web for "parallel prng" or "parallel random number generator". I personally use a library which I already pointed to in https://stackoverflow.com/a/38263032/721644 A simple web search reveals another simple possibility in https://jblevins.org/log/openmp . And then there are many larger and more complex libraries.