TABLE OF CONTENTS
ABINIT/multithreaded_getghc [ Functions ]
NAME
multithreaded_getghc
FUNCTION
COPYRIGHT
Copyright (C) 2016-2018 ABINIT group (JB) 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
cpopt=flag defining the status of cwaveprj%cp(:)=<Proj_i|Cnk> scalars (PAW only) (same meaning as in nonlop.F90 routine) if cpopt=-1, <p_lmn|in> (and derivatives) are computed here (and not saved) if cpopt= 0, <p_lmn|in> are computed here and saved if cpopt= 1, <p_lmn|in> and first derivatives are computed here and saved if cpopt= 2 <p_lmn|in> are already in memory; if cpopt= 3 <p_lmn|in> are already in memory; first derivatives are computed here and saved if cpopt= 4 <p_lmn|in> and first derivatives are already in memory; cwavef(2,npw*my_nspinor*ndat)=planewave coefficients of wavefunction. gs_ham <type(gs_hamiltonian_type)>=all data for the Hamiltonian to be applied lambda=factor to be used when computing <G|H-lambda.S|C> - only for sij_opt=-1 Typically lambda is the eigenvalue (or its guess) mpi_enreg=informations about MPI parallelization ndat=number of FFT to do in parallel prtvol=control print volume and debugging output sij_opt= -PAW ONLY- if 0, only matrix elements <G|H|C> have to be computed (S=overlap) if 1, matrix elements <G|S|C> have to be computed in gsc in addition to ghc if -1, matrix elements <G|H-lambda.S|C> have to be computed in ghc (gsc not used) tim_getghc=timing code of the calling subroutine(can be set to 0 if not attributed) type_calc= option governing which part of Hamitonian is to be applied: 1: local part only 2: non-local+kinetic only (added to the exixting Hamiltonian) 3: local + kinetic only (added to the existing Hamiltonian) ===== Optional inputs ===== [kg_fft_k(3,:)]=optional, (k+G) vector coordinates to be used for the FFT tranformation instead of the one contained in gs_ham datastructure. Typically used for real WF (in parallel) which are FFT-transformed 2 by 2. [kg_fft_kp(3,:)]=optional, (k^prime+G) vector coordinates to be used for the FFT tranformation [select_k]=optional, option governing the choice of k points to be used. gs_ham datastructure contains quantities needed to apply Hamiltonian in reciprocal space between 2 kpoints, k and k^prime (equal in most cases); if select_k=1, <k^prime|H|k> is applied [default] if select_k=2, <k|H|k^prime> is applied if select_k=3, <k|H|k> is applied if select_k=4, <k^prime|H|k^prime> is applied
OUTPUT
ghc(2,npw*my_nspinor*ndat)=matrix elements <G|H|C> (if sij_opt>=0) or <G|H-lambda.S|C> (if sij_opt=-1) gvnlc(2,npw*my_nspinor*ndat)=matrix elements <G|Vnonlocal|C> (if sij_opt>=0) or <G|Vnonlocal-lambda.S|C> (if sij_opt=-1) if (sij_opt=1) gsc(2,npw*my_nspinor*ndat)=matrix elements <G|S|C> (S=overlap).
SIDE EFFECTS
cwaveprj(natom,my_nspinor*(1+cpopt)*ndat)= wave function projected on nl projectors (PAW only)
PARENTS
m_lobpcgwf,prep_getghc
CHILDREN
getghc,mkl_set_num_threads,omp_set_nested
SOURCE
73 #if defined HAVE_CONFIG_H 74 #include "config.h" 75 #endif 76 77 #include "abi_common.h" 78 79 subroutine multithreaded_getghc(cpopt,cwavef,cwaveprj,ghc,gsc,gs_ham,gvnlc,lambda,mpi_enreg,ndat,& 80 & prtvol,sij_opt,tim_getghc,type_calc,& 81 & kg_fft_k,kg_fft_kp,select_k) ! optional arguments 82 83 use defs_basis 84 use defs_abitypes 85 use m_errors 86 use m_profiling_abi 87 use m_xmpi 88 89 use m_pawcprj, only : pawcprj_type,pawcprj_alloc,pawcprj_free,pawcprj_getdim 90 use m_bandfft_kpt, only : bandfft_kpt,bandfft_kpt_get_ikpt 91 use m_hamiltonian, only : gs_hamiltonian_type,KPRIME_H_K,K_H_KPRIME,K_H_K,KPRIME_H_KPRIME 92 use m_fock, only : fock_type,fock_get_getghc_call 93 94 #ifdef HAVE_OPENMP 95 use omp_lib 96 #endif 97 98 !This section has been created automatically by the script Abilint (TD). 99 !Do not modify the following lines by hand. 100 #undef ABI_FUNC 101 #define ABI_FUNC 'multithreaded_getghc' 102 use interfaces_66_wfs, except_this_one => multithreaded_getghc 103 !End of the abilint section 104 105 implicit none 106 107 !Arguments ------------------------------------ 108 !scalars 109 integer,intent(in) :: cpopt,ndat, prtvol 110 integer,intent(in) :: sij_opt,tim_getghc,type_calc 111 integer,intent(in),optional :: select_k 112 real(dp),intent(in) :: lambda 113 type(MPI_type),intent(in) :: mpi_enreg 114 type(gs_hamiltonian_type),intent(inout),target :: gs_ham 115 !arrays 116 integer,intent(in),optional,target :: kg_fft_k(:,:),kg_fft_kp(:,:) 117 real(dp),intent(out),target :: gsc(:,:) 118 real(dp),intent(inout) :: cwavef(:,:) 119 real(dp),intent(out) :: ghc(:,:),gvnlc(:,:) 120 type(pawcprj_type),intent(inout),target :: cwaveprj(:,:) 121 122 !Local variables------------------------------- 123 !scalars 124 integer :: firstelt, lastelt 125 integer :: nthreads 126 integer :: ithread 127 integer :: chunk 128 integer :: residuchunk 129 integer :: firstband 130 integer :: lastband 131 integer :: spacedim 132 #ifdef HAVE_OPENMP 133 logical :: is_nested 134 #endif 135 136 integer :: select_k_default 137 138 ! ************************************************************************* 139 140 select_k_default = 1; if ( present(select_k) ) select_k_default = select_k 141 142 spacedim = size(cwavef,dim=2)/ndat 143 144 !$omp parallel default (none) private(ithread,nthreads,chunk,firstband,lastband,residuchunk,firstelt,lastelt, is_nested), & 145 !$omp& shared(cwavef,ghc,gsc, gvnlc,spacedim,ndat,kg_fft_k,kg_fft_kp,gs_ham,cwaveprj,mpi_enreg), & 146 !$omp& firstprivate(cpopt,lambda,prtvol,sij_opt,tim_getghc,type_calc,select_k_default) 147 #ifdef HAVE_OPENMP 148 ithread = omp_get_thread_num() 149 nthreads = omp_get_num_threads() 150 is_nested = omp_get_nested() 151 call omp_set_nested(.false.) 152 #ifdef HAVE_LINALG_MKL_THREADS 153 call mkl_set_num_threads(1) 154 #endif 155 #else 156 ithread = 0 157 nthreads = 1 158 #endif 159 chunk = ndat/nthreads ! Divide by 2 to construct chunk of even number of bands 160 residuchunk = ndat - nthreads*chunk 161 if ( ithread < nthreads-residuchunk ) then 162 firstband = ithread*chunk+1 163 lastband = (ithread+1)*chunk 164 else 165 firstband = (nthreads-residuchunk)*chunk + ( ithread -(nthreads-residuchunk) )*(chunk+1) +1 166 lastband = firstband+chunk 167 end if 168 169 if ( lastband /= 0 ) then 170 firstelt = (firstband-1)*spacedim+1 171 lastelt = lastband*spacedim 172 ! Don't know how to manage optional arguments .... :( 173 if ( present(kg_fft_k) ) then 174 if (present(kg_fft_kp)) then 175 call getghc(cpopt,cwavef(:,firstelt:lastelt),cwaveprj,ghc(:,firstelt:lastelt),gsc(:,firstelt:lastelt*gs_ham%usepaw),& 176 gs_ham,gvnlc(:,firstelt:lastelt),lambda, mpi_enreg,lastband-firstband+1,prtvol,sij_opt,tim_getghc,type_calc,& 177 select_k=select_k_default,kg_fft_k=kg_fft_k,kg_fft_kp=kg_fft_kp) 178 else 179 call getghc(cpopt,cwavef(:,firstelt:lastelt),cwaveprj,ghc(:,firstelt:lastelt),gsc(:,firstelt:lastelt*gs_ham%usepaw),& 180 gs_ham,gvnlc(:,firstelt:lastelt),lambda, mpi_enreg,lastband-firstband+1,prtvol,sij_opt,tim_getghc,type_calc,& 181 select_k=select_k_default,kg_fft_k=kg_fft_k) 182 end if 183 else 184 if (present(kg_fft_kp)) then 185 call getghc(cpopt,cwavef(:,firstelt:lastelt),cwaveprj,ghc(:,firstelt:lastelt),gsc(:,firstelt:lastelt*gs_ham%usepaw),& 186 gs_ham,gvnlc(:,firstelt:lastelt),lambda, mpi_enreg,lastband-firstband+1,prtvol,sij_opt,tim_getghc,type_calc,& 187 select_k=select_k_default,kg_fft_kp=kg_fft_kp) 188 else 189 call getghc(cpopt,cwavef(:,firstelt:lastelt),cwaveprj,ghc(:,firstelt:lastelt),gsc(:,firstelt:lastelt*gs_ham%usepaw),& 190 gs_ham,gvnlc(:,firstelt:lastelt),lambda, mpi_enreg,lastband-firstband+1,prtvol,sij_opt,tim_getghc,type_calc,& 191 select_k=select_k_default) 192 end if 193 end if 194 end if 195 #ifdef HAVE_OPENMP 196 call omp_set_nested(is_nested) 197 #ifdef HAVE_LINALG_MKL_THREADS 198 call mkl_set_num_threads(nthreads) 199 #endif 200 #endif 201 !$omp end parallel 202 203 end subroutine multithreaded_getghc