TABLE OF CONTENTS


ABINIT/m_orbmag [ Modules ]

[ Top ] [ Modules ]

NAME

  m_orbmag

FUNCTION

  This module contains the declaration of data types and methods
  used to handle orbital magnetization

COPYRIGHT

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

INPUTS

OUTPUT

NOTES

PARENTS

CHILDREN

SOURCE

28 #if defined HAVE_CONFIG_H
29 #include "config.h"
30 #endif
31 
32 #include "abi_common.h"
33 
34 module m_orbmag
35 
36  use defs_basis
37  use m_profiling_abi
38  use m_errors
39 
40  use m_pawcprj, only : pawcprj_type, pawcprj_free
41 
42  implicit none
43 
44  private

defs_datatypes/orbmag_type [ Types ]

[ Top ] [ defs_datatypes ] [ Types ]

NAME

 orbmag_type

FUNCTION

 variables used in orbital magnetism calculation

SOURCE

 57  type, public :: orbmag_type
 58 
 59 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
 60 ! declared in another part of ABINIT, that might need to take into account your modification.
 61 
 62 ! Integer variables
 63   integer :: orbmag              ! value of orbmag input variable in use
 64   integer :: fmkmem              ! number of k-points in the FBZ per cpu
 65   integer :: fmkmem_max          ! max of fmkmem
 66   integer :: fnkpt               ! number of k-points in the FBZ
 67   integer :: lmax
 68   integer :: lmnmax
 69   integer :: lmn2max
 70   integer :: mkmem_max           ! max of mkmem
 71   integer :: natom               ! number of atoms in unit cell
 72   integer :: my_natom            ! number of atoms treated by current proc
 73   integer :: mband_occ           ! max number of occupied bands (over spin)
 74                                  ! this number must be the same for every k
 75   integer :: nspinor             ! nspinor input from data set
 76   integer :: nsym
 77   integer :: usepaw              ! 1 if a PAW calculation, 0 else
 78 
 79 ! Real(dp) scalars
 80   real(dp) :: sdeg               ! spin degeneracy: sdeg = 2 if nsppol = 1
 81 
 82   ! Real(dp) arrays
 83   real(dp) :: chern(2,3)           ! result of chern number calculation
 84   
 85   real(dp) :: dkvecs(3,3)        ! dkvec(:,idir) = vector between a k-poinit
 86                                  ! and its nearest neighbour along idir
 87 ! Integer pointers
 88   integer, allocatable :: atom_indsym(:,:,:) ! atom_indsym(4,nsym,natom)
 89                                          ! this is data on how the symmetries map the atoms in the cell
 90                                          ! see symatm.F90 for full description
 91   integer, allocatable :: cgindex(:,:)    ! cgindex(nkpt,nsppol)
 92                                       ! for each k-point, stores the location
 93                                       ! of the WF in the cg array
 94   integer, allocatable :: cprjindex(:,:)  ! cprjindex(nkpt,nsppol)
 95                                       ! for each k-point, stores the location
 96                                       ! of the cprj in the cprj array (used only
 97                                       ! for PAW calculations)
 98   integer, allocatable :: fkgindex(:)     ! same as kgindex, but defined
 99                                       ! for the FBZ and intended to use
