TABLE OF CONTENTS


ABINIT/m_efield [ Modules ]

[ Top ] [ Modules ]

NAME

  m_efield

FUNCTION

  This module contains the declaration of data types and methods
  used to handle electric fields
  Imported object from defs_datatypes

COPYRIGHT

 Copyright (C) 2011-2018 ABINIT group (MJV)
 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

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

defs_datatypes/efield_type [ Types ]

[ Top ] [ defs_datatypes ] [ Types ]

NAME

 efield_type

FUNCTION

 First-principles calculations in a finite electric field

SOURCE

 58  type, public :: efield_type
 59 
 60 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
 61 ! declared in another part of ABINIT, that might need to take into account your modification.
 62 
 63 ! Integer variables
 64   integer :: berryopt            ! value of berryopt in use
 65   integer :: fmkmem              ! number of k-points in the FBZ per cpu
 66   integer :: fmkmem_max          ! max of fmkmem
 67   integer :: fnkpt               ! number of k-points in the FBZ
 68   integer :: has_epawf3          ! 2 if epawf3 computed, 1 if only allocated, zero else
 69   integer :: has_epaws3          ! 2 if epaws3 computed, 1 if only allocated, zero else
 70   integer :: has_expibi          ! 2 if expibi computed, 1 if only allocated, zero else
 71   integer :: has_rij             ! 2 if paw_rij computed, 1 if only allocated, zero else
 72   integer :: has_qijb            ! 2 if paw_qijb computed, 1 if only allocated, zero else
 73   integer :: lmax
 74   integer :: lmnmax
 75   integer :: lmn2max
 76   integer :: maxnstr             ! max number of strings along idir=1,2,3
 77   integer :: maxnkstr            ! max number of k-points per string
 78   integer :: mkmem_max           ! max of mkmem
 79   integer :: natom               ! number of atoms in unit cell
 80   integer :: my_natom            ! number of atoms treated by current proc
 81   integer :: mband_occ           ! max number of occupied bands (over spin)
 82                                  ! this number must be the same for every k
 83   integer :: nspinor             ! nspinor input from data set
 84   integer :: nsym
 85   integer :: usecprj             ! 1 if efield%cprj allocated (see below), 0 else
 86   integer :: usepaw              ! 1 if a PAW calculation, 0 else
 87 
 88 ! Integer arrays
 89 
 90   integer :: nstr(3)             ! nstr(idir) = number of strings along idir
 91   integer :: nkstr(3)            ! nkstr(idir) = number of k-points per string
 92 
 93 ! Real(dp) scalars
 94   real(dp) :: sdeg               ! spin degeneracy: sdeg = 2 if nsppol = 1
 95                                  !                         1 if nsppol = 2
 96 
 97 ! Real(dp) arrays
 98   real(dp) :: dkvecs(3,3)        ! dkvec(:,idir) = vector between a k-poinit
 99                                  ! and its nearest neighbour along idir
