The most typical reason for this is that not enough space is available on the stack to hold the private copy of numstrain
. Compute and compare the following two values:
- the size of the array in bytes
- the stack size limit
There are two kinds of stack size limits. The stack size of the main thread is controlled by things like process limits on Unix systems (use ulimit -s
to check and modify this limit) or is fixed at link time on Windows (recompilation or binary edit of the executable is necessary in order to change the limit). The stack size of the additional OpenMP threads is controlled by environment variables like the standard OMP_STACKSIZE
, or the implementation-specific GOMP_STACKSIZE
(GNU/GCC OpenMP) and KMP_STACKSIZE
(Intel OpenMP).
Note that most Fortran OpenMP implementations always put private arrays on the stack, no matter if you enable compiler options that allocate large arrays on the heap (tested with GNU's gfortran
and Intel's ifort
).
If you comment out the PRINT
statement, you effectively remove the reference to numstrain
and the compiler is free to optimise it out, e.g. it could simply not make a private copy of numstrain
, thus the stack limit is not exceeded.
After the additional information that you've provided one can conclude, that stack size is not the culprit. When dealing with private
ALLOCATABLE
arrays, you should know that:
- private copies of unallocated arrays remain unallocated;
- private copies of allocated arrays are allocated with the same bounds.
If you do not use numstrain
outside of the parallel region, it is fine to do what you've done in your first case, but with some modifications:
integer, allocatable :: numstrain(:)
integer :: allocate_status
integer, parameter :: n = 1000000
interface
subroutine anotherSubroutine(numstrain)
integer, allocatable :: numstrain(:)
end subroutine anotherSubroutine
end interface
!$OMP PARALLEL NUM_THREADS(4) PRIVATE(numstrain, allocate_status)
allocate (numstrain(n), stat = allocate_status)
!$OMP DO
do irep = 1, nrep
call anotherSubroutine(numstrain)
end do
!$OMP END DO
deallocate (numstrain)
!$OMP END PARALLEL
If you also use numstrain
outside of the parallel region, then the allocation and deallocation go outside:
allocate (numstrain(n), stat = allocate_status)
!$OMP PARALLEL DO NUM_THREADS(4) PRIVATE(numstrain)
do irep = 1, nrep
call anotherSubroutine(numstrain)
end do
!$OMP END PARALLEL DO
deallocate (numstrain)
You should also know that when you call a routine that takes an ALLOCATABLE
array as argument, you have to provide an explicit interface for that routine. You can either write an INTERFACE
block or you can put the called routine in a module and then USE
that module - both cases would provide the explicit interface. If you do not provide the explicit interface, the compiler would not pass the array correctly and the subroutine would fail to access its content.