1! 2! Copyright (C) 2012 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8!--------------------------------------------------------------------------------- 9MODULE bp 10!--------------------------------------------------------------------------------- 11 !! Contains the variables needed for the Berry phase polarization calculation 12 ! 13 USE kinds, ONLY: DP 14 USE becmod, ONLY: bec_type 15 ! 16 SAVE 17 PRIVATE 18 PUBLIC :: lberry, lelfield, lorbm, gdir, nppstr, nberrycyc, evcel, evcelp, evcelm, & 19 fact_hepsi, bec_evcel, mapgp_global, mapgm_global, nppstr_3d, ion_pol, & 20 el_pol, fc_pol, l_el_pol_old, el_pol_old, el_pol_acc, nx_el, l3dstring, & 21 efield, efield_cart, efield_cry, transform_el, mapg_owner, phase_control 22 PUBLIC :: allocate_bp_efield, deallocate_bp_efield, bp_global_map 23 PUBLIC :: pdl_tot 24 ! 25 ! 26 LOGICAL :: lberry = .FALSE. 27 !! if .TRUE. calculate polarization using Berry phase 28 LOGICAL :: lelfield = .FALSE. 29 !! if .TRUE. finite electric field using Berry phase 30 LOGICAL :: lorbm = .FALSE. 31 !! if .TRUE. calculate orbital magnetization (Kubo terms) 32 INTEGER :: gdir 33 !! G-vector for polarization calculation 34 INTEGER :: nppstr 35 !! number of k-points (parallel vector) 36 INTEGER :: nberrycyc 37 !! number of cycles for convergence in electric field 38 !! without changing the selfconsistent charge 39 REAL(DP) :: efield 40 !! electric field intensity in a.u. 41 COMPLEX(DP), ALLOCATABLE , TARGET :: evcel(:,:) 42 !! wavefunctions for calculating the electric field operator 43 COMPLEX(DP), ALLOCATABLE , TARGET :: evcelm(:,:,:) 44 !! wavefunctions for storing projectors for electric field operator 45 COMPLEX(DP), ALLOCATABLE , TARGET :: evcelp(:,:,:) 46 !! wavefunctions for storing projectors for electric field operator 47 COMPLEX(DP), ALLOCATABLE, TARGET :: fact_hepsi(:,:) 48 !! factors for hermitean electric field operators 49 !COMPLEX(DP), ALLOCATABLE, TARGET :: bec_evcel(:,:) 50 !! for storing bec's factors with evcel 51 TYPE(bec_type) :: bec_evcel 52 INTEGER, ALLOCATABLE, TARGET :: mapgp_global(:,:) 53 !! map for G'= G+1 correspondence 54 INTEGER, ALLOCATABLE, TARGET :: mapgm_global(:,:) 55 !! map for G'= G-1 correspondence 56 REAL(DP) :: ion_pol(3) 57 !! the ionic polarization 58 REAL(DP) :: el_pol(3) 59 !! the electronic polarization 60 REAL(DP) :: fc_pol(3) 61 !! the prefactor for the electronic polarization 62 LOGICAL :: l_el_pol_old 63 !! if true there is already stored a n older value for the polarization 64 !! neeeded for having correct polarization during MD 65 REAL(DP) :: el_pol_old(3) 66 !! the old electronic polarization 67 REAL(DP) :: el_pol_acc(3) 68 !! accumulator for the electronic polarization 69 ! 70 INTEGER :: nppstr_3d(3) 71 !! number of element of strings along the reciprocal directions 72 INTEGER, ALLOCATABLE :: nx_el(:,:) 73 !! index for string to k-point map, (nks*nspin,dir=3) 74 LOGICAL :: l3dstring 75 !! if true strings are on the 3 three directions 76 REAL(DP) :: efield_cart(3) 77 !! electric field vector in cartesian units 78 REAL(DP) :: efield_cry(3) 79 !! electric field vector in crystal units 80 REAL(DP) :: transform_el(3,3) 81 !! transformation matrix from cartesian coordinates to normed reciprocal 82 !! space 83 INTEGER, ALLOCATABLE :: mapg_owner(:,:) 84 REAL(DP) :: pdl_tot 85 !! the total phase calculated from bp_c_phase 86 INTEGER :: phase_control 87 !! 0 no control, 1 write, 2 read 88 ! 89 CONTAINS 90 ! 91 !------------------------------------------------------------------------- 92 SUBROUTINE allocate_bp_efield() 93 !----------------------------------------------------------------------- 94 !! Allocate memory for the Berry's phase electric field. 95 ! NOTICE: should be allocated ONLY in parallel case, for gdir=1 or 2 96 ! 97 USE gvect, ONLY : ngm_g 98 ! 99 IMPLICIT NONE 100 ! 101 IF ( lberry .OR. lelfield .OR. lorbm ) THEN 102 ALLOCATE(mapgp_global(ngm_g,3)) 103 ALLOCATE(mapgm_global(ngm_g,3)) 104 ALLOCATE(mapg_owner(2,ngm_g)) 105 ENDIF 106 ! 107 l_el_pol_old = .FALSE. 108 el_pol_acc = 0.d0 109 ! 110 ! 111 RETURN 112 ! 113 END SUBROUTINE allocate_bp_efield 114 ! 115 ! 116 !-------------------------------------------------------------------------- 117 SUBROUTINE deallocate_bp_efield 118 !------------------------------------------------------------------------ 119 !! Deallocate memory used in Berry's phase electric field calculation 120 ! 121 IMPLICIT NONE 122 ! 123 IF ( lberry .OR. lelfield .OR. lorbm ) THEN 124 IF ( ALLOCATED(mapgp_global) ) DEALLOCATE( mapgp_global ) 125 IF ( ALLOCATED(mapgm_global) ) DEALLOCATE( mapgm_global ) 126 IF ( ALLOCATED(nx_el) ) DEALLOCATE( nx_el ) 127 IF ( ALLOCATED(mapg_owner) ) DEALLOCATE( mapg_owner ) 128 ENDIF 129 ! 130 ! 131 RETURN 132 ! 133 END SUBROUTINE deallocate_bp_efield 134 ! 135 ! 136 !------------------------------------------------------------------------- 137 SUBROUTINE bp_global_map 138 !----------------------------------------------------------------------- 139 !! This subroutine sets up the global correspondence map G+1 and G-1. 140 ! 141 USE mp, ONLY : mp_sum 142 USE mp_images, ONLY : me_image, intra_image_comm 143 USE gvect, ONLY : ngm_g, g, ngm, ig_l2g 144 USE fft_base, ONLY : dfftp 145 USE cell_base, ONLY : at 146 ! 147 IMPLICIT NONE 148 ! 149 INTEGER :: ig, mk1,mk2,mk3, idir, imk(3) 150 INTEGER, ALLOCATABLE :: ln_g(:,:,:) 151 INTEGER, ALLOCATABLE :: g_ln(:,:) 152 ! 153 IF ( .NOT.lberry .AND. .NOT. lelfield .AND. .NOT. lorbm ) RETURN 154 ! set up correspondence ln_g ix,iy,iz ---> global g index in 155 ! (for now...) coarse grid 156 ! and inverse realtion global g (coarse) to ix,iy,iz 157 ! 158 ALLOCATE( ln_g(-dfftp%nr1:dfftp%nr1,-dfftp%nr2:dfftp%nr2,-dfftp%nr3:dfftp%nr3) ) 159 ALLOCATE( g_ln(3,ngm_g) ) 160 ! 161 ln_g(:,:,:)=0!it means also not found 162 DO ig=1,ngm 163 mk1 = NINT( g(1,ig)*at(1,1) + g(2,ig)*at(2,1) + g(3,ig)*at(3,1) ) 164 mk2 = NINT( g(1,ig)*at(1,2) + g(2,ig)*at(2,2) + g(3,ig)*at(3,2) ) 165 mk3 = NINT( g(1,ig)*at(1,3) + g(2,ig)*at(2,3) + g(3,ig)*at(3,3) ) 166 ln_g(mk1,mk2,mk3)=ig_l2g(ig) 167 ENDDO 168 ! 169 CALL mp_sum( ln_g(:,:,:), intra_image_comm ) 170 ! 171 g_ln(:,:) = 0 !it means also not found 172 ! 173 DO ig = 1, ngm 174 mk1 = NINT( g(1,ig)*at(1,1) + g(2,ig)*at(2,1) + g(3,ig)*at(3,1) ) 175 mk2 = NINT( g(1,ig)*at(1,2) + g(2,ig)*at(2,2) + g(3,ig)*at(3,2) ) 176 mk3 = NINT( g(1,ig)*at(1,3) + g(2,ig)*at(2,3) + g(3,ig)*at(3,3) ) 177 g_ln(1,ig_l2g(ig)) = mk1 178 g_ln(2,ig_l2g(ig)) = mk2 179 g_ln(3,ig_l2g(ig)) = mk3 180 ENDDO 181 ! 182 CALL mp_sum( g_ln(:,:), intra_image_comm ) 183 ! 184 !loop on direction 185 DO idir = 1, 3 186 ! for every g on global array find G+1 and G-1 and put on 187 DO ig = 1, ngm_g 188 imk(:) = g_ln(:,ig) 189 imk(idir) = imk(idir)+1 190 ! table array 191 mapgp_global(ig,idir) = ln_g(imk(1),imk(2),imk(3)) 192 imk(idir) = imk(idir)-2 193 mapgm_global(ig,idir) = ln_g(imk(1),imk(2),imk(3)) 194 ENDDO 195 ! 196 ENDDO 197 ! 198 mapg_owner = 0 199 DO ig = 1, ngm 200 mapg_owner(1,ig_l2g(ig)) = me_image + 1 201 mapg_owner(2,ig_l2g(ig)) = ig 202 ENDDO 203 ! 204 CALL mp_sum( mapg_owner, intra_image_comm ) 205 ! 206 DEALLOCATE( ln_g, g_ln ) 207 ! 208 ! 209 RETURN 210 ! 211 END SUBROUTINE bp_global_map 212 ! 213END MODULE bp 214