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

COPYRIGHT

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

SOURCE

23 #if defined HAVE_CONFIG_H
24 #include "config.h"
25 #endif
26 
27 #include "abi_common.h"
28 
29 module m_efield
30 
31  use defs_basis
32  use m_abicore
33  use m_errors
34 
35  use m_pawcprj, only : pawcprj_type, pawcprj_free
36 
37  implicit none
38 
39  private

m_efield/destroy_efield [ Functions ]

[ Top ] [ m_efield ] [ Functions ]

NAME

FUNCTION

   deallocate fields in efield structure

INPUTS

OUTPUT

SOURCE

246 subroutine destroy_efield(dtefield)
247 
248 !Arguments ------------------------------------
249 !array
250  type(efield_type),intent(inout) :: dtefield !vz_i
251 
252 ! ************************************************************************
253 
254 ! Integer pointers
255   if(allocated(dtefield%atom_indsym))  then
256     ABI_FREE(dtefield%atom_indsym)
257   end if
258   if(allocated(dtefield%cgindex))  then
259     ABI_FREE(dtefield%cgindex)
260   end if
261   if(allocated(dtefield%cgqindex))  then
262     ABI_FREE(dtefield%cgqindex)
263   end if
264   if(allocated(dtefield%cprjindex))  then
265     ABI_FREE(dtefield%cprjindex)
266   end if
267   if(allocated(dtefield%fkgindex))  then
268     ABI_FREE(dtefield%fkgindex)
269   end if
270   if(allocated(dtefield%idxkstr))  then
271     ABI_FREE(dtefield%idxkstr)
272   end if
273   if(allocated(dtefield%ikpt_dk))  then
274     ABI_FREE(dtefield%ikpt_dk)
275   end if
276   if(allocated(dtefield%indkk_f2ibz))  then
277     ABI_FREE(dtefield%indkk_f2ibz)
278   end if
279   if(allocated(dtefield%i2fbz))  then
280     ABI_FREE(dtefield%i2fbz)
281   end if
282   if(allocated(dtefield%kgindex))  then
283     ABI_FREE(dtefield%kgindex)
284   end if
285   if(allocated(dtefield%lmn_size))  then
286     ABI_FREE(dtefield%lmn_size)
287   end if
288   if(allocated(dtefield%lmn2_size))  then
289     ABI_FREE(dtefield%lmn2_size)
290   end if
291   if(allocated(dtefield%nband_occ))  then
292     ABI_FREE(dtefield%nband_occ)
293   end if
294   if(allocated(dtefield%nneigh))  then
295     ABI_FREE(dtefield%nneigh)
296   end if
297   if(allocated(dtefield%sflag))  then
298     ABI_FREE(dtefield%sflag)
299   end if
300   if(allocated(dtefield%str_neigh))  then
301     ABI_FREE(dtefield%str_neigh)
302   end if
303   if(allocated(dtefield%strg_neigh))  then
304     ABI_FREE(dtefield%strg_neigh)
305   end if
306 
307 ! Real(dp) pointers
308 
309   if(allocated(dtefield%coord_str))  then
310     ABI_FREE(dtefield%coord_str)
311   end if
312   if(allocated(dtefield%epawf3))  then
313     ABI_FREE(dtefield%epawf3)
314   end if
315   if(allocated(dtefield%epaws3))  then
316     ABI_FREE(dtefield%epaws3)
317   end if
318   if(allocated(dtefield%expibi))  then
319     ABI_FREE(dtefield%expibi)
320   end if
321   if(allocated(dtefield%fkptns))  then
322     ABI_FREE(dtefield%fkptns)
323   end if
324   if(allocated(dtefield%qijb_kk))  then
325     ABI_FREE(dtefield%qijb_kk)
326   end if
327   if(allocated(dtefield%rij))  then
328     ABI_FREE(dtefield%rij)
329   end if
330   if(allocated(dtefield%smat))  then
331     ABI_FREE(dtefield%smat)
332   end if
333   if(allocated(dtefield%zarot))  then
334     ABI_FREE(dtefield%zarot)
335   end if
336 
337 ! pointer to cprj
338   if(allocated(dtefield%cprj)) then
339     call pawcprj_free(dtefield%cprj)
340     ABI_FREE(dtefield%cprj)
341   end if
342 
343 end subroutine destroy_efield

m_efield/efield_type [ Types ]

[ Top ] [ m_efield ] [ Types ]

NAME

 efield_type

FUNCTION

 First-principles calculations in a finite electric field

SOURCE

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