TABLE OF CONTENTS


ABINIT/m_vcoul/m_cutoff_sphere [ Modules ]

[ Top ] [ Modules ]

NAME

  m_cutoff_sphere

FUNCTION

SOURCE

 9 #if defined HAVE_CONFIG_H
10 #include "config.h"
11 #endif
12 
13 
14 module m_cutoff_sphere
15 
16  use defs_basis
17  use m_abicore
18  use m_geometry,        only : normv
19 
20  implicit none
21 
22  private
23 
24  public ::  cutoff_sphere

m_cutoff/cutoff_sphere [ Functions ]

[ Top ] [ Functions ]

NAME

 cutoff_sphere

FUNCTION

  Calculate the Fourier transform of the Coulomb interaction with a spherical cutoff in real-space.

   $ v_{cut}(G)= \frac{4\pi}{|q+G|^2} [ 1-cos(|q+G|*R_cut) ] $  (1)

  For |q|<small and G=0 we use 2pi.R_cut^2, namely we consider the limit q-->0 of Eq. (1)

INPUTS

  qpt(3)=q-point where the cutoff Coulomb is required.
  ngvec=Number of G vectors
  gvec(3,ngvec)=G vectors in reduced coordinates.
  gmet(3,3)=Metric in reciprocal space.
  rcut=Cutoff radius of the sphere.

OUTPUT

  vc_cut(ngvec)=Fourier components of the effective Coulomb interaction.

SOURCE

 53 subroutine cutoff_sphere(qpt, ngvec, gvec, gmet, rcut, vc_cut)
 54 
 55 !Arguments ------------------------------------
 56 !scalars
 57  integer,intent(in) :: ngvec
 58  real(dp),intent(in) :: rcut
 59 !arrays
 60  integer,intent(in) :: gvec(3,ngvec)
 61  real(dp),intent(in) :: gmet(3,3),qpt(3)
 62  real(dp),intent(out) :: vc_cut(ngvec)
 63 
 64 !Local variables-------------------------------
 65 !scalars
 66  integer :: ig,igs
 67  real(dp) :: qpg
 68  logical :: ltest
 69 
 70 !************************************************************************
 71 
 72  !------------------------------------------------------------
 73  ! Code modification by Bruno Rousseau, Montreal, 06/11/2013
 74  !------------------------------------------------------------
 75  !
 76  ! In order for the code below to run in parallel using MPI
 77  ! within the gwls code (by Laflamme-Jansen, Cote and Rousseau)
 78  ! the ABI_CHECK below must be removed.
 79  !
 80  ! The ABI_CHECK below will fail if G-vectors are distributed on
 81  ! many processors, as only the master process has G=(0,0,0).
 82  !
 83  ! The test does not seem necessary; IF a process has G=(0,0,0)
 84  ! (namely, the master process), then it will be the first G vector.
 85  ! If a process does not have G=(0,0,0) (ie, all the other processes),
 86  ! then they don't need to worry about the G-> 0 limit.
 87  !
 88 
 89  ltest = ALL(gvec(:,1) == 0)
 90  !ABI_CHECK(ltest,'The first G vector should be Gamma')
 91 
 92  igs=1
 93  if (ltest .and. normv(qpt, gmet, 'G') < tol4) then ! For small q and G=0, use the limit q-->0.
 94    vc_cut(1)=two_pi*rcut**2
 95    igs=2
 96  end if
 97  do ig=igs,ngvec
 98    qpg = normv(qpt(:) + gvec(:,ig), gmet, 'G')
 99    vc_cut(ig) = four_pi* (one-COS(rcut*qpg)) / qpg**2
100  end do
101 
102 end subroutine cutoff_sphere