TABLE OF CONTENTS


ABINIT/m_ab7_kpoints [ Modules ]

[ Top ] [ Modules ]

NAME

FUNCTION

COPYRIGHT

  Copyright (C) 2008-2022 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 .

SOURCE

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

SOURCE

252   subroutine kpoints_binding_auto_k_1(symid, nkpt, kptrlatt, kptrlen, &
253        & nshiftk, shiftk, errno)
254 
255     integer, intent(in)  :: symid
256     integer, intent(out) :: errno
257     integer, intent(out) :: nkpt
258     real(dp), intent(inout) :: kptrlen
259     integer, intent(out) :: nshiftk
260     real(dp), intent(out) :: shiftk(3,MAX_NSHIFTK)
261     integer, intent(out) :: kptrlatt(3,3)
262 
263     type(symmetry_type), pointer  :: sym
264     real(dp), allocatable :: kpt(:,:), wkpt(:)
265 
266     if (AB_DBG) write(std_err,*) "AB symmetry: call get auto k grid1."
267 
268     errno = AB7_NO_ERROR
269     call symmetry_get_from_id(sym, symid, errno)
270     if (errno /= AB7_NO_ERROR) return
271 
272     !  The parameters of the k lattice are not known, compute
273     !  kptrlatt, nshiftk, shiftk.
274     call testkgrid(sym%bravais,6,kptrlatt,kptrlen,&
275          & AB7_MAX_SYMMETRIES,nshiftk,sym%nSym,0,sym%rprimd,&
276          & shiftk,sym%symAfm,sym%sym,sym%vacuum)
277     if (AB_DBG) write(std_err,*) "AB symmetry: testkgrid -> kptrlatt=", kptrlatt
278 
279     nkpt=0
280     ABI_MALLOC(kpt,(3, nkpt))
281     ABI_MALLOC(wkpt,(nkpt))
282 
283     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
284          & AB7_MAX_SYMMETRIES, 0, nkpt, nshiftk, sym%nSym, &
285          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
286          & sym%vacuum, wkpt)
287     if (AB_DBG) write(std_err,*) "AB symmetry: getkgrid -> nkpt=", nkpt
288 
289     ABI_FREE(kpt)
290     ABI_FREE(wkpt)
291 
292   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

SOURCE

311   subroutine kpoints_binding_auto_k_2(symid, nkpt, kpt, wkpt, kptrlatt, kptrlen, &
312        & nshiftk, shiftk, errno)
313 
314     integer, intent(in)  :: symid
315     integer, intent(out) :: errno
316     integer, intent(in) :: nkpt
317     real(dp), intent(out) :: kpt(3,nkpt), wkpt(nkpt)
318     real(dp), intent(inout) :: kptrlen
319     integer, intent(inout) :: nshiftk
320     real(dp), intent(inout) :: shiftk(3,MAX_NSHIFTK)
321     integer, intent(inout) :: kptrlatt(3,3)
322 
323     type(symmetry_type), pointer  :: sym
324     integer :: nkpt_
325 
326     if (AB_DBG) write(std_err,*) "AB symmetry: call get auto k grid2."
327 
328     errno = AB7_NO_ERROR
329     call symmetry_get_from_id(sym, symid, errno)
330     if (errno /= AB7_NO_ERROR) return
331 
332     ! Then, we call it again to get the actual values for the k points.
333     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
334          & AB7_MAX_SYMMETRIES, nkpt, nkpt_, nshiftk, sym%nSym, &
335          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
336          & sym%vacuum, wkpt)
337   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

SOURCE

109   subroutine kpoints_binding_mp_k_1(symid, nkpt, ngkpt, &
110        & kptrlatt, kptrlen, nshiftk, shiftk, errno)
111 
112     integer, intent(in)  :: symid
113     integer, intent(out) :: errno
114     integer, intent(in) :: ngkpt(3)
115     integer, intent(inout) :: nshiftk
116     real(dp), intent(inout) :: shiftk(3,MAX_NSHIFTK)
117     real(dp), intent(out) :: kptrlen
118     integer, intent(out) :: kptrlatt(3,3)
119     integer, intent(out) :: nkpt
120 
121     type(symmetry_type), pointer  :: sym
122     real(dp), allocatable :: kpt(:,:), wkpt(:)
123 
124     if (AB_DBG) write(std_err,*) "AB symmetry: call get k grid1."
125 
126     errno = AB7_NO_ERROR
127     call symmetry_get_from_id(sym, symid, errno)
128     if (errno /= AB7_NO_ERROR) return
129 
130     ! First, compute the number of kpoints
131     kptrlatt(:,:) = 0
132     kptrlatt(1,1) = ngkpt(1)
133     kptrlatt(2,2) = ngkpt(2)
134     kptrlatt(3,3) = ngkpt(3)
135     kptrlen = 20.
136 
137     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
138          & AB7_MAX_SYMMETRIES, 0, nkpt, nshiftk, sym%nSym, &
139          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
140          & sym%vacuum, wkpt)
141   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

SOURCE

