TABLE OF CONTENTS


ABINIT/setmqgrid [ Functions ]

[ Top ] [ Functions ]

NAME

  setmqgrid

FUNCTION

  Sets the number of points needed to represent the pseudopotentials in
  reciprocal space for a specified resolution.

COPYRIGHT

  Copyright (C) 2011-2018 ABINIT group (DW)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  ecut=cutoff energy for the wavefunctions
  ecutdg=cutoff energy for the fine grid in case usepaw==1
  gprimd=primitive translation vectors for reciprocal space
  nptsgvec=number of points along the smallest primitive translation vector
    of the reciprocal space
  usepaw=1 if PAW is used, 0 otherwise

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      m_psps,memory_eval

CHILDREN

SOURCE

 37 #if defined HAVE_CONFIG_H
 38 #include "config.h"
 39 #endif
 40 
 41 #include "abi_common.h"
 42 
 43 
 44 subroutine setmqgrid(mqgrid,mqgriddg,ecut,ecutdg,gprimd,nptsgvec,usepaw)
 45 
 46  use defs_basis
 47  use m_profiling_abi
 48  use m_errors
 49 
 50 !This section has been created automatically by the script Abilint (TD).
 51 !Do not modify the following lines by hand.
 52 #undef ABI_FUNC
 53 #define ABI_FUNC 'setmqgrid'
 54 !End of the abilint section
 55 
 56  implicit none
 57 
 58 !Arguments ------------------------------------
 59  integer , intent(inout)  :: mqgrid,mqgriddg
 60  integer , intent(in)  :: nptsgvec,usepaw
 61  real(dp), intent(in) :: ecut,ecutdg
 62  real(dp), intent(in) :: gprimd(3,3)
 63 
 64 !Local variables-------------------------------
 65  integer :: mqgrid2,mqgriddg2
 66  real(dp) :: gmax,gmaxdg,gvecnorm
 67  character(len=500) :: message  
 68  
 69 ! *************************************************************************
 70  
 71  gvecnorm=sqrt(min(dot_product(gprimd(:,1),gprimd(:,1)), &
 72 & dot_product(gprimd(:,2),gprimd(:,2)), &
 73 & dot_product(gprimd(:,3),gprimd(:,3))))
 74  gmax=one/(sqrt2*pi)*sqrt(ecut)
 75 
 76  if (mqgrid == 0) then
 77    mqgrid2=ceiling(gmax/gvecnorm*nptsgvec)
 78    mqgrid=max(mqgrid2,3001)
 79    write(message, '(5a,i0,a)' )&
 80 &   'The number of points "mqgrid" in reciprocal space used for the',ch10,&
 81 &   'description of the pseudopotentials has been set automatically',ch10,&
 82 &   'by abinit to: ',mqgrid,'.'
 83    !MSG_COMMENT(message)
 84  else
 85    mqgrid2=ceiling(gmax/gvecnorm*nptsgvec)
 86    if (mqgrid2>mqgrid) then
 87      write(message, '(3a,i8,3a,i8,3a)' )&
 88 &     'The number of points "mqgrid" in reciprocal space used for the',ch10,&
 89 &     'description of the pseudopotentials is : ',mqgrid,'.',ch10,&
 90 &     'It would be better to increase it to at least ',mqgrid2,', or',ch10,&
 91 &     'let abinit choose it automatically by setting mqgrid = 0.'
 92      MSG_WARNING(message)
 93    end if
 94  end if
 95 
 96  if (usepaw==1) then
 97    if(ecutdg<tol6)then
 98      write(message,'(a)')&
 99 &     'The value of (paw)ecutdg is zero or negative, which is forbidden.'
100      MSG_ERROR(message)
101    end if
102    gmaxdg=one/(sqrt2*pi)*sqrt(ecutdg)
103    if (mqgriddg == 0) then
104      mqgriddg2=ceiling(gmaxdg/gvecnorm*nptsgvec)
105      mqgriddg=max(mqgriddg2,3001)
106      write(message, '(5a,i0,a)' )&
107 &     'The number of points "mqgriddg" in reciprocal space used for the',ch10,&
108 &     'description of the pseudopotentials has been set automatically',ch10,&
109 &     'by abinit to: ',mqgriddg,'.'
110      !MSG_COMMENT(message)
111    else
112      mqgriddg2=ceiling(gmax/gvecnorm*nptsgvec)
113      if (mqgriddg2>mqgriddg) then
114        write(message, '(3a,i8,3a,i8,3a)' )&
115 &       'The number of points "mqgriddg" in reciprocal space used for the',ch10,&
116 &       'description of the pseudopotentials (fine grid) is :',mqgriddg,'.',ch10,&
117 &       'It would be better to increase it to at least ',mqgriddg2,', or',ch10,&
118 &       'let abinit choose it automatically by setting mqgrid = 0.'
119        MSG_WARNING(message)
120      end if
121    end if
122  end if
123 
124 end subroutine setmqgrid