TABLE OF CONTENTS


ABINIT/alloc_hamilt_gpu [ Functions ]

[ Top ] [ Functions ]

NAME

 alloc_hamilt_gpu

FUNCTION

 allocate several memory pieces on a GPU device for the application
 of Hamiltonian using a GPU

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

 80 subroutine alloc_hamilt_gpu(atindx1,dtset,gprimd,mpi_enreg,nattyp,npwarr,option,psps,use_gpu_cuda)
 81 
 82 
 83 !This section has been created automatically by the script Abilint (TD).
 84 !Do not modify the following lines by hand.
 85 #undef ABI_FUNC
 86 #define ABI_FUNC 'alloc_hamilt_gpu'
 87 !End of the abilint section
 88 
 89  implicit none
 90 
 91 !Arguments ------------------------------------
 92 !scalars
 93  integer,intent(in) :: option,use_gpu_cuda
 94  type(dataset_type),intent(in) :: dtset
 95  type(MPI_type),intent(in) :: mpi_enreg
 96  type(pseudopotential_type),intent(in) :: psps
 97 !arrays
 98  integer,intent(in) :: atindx1(dtset%natom),nattyp(dtset%ntypat),npwarr(dtset%nkpt)
 99  real(dp),intent(in) :: gprimd(3,3)
100 
101 !Local variables-------------------------------
102 !scalars
103 #if defined HAVE_GPU_CUDA
104  integer :: dimekb1_max,dimekb2_max,dimffnl_max,ierr,ikpt,npw_max_loc,npw_max_nonloc
105  integer ::npwarr_tmp(dtset%nkpt)
106 #endif
107 
108 ! *************************************************************************
109 
110  if (use_gpu_cuda==0) return
111 
112 #if defined HAVE_GPU_CUDA
113 !=== Local Hamiltonian ===
114  if (option==0.or.option==2) then
115 !  Compute max of total planes waves
116    npw_max_loc=0
117    if(mpi_enreg%paral_kgb==1) then
118      npwarr_tmp=npwarr
119      call xmpi_sum(npwarr_tmp,mpi_enreg%comm_bandfft,ierr)
120      npw_max_loc =maxval(npwarr_tmp)
121    else
122      npw_max_loc=dtset%mpw
123    end if
124 !  Initialize gpu data needed in fourwf
125 !  ndat=bandpp when paral_kgb=1
126    if(mpi_enreg%paral_kgb==1) then
127      call alloc_gpu_fourwf(dtset%ngfft,dtset%bandpp,npw_max_loc,npw_max_loc)
128    else
129      call alloc_gpu_fourwf(dtset%ngfft,1,npw_max_loc,npw_max_loc)
130    end if
131  end if
132 !=== Nonlocal Hamiltonian ===
133  if (option==1.or.option==2) then
134 !  Compute max of total planes waves
135    npw_max_nonloc=0
136    if(mpi_enreg%paral_kgb==1) then
137      npwarr_tmp=npwarr
138      call xmpi_sum(npwarr_tmp,mpi_enreg%comm_band,ierr)
139      npw_max_nonloc =maxval(npwarr_tmp)
140    else
141      npw_max_nonloc=dtset%mpw
142    end if
143 !  Initialize all gpu data needed in nonlop
144    dimffnl_max=4
145 !  if (abs(dtset%berryopt) == 5) dimffnl_max=4
146    dimekb1_max=psps%dimekb
147    if (dtset%usepaw==0) dimekb2_max=psps%ntypat
148    if (dtset%usepaw==1) dimekb2_max=dtset%natom
149    call alloc_nonlop_gpu(npw_max_nonloc,npw_max_nonloc,dtset%nspinor,dtset%natom,&
150 &   dtset%ntypat,psps%lmnmax,psps%indlmn,nattyp,atindx1,gprimd,&
151 &   dimffnl_max,dimffnl_max,dimekb1_max,dimekb2_max)
152  end if
153  call xmpi_barrier(mpi_enreg%comm_cell)
154 #else
155  if (.false.) then
156    write(std_out,*) atindx1(1),dtset%natom,gprimd(1,1),mpi_enreg%me,nattyp(1),option
157  end if
158 #endif
159 
160 end subroutine alloc_hamilt_gpu

ABINIT/dealloc_hamilt_gpu [ Functions ]

[ Top ] [ Functions ]

NAME

 dealloc_hamilt_gpu

FUNCTION

 deallocate several memory pieces on a GPU device used for the application
 of Hamiltonian using a GPU

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

188 subroutine dealloc_hamilt_gpu(option,use_gpu_cuda)
189 
190 
191 !This section has been created automatically by the script Abilint (TD).
192 !Do not modify the following lines by hand.
193 #undef ABI_FUNC
194 #define ABI_FUNC 'dealloc_hamilt_gpu'
195 !End of the abilint section
196 
197  implicit none
198 
199 !Arguments ------------------------------------
200 !scalars
201  integer,intent(in) :: option,use_gpu_cuda
202 !arrays
203 
204 !Local variables-------------------------------
205 
206 ! *************************************************************************
207 
208  if (use_gpu_cuda==0) return
209 
210 #if defined HAVE_GPU_CUDA
211  if (option==0.or.option==2) then
212    call free_gpu_fourwf()
213  end if
214  if (option==1.or.option==2) then
215    call free_nonlop_gpu()
216  end if
217 #else
218  if (.false.) then
219    write(std_out,*) option
220  end if
221 #endif
222 
223 end subroutine dealloc_hamilt_gpu

ABINIT/m_alloc_hamilt_gpu [ Modules ]

[ Top ] [ Modules ]

NAME

  m_alloc_hamilt_gpu

FUNCTION

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 .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_alloc_hamilt_gpu
27 
28  use defs_basis
29  use defs_datatypes
30  use defs_abitypes
31  use m_abicore
32  use m_xmpi
33 #if defined HAVE_GPU_CUDA
34  use m_initcuda
35 #endif
36 
37  implicit none
38 
39  private