I recently read (here) that pure subroutines can allow better parallization optimizations. Assuming this is true, is there a way that I can make the following routine pure?
subroutine diff_stag(operator,dfdh,f,T,dir,pad,gt)
implicit none
procedure(stencils_stag) :: operator
type(realField),intent(inout) :: dfdh
type(realField),intent(in) :: f
type(triDiag),intent(in) :: T
integer,intent(in) :: dir,pad,gt
integer :: i,j,k
select case (dir)
case (1)
!$OMP PARALLEL DO SHARED(T,f,gt)
do k=1+pad,f%s(3)-pad; do j=1+pad,f%s(2)-pad
call operator(dfdh%f(:,j,k),f%f(:,j,k),T,f%s(dir),dfdh%s(dir),gt)
enddo; enddo
!$OMP END PARALLEL DO
case (2)
!$OMP PARALLEL DO SHARED(T,f,gt)
do k=1+pad,f%s(3)-pad; do i=1+pad,f%s(1)-pad
call operator(dfdh%f(i,:,k),f%f(i,:,k),T,f%s(dir),dfdh%s(dir),gt)
enddo; enddo
!$OMP END PARALLEL DO
case (3)
!$OMP PARALLEL DO SHARED(T,f,gt)
do j=1+pad,f%s(2)-pad; do i=1+pad,f%s(1)-pad
call operator(dfdh%f(i,j,:),f%f(i,j,:),T,f%s(dir),dfdh%s(dir),gt)
enddo; enddo
!$OMP END PARALLEL DO
case default
stop 'Error: dir must = 1,2,3 in delGen_T in ops_del.f90.'
end select
end subroutine
The problem, I believe, is that the select case introduces a side effect, which is unallowed.
Is there a way that I can slice the fields f%f(i,j,k)
and dfdh%f(i,j,k)
so that the select case is not needed?
Any help is greatly appreciated!