TABLE OF CONTENTS


ABINIT/m_lgroup [ Modules ]

[ Top ] [ Modules ]

NAME

 m_lgroup

FUNCTION

 The little group of a q-point is defined as the subset of the space group that preserves q,
 modulo a G0 vector (also called umklapp vector). Namely:

    Sq = q + G0

 where S is an operation in reciprocal space (symrec)
 If time reversal symmetry can be used, it is possible to enlarge the little group by
 including the operations such as:

   -Sq = q + G0

 The operations of little group define an irriducible wedge, IBZ(q), that is usually larger
 than the irredubile zone defined by the point group of the crystal. The two zones coincide when q=0

COPYRIGHT

  Copyright (C) 2008-2018 ABINIT group (MG)
  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

CHILDREN

SOURCE

33 #if defined HAVE_CONFIG_H
34 #include "config.h"
35 #endif
36 
37 #include "abi_common.h"
38 
39 module m_lgroup
40 
41  use defs_basis
42  use m_errors
43  use m_abicore
44  use m_crystal
45  use m_symkpt
46 
47  use m_fstrings,   only : ftoa, ktoa, sjoin
48  use m_symtk,      only : chkgrp, littlegroup_q
49 
50  implicit none
51 
52  private

m_lgroup/lgroup_print [ Functions ]

[ Top ] [ m_lgroup ] [ Functions ]

NAME

 lgroup_print

FUNCTION

  Print the object

INPUTS

  [title]=String to be printed as header for additional info.
  [unit]=Unit number for output
  [prtvol]=Verbosity level

OUTPUT

  Only printing

PARENTS

CHILDREN

SOURCE

236 subroutine lgroup_print(self, title, unit, prtvol)
237 
238 
239 !This section has been created automatically by the script Abilint (TD).
240 !Do not modify the following lines by hand.
241 #undef ABI_FUNC
242 #define ABI_FUNC 'lgroup_print'
243 !End of the abilint section
244 
245  implicit none
246 
247 !Arguments ------------------------------------
248  integer,optional,intent(in) :: unit, prtvol
249  character(len=*),optional,intent(in) :: title
250  type(lgroup_t),intent(in) :: self
251 
252 !Local variables-------------------------------
253 !scalars
254  integer :: my_prtvol,my_unt,ik
255  character(len=500) :: msg
256 ! *************************************************************************
257 
258  my_unt = std_out; if (present(unit)) my_unt = unit
259  my_prtvol = 0; if (present(prtvol)) my_prtvol = prtvol
260 
261  msg = ' ==== Info on the <lgroup_t> object ==== '
262  if (present(title)) msg = ' ==== '//trim(adjustl(title))//' ==== '
263  call wrtout(my_unt, msg)
264 
265  write(msg, '(3a, 2(a, i0), a)') &
266   ' Little group point: ................... ', ktoa(self%point), ch10, &
267   ' Number of points in IBZ(p) ............ ', self%nibz, ch10, &
268   ' Time-reversal flag (0: No, 1: Yes) .... ', self%timrev
269  call wrtout(my_unt, msg)
270 
271  if (my_prtvol /= 0) then
272    do ik=1,self%nibz
273       call wrtout(my_unt, sjoin(ktoa(self%ibz(:,ik)), ftoa(self%weights(ik))))
274    end do
275  end if
276 
277 end subroutine lgroup_print

m_sigmaph/lgroup_free [ Functions ]

[ Top ] [ m_sigmaph ] [ Functions ]

NAME

  lgroup_free

FUNCTION

  Free memory

PARENTS

      m_sigmaph

CHILDREN

SOURCE

294 subroutine lgroup_free(self)
295 
296 
297 !This section has been created automatically by the script Abilint (TD).
298 !Do not modify the following lines by hand.
299 #undef ABI_FUNC
300 #define ABI_FUNC 'lgroup_free'
301 !End of the abilint section
302 
303  implicit none
304 
305 !Arguments ------------------------------------
306  type(lgroup_t),intent(inout) :: self
307 
308 ! *************************************************************************
309 
310  ! integer
311  if (allocated(self%symtab)) then
312    ABI_FREE(self%symtab)
313  end if
314  ! real
315  if (allocated(self%ibz)) then
316    ABI_FREE(self%ibz)
317  end if
318  if (allocated(self%weights)) then
319    ABI_FREE(self%weights)
320  end if
321 
322 end subroutine lgroup_free

m_sigmaph/lgroup_new [ Functions ]

[ Top ] [ m_sigmaph ] [ Functions ]

NAME

  lgroup_new

FUNCTION

  Build the little group of the k-point.

INPUTS

  cryst(crystal_t)=Crystalline structure
  kpoint(3)=External k-point defining the little-group
  timrev=1 if time-reversal symmetry can be used, 0 otherwise.
  nkbz=Number of k-points in the BZ.
  kbz(3,nkbz)=K-points in the BZ.
  nkibz=Number of k-points in the IBZ
  kibz(3,nkibz)=Irreducible zone.

PARENTS

CHILDREN

SOURCE

