TABLE OF CONTENTS
ABINIT/getmpw [ Functions ]
NAME
getmpw
FUNCTION
From input ecut, combined with ucvol and gmet, compute recommended mpw mpw is the maximum number of plane-waves in the wave-function basis for one processor of the WF group
COPYRIGHT
Copyright (C) 1998-2017 ABINIT group (DCA, XG, GMR) 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
ecut=plane-wave cutoff energy in Hartrees exchn2n3d=if n1, n2 and n3 are exchanged gmet(3,3)=reciprocal space metric (bohr**-2). istwfk(nkpt)=input option parameter that describes the storage of wfs kptns(3,nkpt)=real(dp) array for k points (normalisation is already taken into account) mpi_enreg=information about MPI parallelization nkpt=integer number of k points in the calculation
OUTPUT
mpw=maximal number of plane waves over all k points of the processor (for one processor of the WF group)
PARENTS
dfpt_looppert,memory_eval,mpi_setup,scfcv
CHILDREN
kpgsph,wrtout
SOURCE
39 #if defined HAVE_CONFIG_H 40 #include "config.h" 41 #endif 42 43 #include "abi_common.h" 44 45 subroutine getmpw(ecut,exchn2n3d,gmet,istwfk,kptns,mpi_enreg,mpw,nkpt) 46 47 use defs_basis 48 use defs_abitypes 49 use m_profiling_abi 50 use m_errors 51 52 use m_fftcore, only : kpgsph 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'getmpw' 58 use interfaces_14_hidewrite 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 integer,intent(in) :: exchn2n3d,nkpt 66 integer,intent(out) :: mpw 67 real(dp),intent(in) :: ecut 68 type(MPI_type),intent(inout) :: mpi_enreg 69 !arrays 70 integer,intent(in) :: istwfk(nkpt) 71 real(dp),intent(in) :: gmet(3,3),kptns(3,nkpt) 72 73 !Local variables------------------------------- 74 !scalars 75 integer :: ikpt,istwf_k,npw 76 ! integer :: npwwrk,pad=50 77 ! real(dp) :: scale=1.3_dp 78 character(len=500) :: message 79 !arrays 80 integer,allocatable :: kg(:,:) 81 real(dp) :: kpoint(3) 82 83 ! ************************************************************************* 84 85 !An upper bound for mpw, might be obtained as follows 86 !the average number of plane-waves in the cutoff sphere is 87 !npwave = (2*ecut)**(3/2) * ucvol / (6*pi**2) 88 !the upper bound is calculated as 89 !npwwrk = int(scale * npwave) + pad 90 !rescale so an upper bound 91 !npwave=nint(ucvol*(2.0_dp*ecut)**1.5_dp/(6.0_dp*pi**2)) 92 !npwwrk=nint(dble(npwave)*scale)+pad 93 94 ABI_ALLOCATE(kg,(3,100)) 95 96 !set mpw to zero, as needed for only counting in kpgsph 97 mpw = 0 98 99 !Might be parallelized over k points ? ! 100 do ikpt = 1,nkpt 101 ! Do computation of G sphere, returning npw 102 kpoint(:)=kptns(:,ikpt) 103 istwf_k=istwfk(ikpt) 104 call kpgsph(ecut,exchn2n3d,gmet,0,ikpt,istwf_k,kg,kpoint,0,mpi_enreg,0,npw) 105 mpw = max(npw,mpw) 106 end do 107 108 write(message,'(a,i0)') ' getmpw: optimal value of mpw= ',mpw 109 call wrtout(std_out,message,'COLL') 110 111 ABI_DEALLOCATE(kg) 112 113 end subroutine getmpw