100 !! if berryopt=4,6,7, dtset%red_efieldbar and dtset%red_dfield will be initialized in 67_common/initberry.F90
101   real(dp) :: efield_dot(3)      ! reciprocal lattice coordinates of the
102                                  ! electric field
103   real(dp) :: red_ptot1(3)       ! reduced total polarization
104   real(dp) :: efield2(3)         ! unreduced electric field, only used when berryopt == 14 in order to save real electric field for print out.
105 
106   real(dp) :: gmet_str(2,2,3)    ! gmet_str(:,:,idir) is the metric of the metric of
107                                  ! the space of strings of direction idir
108 ! Integer pointers
109   integer, allocatable :: atom_indsym(:,:,:) ! atom_indsym(4,nsym,natom)
110                                          ! this is data on how the symmetries map the atoms in the cell
111                                          ! see symatm.F90 for full description
112   integer, allocatable :: cgindex(:,:)    ! cgindex(nkpt,nsppol)
113                                       ! for each k-point, stores the location
114                                       ! of the WF in the cg array
115   integer, allocatable :: cgqindex(:,:,:) ! cgqindex(3,6,nkpt*nsppol)
116                                       ! for each k-point, stores the location
117                                       ! of the WF in the cgq and pwnsfacq
118                                       ! arrays
119                                       ! (see vtorho.f and initberry.f)
120   integer, allocatable :: cprjindex(:,:)  ! cprjindex(nkpt,nsppol)
121                                       ! for each k-point, stores the location
122                                       ! of the cprj in the cprj array (used only
123                                       ! for PAW calculations)
124   integer, allocatable :: fkgindex(:)     ! same as kgindex, but defined
125                                       ! for the FBZ and intended to use
126                                       ! with pwindf
127   integer, allocatable :: idxkstr(:,:,:)  ! idxkstr(maxnkstr,maxnstr,3)
128                                       ! idxkstr(ikstr,istr,idir) index (ikpt) of
129                                       ! k-point ikstr on string istr along idir
130   integer, allocatable :: ikpt_dk(:,:,:)  ! ikpt_dk(nkpt,2,3)
131                                       ! ikpt_dp(ikpt,ii,idir) = index of the
132                                       ! k-point at k+dk (ii=1) and k-dk (ii=2)
133   integer, allocatable :: indkk_f2ibz(:,:)   ! indkk_f2ibz(1:dtefield%fnkpt,1:6)
134                                          ! information needed to fold a
135                                          ! k-point in the FBZ into the IBZ;
136                                          ! the second index (1:6)
137                                          ! is as described in listkk
138   integer, allocatable :: i2fbz(:)           ! i2fbz(1:nkpt) gives index of IBZ
139                                          ! k-points in the FBZ k-point list
140 
141   integer, allocatable :: kgindex(:)      ! kgind(nkpt)
142                                       ! kgind(ikpt) = ikg
143 
144   integer, allocatable :: lmn_size(:)        ! lmn_size(ntypat)
145   integer, allocatable :: lmn2_size(:)       ! lmn2_size(ntypat)
146 
147   integer, allocatable :: nband_occ(:)       ! nband_occ(nsppol) = actual number of occupied bands
148                                              !  can be different for spin up and down!!!
149   integer, allocatable :: nneigh(:)          ! nneigh(nkpt)
150                                          ! for each k-point, nneigh stores
151                                          ! the number of its nearest neighbours
152                                          ! that are not related by symmetry
153   integer, allocatable :: sflag(:,:,:,:)  ! sflag(mband_occ,nkpt*nsppol,2,3)
154                                       ! sflag = 0 : compute the whole row of
155                                       !             smat
156                                       ! sflag = 1 : the row is up to date
157 
158   integer, allocatable :: str_neigh(:,:,:)
159   integer, allocatable :: strg_neigh(:,:,:,:)
160 ! str_neigh(ineigh, istr, idir) is the index ineigh-th neighbour of the istr-th string in
161 ! the direction idir
162 ! str_neigh(ineigh, istr, :, idir) is a 2-dimensional vector which coordinates are 0 or 1,
163 ! useful only if the k-point mesh isn't a full mesh - if it's a single point, a line or a plane.
164 
165 
166 ! Real(dp) allocatables
167 
168 ! the coordinates of the ineigh-th neighbour of the istr-th string in the direction idir are :
169 ! coord_str(:,str_neigh(ineigh,istr,idir),idir) + real(str_neigh(ineigh, istr, :, idir),dp)
170   real(dp),allocatable :: coord_str(:,:,:)
171 ! coord_str(1:2,istr,idir) are the coordinate of the istr-th string in the direction idir.
172 
173   real(dp),allocatable :: epawf3(:,:,:)
174 ! epawf3(natom,3,3) ! F3-type force term (derivatives of projectors with respect to ion posiion)
175 ! that arises in force for finite electric field with PAW
176 ! epawf3(iatom,idir,fdir) is derivative of polarization component idir with respect to iatom
177 ! displaced in direction fdir
178 ! see equation 32 of Torrent et al. CMS 42, 337 (2008)
179 
180   real(dp),allocatable :: epaws3(:,:,:)
181 ! epaws3(natom,3,6) ! F3-type stress term (derivatives of projectors with respect to strain)
182 ! that arises in stress for finite electric field with PAW
183 ! epaws3(iatom,idir,strain) is derivative of polarization component idir with respect to strain
184 ! component for atom iatom (note that these are on-site terms)
185 ! see equation D.7 of Torrent et al. CMS 42, 337 (2008)
186 
187   real(dp), allocatable :: expibi(:,:,:)
188 ! expibi(2,my_natom,3)
189 ! used for PAW field calculations (distributed over atomic sites)
190 ! stores the on-site phase factors arising from
191 ! $\langle\phi_{i,k}|\phi_{j,k+\sigma_k k_k}\rangle$
192 ! where $\sigma = \pm 1$. These overlaps arise in various Berry
193 ! phase calculations of electric and magnetic polarization. The on-site
194 ! phase factor is $\exp[-i\sigma_k k_k)\cdot I]$ where
195 ! $I$ is the nuclear position. Only the following
196 ! are computed and saved, in the given order:
197 ! 1)    -k_1
198 ! 2)    -k_2
199 ! 3)    -k_3
200 
201   real(dp), allocatable :: fkptns(:,:)       ! fkptns(3,1:dtefield%fnkpt)
202                                          ! k-points in FBZ
203 
204   real(dp), allocatable :: qijb_kk(:,:,:,:)
205 ! qijb_kk(2,lmn2max,natom,3)
206 ! on-site part of <u_nk|u_mk+b> matrix elements, relevant for PAW only
207 ! vector b described by idir (1,2,3), forward direction; value for
208 ! reverse direction (ifor = 2 in berryphase_new and cgwf) obtained by
209 ! complex conjugation
210 
211 ! pointer to on-site dipole moment
212   real(dp),allocatable :: rij(:,:,:) ! rij(lmn2_size_max,ntypat,3)
213  ! gives <r-R> at each atom in each of 3 directions
214  ! these are used only in the PAW case with electric field
215 
216   real(dp), allocatable :: smat(:,:,:,:,:,:)
217 ! smat(2,mband_occ,mband_occ,nkpt*nsppol,2,3)
218 ! Overlap matrix for every k-point. In an electric field calculation,
219 ! smat is updated at every iteration.
220 
221   real(dp), allocatable :: zarot(:,:,:,:)
222    !  zarot(l_size_max,l_size_max,l_max,nsym)
223    !  Coeffs of the transformation of real spherical
224    !  harmonics under the symmetry operations. These are needed when the
225    ! cprj's need to be computed in the full BZ, that is,
226    ! in the PAW case with kptopt /= 3.
227 
228 ! pointer to cprj
229    type(pawcprj_type),allocatable :: cprj(:,:)
230 ! used with finite efield and PAW
231 
232  end type efield_type
233 
234  ! Bound methods:
235  public :: destroy_efield

