I just found and corrected OpenMP errors in GRHydro's prim2con. I believe they have been there for some weeks.
When using OpenMP to parallelise loops, one has to take care that there are two kinds of variables: serial variables (that exist once per MPI process) and parallel variables (that exist once per thread). Since OpenMP correctness depends on everybody who modifies the code, let me briefly explain the issue below.
Variables that are declared in Fortran in the usual way are serial by default; if they are set within a loop, there is a conflict if different threads set the variable to different values at the same time. To parallelise variables, they have to be listed in a "omp private" statement at the beginning of the loop. This is inconvenient.
The wrong code looked like this:
!$OMP PARALLEL DO do k = GRHydro_stencil,cctk_lsh(3)-GRHydro_stencil+1 do j = GRHydro_stencil,cctk_lsh(2)-GRHydro_stencil+1 do i = GRHydro_stencil,cctk_lsh(1)-GRHydro_stencil+1
det = SPATIAL_DETERMINANT(gxx(i,j,k),gxy(i,j,k),gxz(i,j,k), \ gyy(i,j,k),gyz(i,j,k),gzz(i,j,k))
The variable "det" is shared between all threads, and so different threads will assign different values to it at the same time, leading to essentially random results at run time. The correct code needs to look like this:
!$OMP PARALLEL DO PRIVATE(det) do k = GRHydro_stencil,cctk_lsh(3)-GRHydro_stencil+1 do j = GRHydro_stencil,cctk_lsh(2)-GRHydro_stencil+1 do i = GRHydro_stencil,cctk_lsh(1)-GRHydro_stencil+1
det = SPATIAL_DETERMINANT(gxx(i,j,k),gxy(i,j,k),gxz(i,j,k), \ gyy(i,j,k),gyz(i,j,k),gzz(i,j,k))
Note the variable "det" is declared private.
-erik
PS: There may be other, similar errors in GRHydro. It may be time to have the code peer-reviewed after the significant changes since the release. We should also begin to use OpenMP in our test cases.
users@lists.einsteintoolkit.org