TABLE OF CONTENTS
ABINIT/alloc_hamilt_gpu [ Functions ]
NAME
alloc_hamilt_gpu
FUNCTION
allocate several memory pieces on a GPU device for the application of Hamiltonian using a GPU
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (MT, FDahm) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
atindx1(natom)=index table for atoms, inverse of atindx (see gstate.f) dtset <type(dataset_type)>=all input variables for this dataset gprimd(3,3)=dimensional reciprocal space primitive translations mpi_enreg=information about MPI parallelization nattyp(ntypat)= # atoms of each type. option=0: allocate data for local operator (FFT) 1: allocate data for nonlocal operator 2: allocate both psps <type(pseudopotential_type)>=variables related to pseudopotentials use_gpu_cuda= 0 or 1 to know if we use cuda
OUTPUT
(no output - only allocation on GPU)
PARENTS
gstate,respfn
CHILDREN
free_gpu_fourwf,free_nonlop_gpu
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 47 subroutine alloc_hamilt_gpu(atindx1,dtset,gprimd,mpi_enreg,nattyp,npwarr,option,psps,use_gpu_cuda) 48 49 use defs_basis 50 use defs_datatypes 51 use defs_abitypes 52 use m_profiling_abi 53 use m_xmpi 54 #if defined HAVE_GPU_CUDA 55 use m_initcuda 56 #endif 57 58 !This section has been created automatically by the script Abilint (TD). 59 !Do not modify the following lines by hand. 60 #undef ABI_FUNC 61 #define ABI_FUNC 'alloc_hamilt_gpu' 62 !End of the abilint section 63 64 implicit none 65 66 !Arguments ------------------------------------ 67 !scalars 68 integer,intent(in) :: option,use_gpu_cuda 69 type(dataset_type),intent(in) :: dtset 70 type(MPI_type),intent(in) :: mpi_enreg 71 type(pseudopotential_type),intent(in) :: psps 72 !arrays 73 integer,intent(in) :: atindx1(dtset%natom),nattyp(dtset%ntypat),npwarr(dtset%nkpt) 74 real(dp),intent(in) :: gprimd(3,3) 75 76 !Local variables------------------------------- 77 !scalars 78 #if defined HAVE_GPU_CUDA 79 integer :: dimekb1_max,dimekb2_max,dimffnl_max,ierr,ikpt,npw_max_loc,npw_max_nonloc 80 integer ::npwarr_tmp(dtset%nkpt) 81 #endif 82 83 ! ************************************************************************* 84 85 if (use_gpu_cuda==0) return 86 87 #if defined HAVE_GPU_CUDA 88 !=== Local Hamiltonian === 89 if (option==0.or.option==2) then 90 ! Compute max of total planes waves 91 npw_max_loc=0 92 if(mpi_enreg%paral_kgb==1) then 93 npwarr_tmp=npwarr 94 call xmpi_sum(npwarr_tmp,mpi_enreg%comm_bandfft,ierr) 95 npw_max_loc =maxval(npwarr_tmp) 96 else 97 npw_max_loc=dtset%mpw 98 end if 99 ! Initialize gpu date needed in fourwf 100 ! ndat=bandpp when paral_kgb=1 101 if(mpi_enreg%paral_kgb==1) then 102 call alloc_gpu_fourwf(dtset%ngfft,dtset%bandpp,npw_max_loc,npw_max_loc) 103 else 104 call alloc_gpu_fourwf(dtset%ngfft,1,npw_max_loc,npw_max_loc) 105 end if 106 end if 107 !=== Nonlocal Hamiltonian === 108 if (option==1.or.option==2) then 109 ! Compute max of total planes waves 110 npw_max_nonloc=0 111 if(mpi_enreg%paral_kgb==1) then 112 npwarr_tmp=npwarr 113 call xmpi_sum(npwarr_tmp,mpi_enreg%comm_band,ierr) 114 npw_max_nonloc =maxval(npwarr_tmp) 115 else 116 npw_max_nonloc=dtset%mpw 117 end if 118 ! Initialize all gpu data needed in nonlop 119 dimffnl_max=4 120 ! if (abs(dtset%berryopt) == 5) dimffnl_max=4 121 dimekb1_max=psps%dimekb 122 if (dtset%usepaw==0) dimekb2_max=psps%ntypat 123 if (dtset%usepaw==1) dimekb2_max=dtset%natom 124 call alloc_nonlop_gpu(npw_max_nonloc,npw_max_nonloc,dtset%nspinor,dtset%natom,& 125 & dtset%ntypat,psps%lmnmax,psps%indlmn,nattyp,atindx1,gprimd,& 126 & dimffnl_max,dimffnl_max,dimekb1_max,dimekb2_max) 127 end if 128 call xmpi_barrier(mpi_enreg%comm_cell) 129 #else 130 if (.false.) then 131 write(std_out,*) atindx1(1),dtset%natom,gprimd(1,1),mpi_enreg%me,nattyp(1),option 132 end if 133 #endif 134 135 end subroutine alloc_hamilt_gpu
ABINIT/dealloc_hamilt_gpu [ Functions ]
NAME
dealloc_hamilt_gpu
FUNCTION
deallocate several memory pieces on a GPU device used for the application of Hamiltonian using a GPU
COPYRIGHT
Copyright (C) 2000-2018 ABINIT group (MT, FDahm) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
option=0: deallocate data for local operator (FFT) 1: deallocate data for nonlocal operator 2: deallocate both use_gpu_cuda= 0 or 1 to know if we use cuda
OUTPUT
(no output - only deallocation on GPU)
PARENTS
gstate,respfn
CHILDREN
free_gpu_fourwf,free_nonlop_gpu
SOURCE
172 subroutine dealloc_hamilt_gpu(option,use_gpu_cuda) 173 174 use m_profiling_abi 175 176 use defs_basis 177 #if defined HAVE_GPU_CUDA 178 use m_initcuda 179 #endif 180 181 !This section has been created automatically by the script Abilint (TD). 182 !Do not modify the following lines by hand. 183 #undef ABI_FUNC 184 #define ABI_FUNC 'dealloc_hamilt_gpu' 185 !End of the abilint section 186 187 implicit none 188 189 !Arguments ------------------------------------ 190 !scalars 191 integer,intent(in) :: option,use_gpu_cuda 192 !arrays 193 194 !Local variables------------------------------- 195 196 ! ************************************************************************* 197 198 if (use_gpu_cuda==0) return 199 200 #if defined HAVE_GPU_CUDA 201 if (option==0.or.option==2) then 202 call free_gpu_fourwf() 203 end if 204 if (option==1.or.option==2) then 205 call free_nonlop_gpu() 206 end if 207 #else 208 if (.false.) then 209 write(std_out,*) option 210 end if 211 #endif 212 213 end subroutine dealloc_hamilt_gpu