# OpenMP inconsistently segfaulting

I parallelized my code with openMP, and now have this bug in my code that's really odd. If it's the first time I run on a computer, it segfaults, but if I run the program again it runs fine. I've never seen anything like this and was curious if anyone has seen this before or has good resources or hints on how to fix this. It's a very simple loop that it segfaults on, and I parallelized every loop with default none. The variable declaration followed by the loop is below. I had previously posted a stripped down version of the loop for readability, but it was requested that I fill it in for the sake of more clarity :

 real(dp)    :: q, oneOverq, L, D, Pressure, u, v, En, rho, oneOverrho
real(dp)    :: normal(Mesh%numDims), wL(Mesh%numFields), cd_dql(Mesh%numFields, Mesh%numEdges)
real(dp)    :: norm_dx1(Mesh%numDims, Mesh%numDims), norm_dx2(Mesh%numDims, Mesh%numDims)
real(dp)    :: L_dx(Mesh%numDims, Mesh%numNodes), D_dx(Mesh%numDims, Mesh%numNodes)
real(dp)    :: p_dql(Mesh%numFields), L_dql(Mesh%numFields, Mesh%numEdges), Res(Mesh%numFields, Mesh%numTri)
INTEGER(i4) :: i, j, k, m, nodecount, cellNum, edgeNum, Node1, Node2, node, cell, cells(3)

!$$OMP PARALLEL DO Reduction(+:L, D, cL_dql, L_dx, cD_dql, D_Dx) default(none) shared(Mesh, W, grad, gamma, dRdX, Jac_2nd) & !$$OMP private(cellnum, edgenum, node1, node2, rho, oneoverrho, u, v, En, pressure, normal, norm_dx1, norm_dx2, p_dql, wL, cell, node, cells)
!DIR$IVDEP:LOOP do cellnum = 1, Mesh%numTri edgeNum = Mesh%airfoilEdges(i) Node1 = Mesh%edgelist(edgenum,1) Node2 = Mesh%edgelist(edgenum,2) cellNum = Mesh%edgeList(edgenum,3) call cell_stencil(Mesh, cellnum, cells) wL = W(:, cellnum) call linear_limited_reconstruction(cellnum, edgenum, Mesh, grad(:,:,cellNum), wL) rho = wL(1) oneOverrho = 1.0_dp/rho u = wL(2)*oneOverrho v = wL(3)*oneOverrho En = wL(4)*oneOverrho call get_pressure(rho, En, u, v, pressure) p_dql(1) = (gamma - 1._dp)*half*( u**2 + v**2 ) p_dql(2) = (gamma - 1._dp)*-u p_dql(3) = (gamma - 1._dp)*-v p_dql(4) = gamma-1._dp call norm_dx(node1, node2, Mesh, normal, norm_dx1, norm_dx2) L = L + Pressure*normal(2) D = D + Pressure*normal(1) do j = 1, 4 if (j == 1) then cell = cellnum else cell = cells(j-1) end if if (cell /= 0) then cL_dql(:,cell) = cL_dql(:,cell) + matmul(P_dql(:),Jac_2nd%wL_rc_dq(:,:,j,edgenum))*normal(2) cD_dql(:,cell) = cD_dql(:,cell) + matmul(P_dql(:),Jac_2nd%wL_rc_dq(:,:,j,edgenum))*normal(1) end if end do L_dx(:, node1) = L_dx(:, node1) + Pressure*norm_dx1(2, :) L_dx(:, node2) = L_dx(:, node2) + Pressure*norm_dx2(2, :) D_dx(:, node1) = D_dx(:, node1) + Pressure*norm_dx1(1, :) D_dx(:, node2) = D_dx(:, node2) + Pressure*norm_dx2(1, :) do j = 1, 6 node = Mesh%NodeStencil(j, cellnum) if (node /= 0) then do k = 1, Mesh%numFields L_dx(:,node) = L_dx(:,node) + P_dql(k)*dRdX%wL_rc_dx(k,:,j,edgenum)*normal(2) D_dx(:,node) = D_dx(:,node) + P_dql(k)*dRdX%wL_rc_dx(k,:,j,edgenum)*normal(1) end do end if end do end do !$OMP END PARALLEL DO


I can put in the full code up in the code block, but for the sake of readability I thought it was best not to. If anyone has any advice I'd greatly appreciate it.

Thanks!

• Often these sorts of heisenbugs happens because of an out-of-bounds memory access. Depending on the operating system's whim of where arrays get allocated into heap memory, that memory access can either segfault or run over into already allocated data, the first of which is obvious and the second of which is not. You can try compiling with -fsanitize=address and see what happens; there are also thread, leak, memory sanitizers, etc. Jul 24 '20 at 15:20
• As with all things, I posted and then found the answer (I think). It looks like default none didn't properly work. I hadn't specified two indices, i and j, that I was using to assign the data to the appropriate arrays, so I guess sometimes it corrected had them as private, and other times not so. It looks like your hunch was right, and I'll definintely use that flag in the future for my memory debugging. Do you have any other suggestions as to what sanitizers to use, etc? thanks!
– EMP
Jul 24 '20 at 15:30
• The address and leak sanitizers have been the most useful to me, but I code more in C/C++ than Fortran and there you have to be more careful about memory leaks. What I would definitely recommend is that you put your project is on github or gitlab, set up an automated test suite for it, and enable address sanitizer for the testing build. Once that's set up, it's easy to add more sanitizers later. Jul 24 '20 at 17:16
• If default (none) is not working properly that is a VERY broken compiler - what is it? To be honest that would be so broken that I can't believe what you attest to be the problem. But to say more we really need to see at the very least "Computational things" and the declaration of the variables. Jul 24 '20 at 18:32
• It's the intel ifort compiler. I used the i and j variables inside the computational section. From my understanding of default(none) it should have given an error since I hadn't specified either public or private. Since I made this change, I haven't had a problem.
– EMP
Jul 24 '20 at 19:07