160   subroutine kpoints_binding_mp_k_2(symid, nkpt, kpt, wkpt, &
161        & kptrlatt, kptrlen, nshiftk, shiftk, errno)
162 
163     integer, intent(in)  :: symid
164     integer, intent(out) :: errno
165     integer, intent(inout) :: nshiftk
166     real(dp), intent(inout) :: shiftk(3,MAX_NSHIFTK)
167     integer, intent(in) :: nkpt
168     real(dp), intent(out) :: kpt(3,nkpt), wkpt(nkpt)
169     real(dp), intent(inout) :: kptrlen
170     integer, intent(inout) :: kptrlatt(3,3)
171 
172     type(symmetry_type), pointer  :: sym
173     integer :: nkpt_
174 
175     if (AB_DBG) write(std_err,*) "AB symmetry: call get k grid2."
176 
177     errno = AB7_NO_ERROR
178     call symmetry_get_from_id(sym, symid, errno)
179     if (errno /= AB7_NO_ERROR) return
180 
181     ! Then, we call it again to get the actual values for the k points.
182     call getkgrid(0, 0, 1, kpt, 1, kptrlatt, kptrlen, &
183          & AB7_MAX_SYMMETRIES, nkpt, nkpt_, nshiftk, sym%nSym, &
184          & sym%rprimd, shiftk, sym%symAfm, sym%sym, &
185          & sym%vacuum, wkpt)
186   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

SOURCE

356   subroutine kpoints_get_auto_k_grid(symid, nkpt, kpt, wkpt, &
357        & kptrlen, errno)
358 
359     integer, intent(in)  :: symid
360     integer, intent(out) :: errno
361     integer, intent(out) :: nkpt
362     real(dp), intent(in) :: kptrlen
363     real(dp), pointer :: kpt(:,:), wkpt(:)
364 
365     real(dp) :: kptrlen_
366     integer :: kptrlatt(3,3)
367     integer :: nshiftk
368     real(dp) :: shiftk(3,MAX_NSHIFTK)
369 
370     if (AB_DBG) write(std_err,*) "AB symmetry: call get auto k grid."
371 
372     kptrlen_ = kptrlen
373     call kpoints_binding_auto_k_1(symid, nkpt, kptrlatt, kptrlen_, &
374        & nshiftk, shiftk, errno)
375     if (errno /= AB7_NO_ERROR) return
376     ABI_MALLOC(kpt,(3, nkpt))
377     ABI_MALLOC(wkpt,(nkpt))
378     call kpoints_binding_auto_k_2(symid, nkpt, kpt, wkpt, kptrlatt, kptrlen_, &
379        & nshiftk, shiftk, errno)
380   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

SOURCE

65   subroutine kpoints_get_irreductible_zone(irrzon, phnons, &
66        & n1, n2, n3, nsppol, nspden, symid, errno)
67 
68     integer, intent(in)   :: symid
69     integer, intent(in)   :: n1, n2, n3, nsppol, nspden
70     integer, intent(out)  :: irrzon(n1*n2*n3,2,(nspden/nsppol)-3*(nspden/4))
71     real(dp), intent(out) :: phnons(2,n1*n2*n3,(nspden/nsppol)-3*(nspden/4))
72     integer, intent(out)  :: errno
73 
74     type(symmetry_type), pointer  :: sym
75 
76     if (AB_DBG) write(std_err,*) "AB kpoints: call get irreductible zone."
77 
78     errno = AB7_NO_ERROR
79     call symmetry_get_from_id(sym, symid, errno)
80     if (errno /= AB7_NO_ERROR) return
81 
82     if (sym%withSpin /= nspden) then
83        errno = AB7_ERROR_ARG
84        return
85     end if
86 
87     call irrzg(irrzon, nspden, nsppol, sym%nSym, n1, n2, n3, phnons, &
88          & sym%symAfm, sym%sym, sym%transNon)
89   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

SOURCE

205   subroutine kpoints_get_mp_k_grid(symid, nkpt, kpt, wkpt, &
206        & ngkpt, nshiftk, shiftk, errno)
207 
208     integer, intent(in)  :: symid
209     integer, intent(out) :: errno
210     integer, intent(in) :: ngkpt(3)
211     integer, intent(in) :: nshiftk
212     real(dp), intent(in) :: shiftk(3, nshiftk)
213     integer, intent(out) :: nkpt
214     real(dp), pointer :: kpt(:,:), wkpt(:)
215 
216     real(dp) :: kptrlen
217     integer :: kptrlatt(3,3)
218     integer :: nshiftk_
219     real(dp) :: shiftk_(3,MAX_NSHIFTK)
220 
221     if (AB_DBG) write(std_err,*) "AB symmetry: call get k grid."
222 
223     nshiftk_ = nshiftk
224     shiftk_(:,1:nshiftk_) = shiftk(:,:)
225 
226     call kpoints_binding_mp_k_1(symid, nkpt, ngkpt, kptrlatt, kptrlen, &
227          & nshiftk_, shiftk_, errno)
228     if (errno /= AB7_NO_ERROR) return
229     ABI_MALLOC(kpt,(3, nkpt))
230     ABI_MALLOC(wkpt,(nkpt))
231     call kpoints_binding_mp_k_2(symid, nkpt, kpt, wkpt, &
232        & kptrlatt, kptrlen, nshiftk_, shiftk_, errno)
233   end subroutine kpoints_get_mp_k_grid