TABLE OF CONTENTS


ABINIT/m_ab7_kpoints [ Modules ]

[ Top ] [ Modules ]

NAME

FUNCTION

COPYRIGHT

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

PARENTS

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 module m_ab7_kpoints
24 
25   use defs_basis
26   use m_ab7_symmetry
27   use m_abicore
28 
29   use m_kpts,      only : getkgrid, testkgrid
30   use m_spacepar,  only : irrzg
31 
32   implicit none
33 
34   private
35 
36   logical, private, parameter :: AB_DBG = .false.
37 
38   public :: kpoints_get_irreductible_zone
39 
40   public :: kpoints_get_mp_k_grid
41   public :: kpoints_get_auto_k_grid
42 
43   public :: kpoints_binding_mp_k_1
44   public :: kpoints_binding_mp_k_2
45   public :: kpoints_binding_auto_k_1
46   public :: kpoints_binding_auto_k_2
47 
48 contains

m_ab7_kpoints/kpoints_binding_auto_k_1 [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

  kpoints_binding_auto_k_1

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      m_ab7_kpoints

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

311   subroutine kpoints_binding_auto_k_1(symid, nkpt, kptrlatt, kptrlen, &
312        & nshiftk, shiftk, errno)
313 
314 
315 !This section has been created automatically by the script Abilint (TD).
316 !Do not modify the following lines by hand.
317 #undef ABI_FUNC
318 #define ABI_FUNC 'kpoints_binding_auto_k_1'
319 !End of the abilint section
320 
321     integer, intent(in)  :: symid
322     integer, intent(out) :: errno
323     integer, intent(out) :: nkpt
324     real(dp), intent(inout) :: kptrlen
325     integer, intent(out) :: nshiftk
326     real(dp), intent(out) :: shiftk(3,210)
327     integer, intent(out) :: kptrlatt(3,3)
328 
329     type(symmetry_type), pointer  :: sym
330     real(dp), allocatable :: kpt(:,:), wkpt(:)
331 
332     if (AB_DBG) write(std_err,*) "AB symmetry: call get auto k grid1."
333 
334     errno = AB7_NO_ERROR
335     call symmetry_get_from_id(sym, symid, errno)
336     if (errno /= AB7_NO_ERROR) return
337 
338     !  The parameters of the k lattice are not known, compute
339     !  kptrlatt, nshiftk, shiftk.
340     call testkgrid(sym%bravais,6,kptrlatt,kptrlen,&
341          & AB7_MAX_SYMMETRIES,nshiftk,sym%nSym,0,sym%rprimd,&
342          & shiftk,sym%symAfm,sym%sym,sym%vacuum)
343     if (AB_DBG) write(std_err,*) "AB symmetry: testkgrid -> kptrlatt=", kptrlatt
344 
345     nkpt=0
346     ABI_ALLOCATE(kpt,(3, nkpt))
347     ABI_ALLOCATE(wkpt,(nkpt))
348 
349     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
350          & AB7_MAX_SYMMETRIES, 0, nkpt, nshiftk, sym%nSym, &
351          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
352          & sym%vacuum, wkpt)
353     if (AB_DBG) write(std_err,*) "AB symmetry: getkgrid -> nkpt=", nkpt
354 
355     ABI_DEALLOCATE(kpt)
356     ABI_DEALLOCATE(wkpt)
357 
358   end subroutine kpoints_binding_auto_k_1

m_ab7_kpoints/kpoints_binding_auto_k_2 [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

  kpoints_binding_auto_k_2

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      m_ab7_kpoints

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

383   subroutine kpoints_binding_auto_k_2(symid, nkpt, kpt, wkpt, kptrlatt, kptrlen, &
384        & nshiftk, shiftk, errno)
385 
386 
387 !This section has been created automatically by the script Abilint (TD).
388 !Do not modify the following lines by hand.
389 #undef ABI_FUNC
390 #define ABI_FUNC 'kpoints_binding_auto_k_2'
391 !End of the abilint section
392 
393     integer, intent(in)  :: symid
394     integer, intent(out) :: errno
395     integer, intent(in) :: nkpt
396     real(dp), intent(out) :: kpt(3,nkpt), wkpt(nkpt)
397     real(dp), intent(inout) :: kptrlen
398     integer, intent(inout) :: nshiftk
399     real(dp), intent(inout) :: shiftk(3,210)
400     integer, intent(inout) :: kptrlatt(3,3)
401 
402     type(symmetry_type), pointer  :: sym
403     integer :: nkpt_
404 
405     if (AB_DBG) write(std_err,*) "AB symmetry: call get auto k grid2."
406 
407     errno = AB7_NO_ERROR
408     call symmetry_get_from_id(sym, symid, errno)
409     if (errno /= AB7_NO_ERROR) return
410 
411     ! Then, we call it again to get the actual values for the k points.
412     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
413          & AB7_MAX_SYMMETRIES, nkpt, nkpt_, nshiftk, sym%nSym, &
414          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
415          & sym%vacuum, wkpt)
416   end subroutine kpoints_binding_auto_k_2

