PGI, OpenMP and namelist in Fortran

Question

I encounter a problem while trying to read a namelist in a Fortran program, using OpenMP and the Portland Group compiler.

What I am trying to do is simple: I call a read_namelist subroutine in a SINGLE region, where I initialize the parameters I want to read from the namelist, and then I open, read, close the namelist. The parameters I'm reading in the namelist are threadprivate, and I spread them to the other threads after the reading.

While it works perfectly well with GNU and Intel compilers, it fails with PGI and I cannot get why. I get no error, but the read parameters are equal to the default parameters, not the ones I read from the namelist.

Here is an example of what I am trying to do:

program read_input
  !$ use OMP_LIB
  use params
  implicit none

  integer :: rank=0, nthreads=1

  !$OMP PARALLEL DEFAULT(PRIVATE)
  !$ rank = OMP_GET_THREAD_NUM()
  !$ nthreads = OMP_GET_NUM_THREADS()   
  !$OMP SINGLE 
  print*, 'There is ', nthreads, ' threads running'

  call read_nml
  !$OMP END SINGLE COPYPRIVATE(nx, ny, nz)

  print*, 'Rank: ', rank
  print*, 'nx, ny, nz: ', nx, ny, nz
  !$OMP END PARALLEL   

contains
  subroutine read_nml
    use params
    implicit none
    namelist /input_params/ nx, ny, nz

    call default_parameters
    print*, 'nx, ny, nz (default): ', nx, ny, nz

    open(unit=1, file='input', status='old')
    read(1, input_params)
    close(1)
    print*, 'nx, ny, nz (read): ', nx, ny, nz

    return
  end subroutine read_nml

  subroutine default_parameters
    use params
    implicit none

    nx = 2; ny = 2; nz = 2

    return
  end subroutine default_parameters
end program read_input

the module params is very simple and contains only:

module params
  integer :: nx, ny, nz
  !$OMP THREADPRIVATE(nx, ny, nz)
end module params

Compiling with pgfortran, here is the output I get (with 2 threads):

 Start program: read_input
 There is             2  threads running
 nx, ny, nz (default):             2            2            2
 Rank:             0
 nx, ny, nz:             2            2            2
 Rank:             1
 nx, ny, nz:             2            2            2

And if I compile the same piece of code with Intel or GNU compilers (still with 2 threads):

 Start program: read_input
 There is            2  threads running
 nx, ny, nz (default):            2           2           2
 nx, ny, nz (read):           10          10          10
 Rank:            0
 nx, ny, nz:           10          10          10
 Rank:            1
 nx, ny, nz:           10          10          10

Any thought or hint will be appreciated!

Answer 1

I don't get the exact reason, but at least I found a workaround, playing with the code.

If the parameters read in the namelist are private in the subroutine, they can be read with no problem; therefore replacing

call read_nml

by

call read_nml(nx, ny, nz)

and the subroutine read_nml by

  subroutine read_nml(nx, ny, nz)
    implicit none
    integer :: nx, ny, nz             
    namelist /input_params/ nx, ny, nz

    call default_parameters
    print*, 'nx, ny, nz (default): ', nx, ny, nz

    open(unit=1, file='input', status='old')
    read(1, input_params)
    close(1)
    print*, 'nx, ny, nz (read): ', nx, ny, nz

    return
  end subroutine read_nml

works well. I guess it is once again a problem of status (private) of the attributes, but I don't get why Intel and GNU compilers handled it with no problem whereas PGI compiler is not able to handle it. Actually, that is why the parameters were threadprivate in my module, to prevent this kind of behaviour. If somebody can give me a better answer than the one I have, I'm still interested!