TABLE OF CONTENTS


ABINIT/bound [ Functions ]

[ Top ] [ Functions ]

NAME

 bound

FUNCTION

 For given kpt, ngfft, and gmet,
  Find distance**2 to boundary point of fft box nearest to kpt
  Find distance**2 to boundary point of fft box farthest to kpt

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG)
 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

  kpt(3)=real input k vector (reduced coordinates)
  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
  gmet(3,3)=reciprocal space metric (currently in Bohr**-2)

OUTPUT

  dsqmax=maximum distance**2 from k to boundary in Bohr**-2.
  dsqmin=minimum distance**2 from k to boundary in Bohr**-2.
  gbound(3)=coords of G on boundary (correspnding to gsqmin)
  plane=which plane min occurs in (1,2, or 3 for G1,etc).

SIDE EFFECTS

NOTES

 Potential trouble: this routine was written assuming kpt lies inside
 first Brillouin zone.  No measure is taken to fold input kpt back
 into first zone.  Given arbitrary kpt, this will cause trouble.

PARENTS

      getcut,getng

CHILDREN

SOURCE

 45 #if defined HAVE_CONFIG_H
 46 #include "config.h"
 47 #endif
 48 
 49 #include "abi_common.h"
 50 
 51 
 52 subroutine bound(dsqmax,dsqmin,gbound,gmet,kpt,ngfft,plane)
 53 
 54  use defs_basis
 55  use m_profiling_abi
 56  use m_errors
 57 
 58 !This section has been created automatically by the script Abilint (TD).
 59 !Do not modify the following lines by hand.
 60 #undef ABI_FUNC
 61 #define ABI_FUNC 'bound'
 62 !End of the abilint section
 63 
 64  implicit none
 65 
 66 !Arguments ------------------------------------
 67 !scalars
 68  integer,intent(out) :: plane
 69  real(dp),intent(out) :: dsqmax,dsqmin
 70 !arrays
 71  integer,intent(in) :: ngfft(18)
 72  integer,intent(out) :: gbound(3)
 73  real(dp),intent(in) :: gmet(3,3),kpt(3)
 74 
 75 !Local variables-------------------------------
 76 !scalars
 77  integer :: i1,i1min,i2,i2min,i3,i3min
 78  real(dp) :: dsm,dsp
 79  character(len=500) :: message
 80 
 81 ! *************************************************************************
 82 
 83 ! dsq(i1,i2,i3)=gmet(1,1)*(kpt(1)+dble(i1))**2&
 84 !& +gmet(2,2)*(kpt(2)+dble(i2))**2&
 85 !& +gmet(3,3)*(kpt(3)+dble(i3))**2&
 86 !& +2._dp*(gmet(1,2)*(kpt(1)+dble(i1))*(kpt(2)+dble(i2))&
 87 !& +gmet(2,3)*(kpt(2)+dble(i2))*(kpt(3)+dble(i3))&
 88 !& +gmet(3,1)*(kpt(3)+dble(i3))*(kpt(1)+dble(i1)))
 89 
 90 !Set plane to impossible value
 91  plane=0
 92 
 93 !look at +/- g1 planes:
 94  dsqmax=zero
 95  dsqmin=dsq(ngfft(1)/2,-ngfft(2)/2,-ngfft(3)/2,gmet,kpt)+0.01_dp
 96  do i2=-ngfft(2)/2,ngfft(2)/2
 97    do i3=-ngfft(3)/2,ngfft(3)/2
 98      dsp = dsq(ngfft(1)/2, i2, i3,gmet,kpt)
 99      dsm = dsq( - ngfft(1)/2, i2, i3,gmet,kpt)
100      if (dsp>dsqmax) dsqmax = dsp
101      if (dsm>dsqmax) dsqmax = dsm
102      if (dsp<dsqmin) then
103        dsqmin = dsp
104        i1min = ngfft(1)/2
105        i2min = i2
106        i3min = i3
107        plane=1
108      end if
109      if (dsm<dsqmin) then
110        dsqmin = dsm
111        i1min =  - ngfft(1)/2
112        i2min = i2
113        i3min = i3
114        plane=1
115      end if
116    end do
117  end do
118 !
119 !+/- g2 planes:
120  do i1=-ngfft(1)/2,ngfft(1)/2
121    do i3=-ngfft(3)/2,ngfft(3)/2
122      dsp = dsq(i1,ngfft(2)/2,i3,gmet,kpt)
123      dsm = dsq(i1,-ngfft(2)/2,i3,gmet,kpt)
124      if (dsp>dsqmax) dsqmax = dsp
125      if (dsm>dsqmax) dsqmax = dsm
126      if (dsp<dsqmin) then
127        dsqmin = dsp
128        i1min = i1
129        i2min = ngfft(2)/2
130        i3min = i3
131        plane=2
132      end if
133      if (dsm<dsqmin) then
134        dsqmin = dsm
135        i1min = i1
136        i2min =  - ngfft(2)/2
137        i3min = i3
138        plane=2
139      end if
140    end do
141  end do
142 !
143 !+/- g3 planes:
144  do i1=-ngfft(1)/2,ngfft(1)/2
145    do i2=-ngfft(2)/2,ngfft(2)/2
146      dsp = dsq(i1,i2,ngfft(3)/2,gmet,kpt)
147      dsm = dsq(i1,i2,-ngfft(3)/2,gmet,kpt)
148      if (dsp>dsqmax) dsqmax = dsp
149      if (dsm>dsqmax) dsqmax = dsm
150      if (dsp<dsqmin) then
151        dsqmin = dsp
152        i1min = i1
153        i2min = i2
154        i3min = ngfft(3)/2
155        plane=3
156      end if
157      if (dsm<dsqmin) then
158        dsqmin = dsm
159        i1min = i1
160        i2min = i2
161        i3min =  - ngfft(3)/2
162        plane=3
163      end if
164    end do
165  end do
166 
167  if (plane==0) then
168 !  Trouble: missed boundary somehow
169    write(message, '(a,a,a,3f9.4,a,3i5,a,a,a,a,a)' )&
170 &   'Trouble finding boundary of G sphere for',ch10,&
171 &   'kpt=',kpt(:),' and ng=',ngfft(1:3),ch10,&
172 &   'Action : check that kpt lies',&
173 &   'reasonably within first Brillouin zone; ',ch10,&
174 &   'else code bug, contact ABINIT group.'
175    MSG_BUG(message)
176  end if
177 
178  gbound(1)=i1min
179  gbound(2)=i2min
180  gbound(3)=i3min
181 
182  contains
183 
184    function dsq(i1,i2,i3,gmet,kpt)
185 
186 
187 !This section has been created automatically by the script Abilint (TD).
188 !Do not modify the following lines by hand.
189 #undef ABI_FUNC
190 #define ABI_FUNC 'dsq'
191 !End of the abilint section
192 
193      integer :: i1,i2,i3
194      real(dp) :: dsq
195      real(dp) :: kpt(3),gmet(3,3)
196 
197      dsq=gmet(1,1)*(kpt(1)+dble(i1))**2&
198 &      +gmet(2,2)*(kpt(2)+dble(i2))**2&
199 &      +gmet(3,3)*(kpt(3)+dble(i3))**2&
200 &      +2._dp*(gmet(1,2)*(kpt(1)+dble(i1))*(kpt(2)+dble(i2))&
201 &      +gmet(2,3)*(kpt(2)+dble(i2))*(kpt(3)+dble(i3))&
202 &      +gmet(3,1)*(kpt(3)+dble(i3))*(kpt(1)+dble(i1)))
203    end function dsq
204 
205 end subroutine bound