m_ab7_kpoints/kpoints_binding_mp_k_1 [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

  kpoints_binding_mp_k_1

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      m_ab7_kpoints

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

130   subroutine kpoints_binding_mp_k_1(symid, nkpt, ngkpt, &
131        & kptrlatt, kptrlen, nshiftk, shiftk, errno)
132 
133 
134 !This section has been created automatically by the script Abilint (TD).
135 !Do not modify the following lines by hand.
136 #undef ABI_FUNC
137 #define ABI_FUNC 'kpoints_binding_mp_k_1'
138 !End of the abilint section
139 
140     integer, intent(in)  :: symid
141     integer, intent(out) :: errno
142     integer, intent(in) :: ngkpt(3)
143     integer, intent(inout) :: nshiftk
144     real(dp), intent(inout) :: shiftk(3,210)
145     real(dp), intent(out) :: kptrlen
146     integer, intent(out) :: kptrlatt(3,3)
147     integer, intent(out) :: nkpt
148 
149     type(symmetry_type), pointer  :: sym
150     real(dp), allocatable :: kpt(:,:), wkpt(:)
151 
152     if (AB_DBG) write(std_err,*) "AB symmetry: call get k grid1."
153 
154     errno = AB7_NO_ERROR
155     call symmetry_get_from_id(sym, symid, errno)
156     if (errno /= AB7_NO_ERROR) return
157 
158     ! First, compute the number of kpoints
159     kptrlatt(:,:) = 0
160     kptrlatt(1,1) = ngkpt(1)
161     kptrlatt(2,2) = ngkpt(2)
162     kptrlatt(3,3) = ngkpt(3)
163     kptrlen = 20.
164 
165     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
166          & AB7_MAX_SYMMETRIES, 0, nkpt, nshiftk, sym%nSym, &
167          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
168          & sym%vacuum, wkpt)
169   end subroutine kpoints_binding_mp_k_1

m_ab7_kpoints/kpoints_binding_mp_k_2 [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

  kpoints_binding_mp_k_2

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      m_ab7_kpoints

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

194   subroutine kpoints_binding_mp_k_2(symid, nkpt, kpt, wkpt, &
195        & kptrlatt, kptrlen, nshiftk, shiftk, errno)
196 
197 
198 !This section has been created automatically by the script Abilint (TD).
199 !Do not modify the following lines by hand.
200 #undef ABI_FUNC
201 #define ABI_FUNC 'kpoints_binding_mp_k_2'
202 !End of the abilint section
203 
204     integer, intent(in)  :: symid
205     integer, intent(out) :: errno
206     integer, intent(inout) :: nshiftk
207     real(dp), intent(inout) :: shiftk(3,210)
208     integer, intent(in) :: nkpt
209     real(dp), intent(out) :: kpt(3,nkpt), wkpt(nkpt)
210     real(dp), intent(inout) :: kptrlen
211     integer, intent(inout) :: kptrlatt(3,3)
212 
213     type(symmetry_type), pointer  :: sym
214     integer :: nkpt_
215 
216     if (AB_DBG) write(std_err,*) "AB symmetry: call get k grid2."
217 
218     errno = AB7_NO_ERROR
219     call symmetry_get_from_id(sym, symid, errno)
220     if (errno /= AB7_NO_ERROR) return
221 
222     ! Then, we call it again to get the actual values for the k points.
223     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
224          & AB7_MAX_SYMMETRIES, nkpt, nkpt_, nshiftk, sym%nSym, &
225          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
226          & sym%vacuum, wkpt)
227   end subroutine kpoints_binding_mp_k_2