100                                       ! with pwindf
101   integer, allocatable :: ikpt_dk(:,:,:)  ! ikpt_dk(nkpt,2,3)
102                                       ! ikpt_dp(ikpt,ii,idir) = index of the
103                                       ! k-point at k+dk (ii=1) and k-dk (ii=2)
104   integer, allocatable :: indkk_f2ibz(:,:)   ! indkk_f2ibz(1:dtorbmag%fnkpt,1:6)
105                                          ! information needed to fold a
106                                          ! k-point in the FBZ into the IBZ;
107                                          ! the second index (1:6)
108                                          ! is as described in listkk
109   integer, allocatable :: i2fbz(:)           ! i2fbz(1:nkpt) gives index of IBZ
110                                          ! k-points in the FBZ k-point list
111 
112   integer, allocatable :: kgindex(:)      ! kgind(nkpt)
113                                       ! kgind(ikpt) = ikg
114 
115   integer, allocatable :: lmn_size(:)        ! lmn_size(ntypat)
116   integer, allocatable :: lmn2_size(:)       ! lmn2_size(ntypat)
117 
118   integer, allocatable :: nband_occ(:)       ! nband_occ(nsppol) = actual number of occupied bands
119                                              !  can be different for spin up and down!!!
120 ! Real(dp) allocatables
121 
122   real(dp), allocatable :: fkptns(:,:)       ! fkptns(3,1:dtorbmag%fnkpt)
123                                          ! k-points in FBZ
124 
125   real(dp), allocatable :: zarot(:,:,:,:)
126    !  zarot(l_size_max,l_size_max,l_max,nsym)
127    !  Coeffs of the transformation of real spherical
128    !  harmonics under the symmetry operations. These are needed when the
129    ! cprj's need to be computed in the full BZ, that is,
130    ! in the PAW case with kptopt /= 3.
131 
132  end type orbmag_type
133 
134  ! Bound methods:
135  public :: destroy_orbmag

m_orbmag/destroy_orbmag [ Functions ]

[ Top ] [ m_orbmag ] [ Functions ]

NAME

FUNCTION

   deallocate fields in orbmag structure

INPUTS

OUTPUT

SOURCE

152 subroutine destroy_orbmag(dtorbmag)
153 
154 
155 !This section has been created automatically by the script Abilint (TD).
156 !Do not modify the following lines by hand.
157 #undef ABI_FUNC
158 #define ABI_FUNC 'destroy_orbmag'
159 !End of the abilint section
160 
161  implicit none
162 
163 !Arguments ------------------------------------
164 !array
165  type(orbmag_type),intent(inout) :: dtorbmag
166 
167 ! ************************************************************************
168 
169 ! Integer pointers
170   if(allocated(dtorbmag%atom_indsym))  then
171     ABI_DEALLOCATE(dtorbmag%atom_indsym)
172   end if
173   if(allocated(dtorbmag%cgindex))  then
174     ABI_DEALLOCATE(dtorbmag%cgindex)
175   end if
176   if(allocated(dtorbmag%cprjindex))  then
177     ABI_DEALLOCATE(dtorbmag%cprjindex)
178   end if
179   if(allocated(dtorbmag%fkgindex))  then
180     ABI_DEALLOCATE(dtorbmag%fkgindex)
181   end if
182   if(allocated(dtorbmag%ikpt_dk))  then
183     ABI_DEALLOCATE(dtorbmag%ikpt_dk)
184   end if
185   if(allocated(dtorbmag%indkk_f2ibz))  then
186     ABI_DEALLOCATE(dtorbmag%indkk_f2ibz)
187   end if
188   if(allocated(dtorbmag%i2fbz))  then
189     ABI_DEALLOCATE(dtorbmag%i2fbz)
190   end if
191   if(allocated(dtorbmag%kgindex))  then
192     ABI_DEALLOCATE(dtorbmag%kgindex)
193   end if
194   if(allocated(dtorbmag%lmn_size))  then
195     ABI_DEALLOCATE(dtorbmag%lmn_size)
196   end if
197   if(allocated(dtorbmag%lmn2_size))  then
198     ABI_DEALLOCATE(dtorbmag%lmn2_size)
199   end if
200   if(allocated(dtorbmag%nband_occ))  then
201     ABI_DEALLOCATE(dtorbmag%nband_occ)
202   end if
203 ! Real(dp) pointers
204 
205   if(allocated(dtorbmag%fkptns))  then
206     ABI_DEALLOCATE(dtorbmag%fkptns)
207   end if
208   if(allocated(dtorbmag%zarot))  then
209     ABI_DEALLOCATE(dtorbmag%zarot)
210   end if
211 
212 end subroutine destroy_orbmag