m_efield/destroy_efield [ Functions ]

[ Top ] [ m_efield ] [ Functions ]

NAME

FUNCTION

   deallocate fields in efield structure

INPUTS

OUTPUT

SOURCE

252 subroutine destroy_efield(dtefield)
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 'destroy_efield'
259 !End of the abilint section
260 
261  implicit none
262 
263 !Arguments ------------------------------------
264 !array
265  type(efield_type),intent(inout) :: dtefield !vz_i
266 
267 ! ************************************************************************
268 
269 ! Integer pointers
270   if(allocated(dtefield%atom_indsym))  then
271     ABI_DEALLOCATE(dtefield%atom_indsym)
272   end if
273   if(allocated(dtefield%cgindex))  then
274     ABI_DEALLOCATE(dtefield%cgindex)
275   end if
276   if(allocated(dtefield%cgqindex))  then
277     ABI_DEALLOCATE(dtefield%cgqindex)
278   end if
279   if(allocated(dtefield%cprjindex))  then
280     ABI_DEALLOCATE(dtefield%cprjindex)
281   end if
282   if(allocated(dtefield%fkgindex))  then
283     ABI_DEALLOCATE(dtefield%fkgindex)
284   end if
285   if(allocated(dtefield%idxkstr))  then
286     ABI_DEALLOCATE(dtefield%idxkstr)
287   end if
288   if(allocated(dtefield%ikpt_dk))  then
289     ABI_DEALLOCATE(dtefield%ikpt_dk)
290   end if
291   if(allocated(dtefield%indkk_f2ibz))  then
292     ABI_DEALLOCATE(dtefield%indkk_f2ibz)
293   end if
294   if(allocated(dtefield%i2fbz))  then
295     ABI_DEALLOCATE(dtefield%i2fbz)
296   end if
297   if(allocated(dtefield%kgindex))  then
298     ABI_DEALLOCATE(dtefield%kgindex)
299   end if
300   if(allocated(dtefield%lmn_size))  then
301     ABI_DEALLOCATE(dtefield%lmn_size)
302   end if
303   if(allocated(dtefield%lmn2_size))  then
304     ABI_DEALLOCATE(dtefield%lmn2_size)
305   end if
306   if(allocated(dtefield%nband_occ))  then
307     ABI_DEALLOCATE(dtefield%nband_occ)
308   end if
309   if(allocated(dtefield%nneigh))  then
310     ABI_DEALLOCATE(dtefield%nneigh)
311   end if
312   if(allocated(dtefield%sflag))  then
313     ABI_DEALLOCATE(dtefield%sflag)
314   end if
315   if(allocated(dtefield%str_neigh))  then
316     ABI_DEALLOCATE(dtefield%str_neigh)
317   end if
318   if(allocated(dtefield%strg_neigh))  then
319     ABI_DEALLOCATE(dtefield%strg_neigh)
320   end if
321 
322 ! Real(dp) pointers
323 
324   if(allocated(dtefield%coord_str))  then
325     ABI_DEALLOCATE(dtefield%coord_str)
326   end if
327   if(allocated(dtefield%epawf3))  then
328     ABI_DEALLOCATE(dtefield%epawf3)
329   end if
330   if(allocated(dtefield%epaws3))  then
331     ABI_DEALLOCATE(dtefield%epaws3)
332   end if
333   if(allocated(dtefield%expibi))  then
334     ABI_DEALLOCATE(dtefield%expibi)
335   end if
336   if(allocated(dtefield%fkptns))  then
337     ABI_DEALLOCATE(dtefield%fkptns)
338   end if
339   if(allocated(dtefield%qijb_kk))  then
340     ABI_DEALLOCATE(dtefield%qijb_kk)
341   end if
342   if(allocated(dtefield%rij))  then
343     ABI_DEALLOCATE(dtefield%rij)
344   end if
345   if(allocated(dtefield%smat))  then
346     ABI_DEALLOCATE(dtefield%smat)
347   end if
348   if(allocated(dtefield%zarot))  then
349     ABI_DEALLOCATE(dtefield%zarot)
350   end if
351 
352 ! pointer to cprj
353   if(allocated(dtefield%cprj)) then
354     call pawcprj_free(dtefield%cprj)
355     ABI_DATATYPE_DEALLOCATE(dtefield%cprj)
356   end if
357 
358 end subroutine destroy_efield