m_ab7_kpoints/kpoints_get_auto_k_grid [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

  kpoints_get_auto_k_grid

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

440   subroutine kpoints_get_auto_k_grid(symid, nkpt, kpt, wkpt, &
441        & kptrlen, errno)
442 
443 
444 !This section has been created automatically by the script Abilint (TD).
445 !Do not modify the following lines by hand.
446 #undef ABI_FUNC
447 #define ABI_FUNC 'kpoints_get_auto_k_grid'
448 !End of the abilint section
449 
450     integer, intent(in)  :: symid
451     integer, intent(out) :: errno
452     integer, intent(out) :: nkpt
453     real(dp), intent(in) :: kptrlen
454     real(dp), pointer :: kpt(:,:), wkpt(:)
455 
456     real(dp) :: kptrlen_
457     integer :: kptrlatt(3,3)
458     integer :: nshiftk
459     real(dp) :: shiftk(3,210)
460 
461     if (AB_DBG) write(std_err,*) "AB symmetry: call get auto k grid."
462 
463     kptrlen_ = kptrlen
464     call kpoints_binding_auto_k_1(symid, nkpt, kptrlatt, kptrlen_, &
465        & nshiftk, shiftk, errno)
466     if (errno /= AB7_NO_ERROR) return
467     ABI_ALLOCATE(kpt,(3, nkpt))
468     ABI_ALLOCATE(wkpt,(nkpt))
469     call kpoints_binding_auto_k_2(symid, nkpt, kpt, wkpt, kptrlatt, kptrlen_, &
470        & nshiftk, shiftk, errno)
471   end subroutine kpoints_get_auto_k_grid

m_ab7_kpoints/kpoints_get_irreductible_zone [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

  kpoints_get_irreductible_zone

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

 73   subroutine kpoints_get_irreductible_zone(irrzon, phnons, &
 74        & n1, n2, n3, nsppol, nspden, symid, errno)
 75 
 76 
 77 !This section has been created automatically by the script Abilint (TD).
 78 !Do not modify the following lines by hand.
 79 #undef ABI_FUNC
 80 #define ABI_FUNC 'kpoints_get_irreductible_zone'
 81 !End of the abilint section
 82 
 83     integer, intent(in)   :: symid
 84     integer, intent(in)   :: n1, n2, n3, nsppol, nspden
 85     integer, intent(out)  :: irrzon(n1*n2*n3,2,(nspden/nsppol)-3*(nspden/4))
 86     real(dp), intent(out) :: phnons(2,n1*n2*n3,(nspden/nsppol)-3*(nspden/4))
 87     integer, intent(out)  :: errno
 88 
 89     type(symmetry_type), pointer  :: sym
 90 
 91     if (AB_DBG) write(std_err,*) "AB kpoints: call get irreductible zone."
 92 
 93     errno = AB7_NO_ERROR
 94     call symmetry_get_from_id(sym, symid, errno)
 95     if (errno /= AB7_NO_ERROR) return
 96 
 97     if (sym%withSpin /= nspden) then
 98        errno = AB7_ERROR_ARG
 99        return
100     end if
101 
102     call irrzg(irrzon, nspden, nsppol, sym%nSym, n1, n2, n3, phnons, &
103          & sym%symAfm, sym%sym, sym%transNon)
104   end subroutine kpoints_get_irreductible_zone

m_ab7_kpoints/kpoints_get_mp_k_grid [ Functions ]

[ Top ] [ m_ab7_kpoints ] [ Functions ]

NAME

 kpoints_get_mp_k_grid

FUNCTION

INPUTS

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

CHILDREN

      kpoints_binding_auto_k_1,kpoints_binding_auto_k_2

SOURCE

251   subroutine kpoints_get_mp_k_grid(symid, nkpt, kpt, wkpt, &
252        & ngkpt, nshiftk, shiftk, errno)
253 
254 
255 !This section has been created automatically by the script Abilint (TD).
256 !Do not modify the following lines by hand.
257 #undef ABI_FUNC
258 #define ABI_FUNC 'kpoints_get_mp_k_grid'
259 !End of the abilint section
260 
261     integer, intent(in)  :: symid
262     integer, intent(out) :: errno
263     integer, intent(in) :: ngkpt(3)
264     integer, intent(in) :: nshiftk
265     real(dp), intent(in) :: shiftk(3, nshiftk)
266     integer, intent(out) :: nkpt
267     real(dp), pointer :: kpt(:,:), wkpt(:)
268 
269     real(dp) :: kptrlen
270     integer :: kptrlatt(3,3)
271     integer :: nshiftk_
272     real(dp) :: shiftk_(3,210)
273 
274     if (AB_DBG) write(std_err,*) "AB symmetry: call get k grid."
275 
276     nshiftk_ = nshiftk
277     shiftk_(:,1:nshiftk_) = shiftk(:,:)
278 
279     call kpoints_binding_mp_k_1(symid, nkpt, ngkpt, kptrlatt, kptrlen, &
280          & nshiftk_, shiftk_, errno)
281     if (errno /= AB7_NO_ERROR) return
282     ABI_ALLOCATE(kpt,(3, nkpt))
283     ABI_ALLOCATE(wkpt,(nkpt))
284     call kpoints_binding_mp_k_2(symid, nkpt, kpt, wkpt, &
285        & kptrlatt, kptrlen, nshiftk_, shiftk_, errno)
286   end subroutine kpoints_get_mp_k_grid