Working on openACC

examples/gpmdk/src/gpmdcov_energandforces.F90

@@ -82,27 +82,29 @@ subroutine get_SKBlock_vect_local(sp,refcoord,coord,lattice_vectors&
     integer, intent(in)                  ::  norbs(:), sp(:), atnum
     real(dp), allocatable                ::  HPPP(:), HPPS(:), HSPS(:), HSSS(:), PPSMPP(:)
     real(dp), allocatable                ::  L(:), M(:), N(:)
-    real(dp), allocatable, save                ::  dr(:), rab(:,:)
-    real(dp), allocatable                :: dr_m(:)
+    real(dp), allocatable                ::  dr(:), dr_m(:), rab(:,:)
     real(dp), allocatable                ::  onsites_m(:)
     real(dp), intent(inout)              ::  blk_out(:,:)
-    real(dp), allocatable, save                ::  blk(:,:)
+    real(dp), allocatable                ::  blk(:,:)
     real(dp), intent(in)                 ::  refcoord(:),coord(:,:), lattice_vectors(:,:)
     real(dp), intent(in)                 ::  onsites(:,:)
     real(dp), intent(in)                 ::  intParams1(:,:,:),intParams2(:,:,:)
-    logical, allocatable,save                 ::  dist_mask(:), onsite_mask(:), calc_mask(:), calcs_mask(:)
-    logical, allocatable,save                 ::  calcsp_mask(:), param_mask(:,:), calc_mask_for_porbs(:)
-    logical, allocatable,save                 ::  sorb_mask(:), pxorb_mask(:), pyorb_mask(:), pzorb_mask(:)
-    logical, allocatable,save                 ::  sorb_at_mask(:), sporb_at_mask(:)
-    integer, allocatable,save                 ::  atomidx(:), atomidx_m(:), orbidx(:), orbidx_m(:), orbidx_sel(:)
-    real(dp), allocatable,save                ::  intParams(:,:)
+    logical, allocatable                 ::  dist_mask(:), onsite_mask(:), calc_mask(:), calcs_mask(:)
+    logical, allocatable                 ::  calcsp_mask(:), param_mask(:,:), calc_mask_for_porbs(:)
+    logical, allocatable                 ::  sorb_mask(:), pxorb_mask(:), pyorb_mask(:), pzorb_mask(:)
+    logical, allocatable                 ::  sorb_at_mask(:), sporb_at_mask(:)
+    integer, allocatable                 ::  atomidx(:), atomidx_m(:), orbidx(:), orbidx_m(:), orbidx_sel(:)
+    real(dp), allocatable                ::  intParams(:,:)
 
     nats = size(coord,dim=2)
     norbsall = sum(norbs)
 
+    allocate(blk(size(blk_out,dim=1),size(blk_out,dim=2)))
+
+    blk(:,:)=0.0_dp
+
     if(allocated(dr))then
        if(nats.ne.size(dr,dim=1))then
-          deallocate(blk)
           deallocate(dr)
           deallocate(atomidx)
           deallocate(orbidx)
@@ -123,7 +125,6 @@ subroutine get_SKBlock_vect_local(sp,refcoord,coord,lattice_vectors&
     endif
     if(.not.allocated(dr))then
        allocate(dr(nats))
-       allocate(blk(4,size(blk_out,dim=2)))
        allocate(atomidx(nats))
        atomidx = (/(i,i=1,nats)/)
        allocate(orbidx(norbsall))
@@ -175,9 +176,7 @@ subroutine get_SKBlock_vect_local(sp,refcoord,coord,lattice_vectors&
           endif
        enddo
     endif
-
-    blk(:,:)=0.0_dp
-
+
     do i = 1,3
        Rab(:,i) = coord(i,:)
        Rab(:,i) = modulo((Rab(:,i) - refcoord(i) + 0.5_dp*lattice_vectors(i,i)),lattice_vectors(i,i)) - 0.5_dp * lattice_vectors(i,i)
@@ -331,7 +330,7 @@ subroutine get_dH_or_dS_vect_local(dx,coords,hindex,spindex,intPairsH,onsitesH,s
     real(dp)                           ::  Rax_m(3), Rax_p(3), Ray_m(3), Ray_p(3)
     real(dp)                           ::  Raz_m(3), Raz_p(3), Rb(3), d, maxblockij
     real(dp), allocatable              ::  Rx(:), Ry(:), Rz(:), blockm(:,:,:)
-    real(dp), allocatable              ::  blockp(:,:,:), dh0(:,:), dH0x(:,:), dH0y(:,:), dH0z(:,:)
+    real(dp), allocatable              ::  blockp(:,:,:), dH0x(:,:), dH0y(:,:), dH0z(:,:)
     real(dp), intent(in)               ::  coords(:,:), dx, lattice_vectors(:,:), onsitesH(:,:)
     real(dp), intent(in)               ::  threshold
     type(bml_matrix_t), intent(inout)  ::  dH0x_bml, dH0y_bml, dH0z_bml
@@ -620,6 +619,7 @@ subroutine gpmdcov_EnergAndForces(charges)
     type(rankReduceData_t) :: mpimax_in(1), mpimax_out(1)
     integer :: k
     logical :: testsmd
+
     call gpmdcov_msMem("gpmdcov","Before gpmd_EnergAndForces",lt%verbose,myRank)
 
     if(.not.allocated(coul_forces)) allocate(coul_forces(3,sy%nats))
@@ -781,7 +781,7 @@ subroutine gpmdcov_EnergAndForces(charges)
         GFSCOUL(:,gpat%sgraph(ipt)%core_halo_index(i)+1) = syprt(ipt)%estr%FSCOUL(:,i)
         SKForce(:,gpat%sgraph(ipt)%core_halo_index(i)+1) = syprt(ipt)%estr%SKForce(:,i)
       enddo
-      
+
       call bml_deallocate(dSx_bml)
       call bml_deallocate(dSy_bml)
       call bml_deallocate(dSz_bml)

examples/gpmdk/src/gpmdcov_response.F90

@@ -86,37 +86,63 @@ subroutine gpmdcov_response_dpdmu(P1_bml,dPdMu,H1_bml,Norbs,beta,Q_bml,evals,mu0
 #endif
 #endif ! USE_OFFLOAD
 
+#ifdef USE_NVTX
+    call nvtxStartRange("Response Kernel",3)
+#endif
+
 #ifdef USE_OFFLOAD
     P1_bml_c_ptr = bml_get_data_ptr_dense(P1_bml)
     P1_bml_ld = bml_get_ld_dense(P1_bml)
 
-    !$omp target enter data map(alloc:p_0(1:HDIM))
-    !$omp target update to(p_0(1:HDIM))
-
     call c_f_pointer(P1_bml_c_ptr,P1_bml_ptr,shape=[P1_bml_ld,HDIM])
 
-#endif
-#ifdef USE_NVTX
-    call nvtxStartRange("Response Kernel",3)
-#endif
+    !call offload_kernel(p_0,P1_bml_ptr,P1_bml_ld,HDIM)
+
+#ifdef USE_OMP
+    !$omp target enter data map(alloc:p_0(1:HDIM))
+    !$omp target update to(p_0(1:HDIM))
+#else
+    !$acc enter data copyin(p_0(1:HDIM))
+#endif    
     do i = 1,m  ! Loop over m recursion steps
-#ifdef USE_OFFLOAD
+#ifdef USE_OMP
        !$omp target teams distribute default(none) &
        !$omp shared(HDIM,P1_bml_ptr,p_0)
-
+#else
+       !$acc parallel loop deviceptr(P1_bml_ptr) present(p_0)
+#endif
        do k = 1,HDIM
+#ifdef USE_OMP
           !$omp parallel do
+#else
+          !$acc loop
+#endif
           do j = 1,HDIM
              P1_bml_ptr(j,k) = 1.D0/(2.D0*p_0(j)*(p_0(j)-1.D0)+1.D0)*((p_0(j) + p_0(k))*P1_bml_ptr(j,k) &
                   & + 2.D0*(P1_bml_ptr(j,k)-(p_0(j) + p_0(k))*P1_bml_ptr(j,k))*1.D0/(2.D0*p_0(k)*(p_0(k)-1.0D0)+1.D0)*p_0(k)*p_0(k))
           enddo
+#ifdef USE_OMP
           !$omp end parallel do
-     enddo
-     !$omp end target teams distribute
-     !$omp target
-     p_0 = 1.D0/(2.D0*(p_0(:)*p_0(:)-p_0(:))+1.D0)*p_0(:)*p_0(:)
-     !$omp end target
 #else
+          !$acc end loop
+#endif
+       enddo
+#ifdef USE_OMP
+       !$omp end target teams distribute
+       !$omp target
+#else
+       !$acc end parallel
+       !$acc kernels present(p_0)
+#endif
+       p_0 = 1.D0/(2.D0*(p_0(:)*p_0(:)-p_0(:))+1.D0)*p_0(:)*p_0(:)
+#ifdef USE_OMP
+       !$omp end target
+#else
+       !$acc end kernels
+#endif
+    enddo
+#else
+    do i = 1,m  ! Loop over m recursion steps
      !$omp parallel do default(none) &
      !$omp private(k,j) &
      !$omp shared(HDIM,P1,p_0)
@@ -130,13 +156,17 @@ subroutine gpmdcov_response_dpdmu(P1_bml,dPdMu,H1_bml,Norbs,beta,Q_bml,evals,mu0
      enddo
      !$omp end parallel do
      p_0 = 1.D0/(2.D0*(p_0(:)*p_0(:)-p_0(:))+1.D0)*p_0(:)*p_0(:)
-#endif
     enddo
+#endif
 #ifdef USE_NVTX
     call nvtxEndRange
 #endif
 #ifdef USE_OFFLOAD
+#ifdef USE_OMP
     !$omp target exit data map(delete:p_0(1:HDIM))
+#else
+    !$acc exit data delete(p_0(1:HDIM))
+#endif
 #else
 
     bml_type = bml_get_type(P1_bml)