122 type (lgroup_t) function lgroup_new(cryst, kpoint, timrev, nkbz, kbz, nkibz, kibz) result(new)
123 
124 
125 !This section has been created automatically by the script Abilint (TD).
126 !Do not modify the following lines by hand.
127 #undef ABI_FUNC
128 #define ABI_FUNC 'lgroup_new'
129 !End of the abilint section
130 
131  implicit none
132 
133 !Arguments ------------------------------------
134 !scalars
135  integer,intent(in) :: timrev,nkibz,nkbz
136  type(crystal_t),intent(in) :: cryst
137 !arrays
138  real(dp),intent(in) :: kpoint(3),kbz(3,nkbz),kibz(3,nkibz)
139 
140 !Local variables ------------------------------
141 !scalars
142  integer,parameter :: iout0=0,my_timrev0=0,chksymbreak0=0,debug=0
143  integer :: otimrev_k,ierr,itim,isym,nsym_lg,ik
144 !arrays
145  integer :: symrec_lg(3,3,2*cryst%nsym),symafm_lg(3,3,2*cryst%nsym)
146  integer,allocatable :: ibz2bz(:)
147  real(dp),allocatable :: wtk(:),wtk_folded(:)
148 
149 ! *************************************************************************
150 
151  ! TODO: Option to exclude umklapp/time-reversal symmetry and kptopt
152  new%point = kpoint
153  new%timrev = timrev
154 
155  ! Determines the symmetry operations by which the k-point is preserved,
156  ABI_MALLOC(new%symtab, (4, 2, cryst%nsym))
157  call littlegroup_q(cryst%nsym, kpoint, new%symtab, cryst%symrec, cryst%symafm, otimrev_k, prtvol=0)
158 
159  ABI_CHECK(new%timrev==1, "timrev == 0 not coded")
160  nsym_lg = 0
161  do itim=1,2
162    do isym=1,cryst%nsym
163      if (cryst%symafm(isym) == -1) cycle
164      if (new%symtab(4, itim, isym) /= 1) cycle ! not \pm Sq = q+g0
165      nsym_lg = nsym_lg + 1
166      symrec_lg(:,:,nsym_lg) = cryst%symrec(:,:,isym) * (-2* itim + 3)
167    end do
168  end do
169 
170  ! Check group closure.
171  symafm_lg = 1
172  call chkgrp(nsym_lg, symafm_lg, symrec_lg, ierr)
173  ABI_CHECK(ierr == 0, "Error in group closure")
174 
175  ! Find the irreducible zone with the little group operations.
176  ! Do not use time-reversal since it has been manually introduced previously
177  ABI_MALLOC(ibz2bz, (nkbz))
178  ABI_MALLOC(wtk_folded, (nkbz))
179  ABI_MALLOC(wtk, (nkbz))
180  wtk = one / nkbz ! Weights sum up to one
181 
182  ! TODO: In principle here we would like to have a set that contains the initial IBZ.
183  call symkpt(chksymbreak0,cryst%gmet,ibz2bz,iout0,kbz,nkbz,new%nibz,&
184    nsym_lg,symrec_lg,my_timrev0,wtk,wtk_folded)
185 
186  ABI_MALLOC(new%ibz, (3, new%nibz))
187  ABI_MALLOC(new%weights, (new%nibz))
188 
189  do ik=1,new%nibz
190    new%weights(ik) = wtk_folded(ibz2bz(ik))
191    new%ibz(:,ik) = kbz(:, ibz2bz(ik))
192  end do
193  ABI_CHECK(sum(new%weights) - one < tol12, sjoin("Weights don't sum up to one but to:", ftoa(sum(new%weights))))
194 
195  ABI_FREE(ibz2bz)
196  ABI_FREE(wtk_folded)
197  ABI_FREE(wtk)
198 
199  ! Debug section.
200  if (debug /= 0) then
201    do ik=1,new%nibz
202      if (ik <= nkibz) then
203        write(std_out,"(a)")sjoin(ktoa(new%ibz(:,ik)), ktoa(new%ibz(:,ik) - kibz(:,ik)), ftoa(new%weights(ik)))
204      else
205        write(std_out,"(a)")sjoin(ktoa(new%ibz(:,ik)), "[---]", ftoa(new%weights(ik)))
206      end if
207    end do
208  end if
209 
210 end function lgroup_new

m_sigmaph/lgroup_t [ Types ]

[ Top ] [ m_sigmaph ] [ Types ]

NAME

 lgroup_t

FUNCTION

  Stores tables associated to the little-group.

SOURCE

64  type, public :: lgroup_t
65 
66    integer :: nibz
67    ! Number of points in the IBZ(q)
68 
69    integer :: timrev
70    ! timrev=1 if time-reversal symmetry can be used, 0 otherwise.
71 
72    real(dp) :: point(3)
73    ! The external q-point.
74 
75    integer,allocatable :: symtab(:,:,:)
76    ! symtab(4, 2, cryst%nsym)
77    ! nsym is the **total** number of symmetries of the system as given by cryst%nsym
78    ! three first numbers define the G vector;
79    ! fourth number is zero if the q-vector is not preserved, is 1 otherwise.
80    ! second index is one without time-reversal symmetry, two with time-reversal symmetry
81 
82    real(dp),allocatable :: ibz(:,:)
83    ! ibz(3, nibz)
84    ! K-points in the IBZ(q)
85 
86    real(dp),allocatable :: weights(:)
87    ! weights(nibz)
88    ! Weights in the IBZ(q), normalized to 1
89 
90  end type lgroup_t