TABLE OF CONTENTS
ABINIT/completeperts [ Functions ]
NAME
completeperts
FUNCTION
Complete perturbations wrt atoms and reduced directions for a fixed qpoint. Normally there is a test in read_gkk which guarantees that enough irreducible perturbations are present to generate everything. h1_mat_el is first squared, making a (ipert,jpert) matrix which has the same symmetry properties as the dynamical matrix.
COPYRIGHT
Copyright (C) 2004-2018 ABINIT group (MVer) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
Cryst<crystal_t>=Info on the unit cell and symmetries. nbranch=Number of phonon branches. nFSband=Number of bands in H1 matrix elements. nkpt=Number of k-points in matrix elements. nsppol=Number of independent spin polarizations. gkk_flag = flags for presence of gkk matrix elements h1_mat_el = irreducible matrix elements to be completed and squared qpt = qpoint symq = flags for symmetry elements conserving the present qpoint tnons = translation vectors associated with symops
OUTPUT
h1_mat_el_sq = irreducible matrix elements squared and completed gkk_flag = changed on output
PARENTS
read_gkk
CHILDREN
d2sym3
SOURCE
45 #if defined HAVE_CONFIG_H 46 #include "config.h" 47 #endif 48 49 #include "abi_common.h" 50 51 subroutine completeperts(Cryst,nbranch,nFSband,nkpt,nsppol,gkk_flag,h1_mat_el,h1_mat_el_sq,& 52 & qpt,symq,qtimrev) 53 54 use defs_basis 55 use defs_elphon 56 use m_profiling_abi 57 use m_errors 58 59 use m_dynmat, only : d2sym3 60 use m_crystal, only : crystal_t 61 62 !This section has been created automatically by the script Abilint (TD). 63 !Do not modify the following lines by hand. 64 #undef ABI_FUNC 65 #define ABI_FUNC 'completeperts' 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: qtimrev,nbranch,nFSband,nkpt,nsppol 73 type(crystal_t),intent(in) :: Cryst 74 !arrays 75 integer,intent(in) :: symq(4,2,Cryst%nsym) 76 integer,intent(inout) :: gkk_flag(nbranch,nbranch,nkpt,nsppol) 77 real(dp),intent(in) :: qpt(3) 78 real(dp),intent(in) :: h1_mat_el(2,nFSband**2,nbranch,nkpt,nsppol) 79 real(dp),intent(out) :: h1_mat_el_sq(2,nFSband**2,nbranch**2,nkpt,nsppol) 80 81 !Local variables------------------------------- 82 !scalars 83 integer :: ikpt_phon,iatom1,iatom2,ibb,idir1,idir2,ipert1,ipert2,isppol,mpert,natom 84 real(dp) :: im1,im2,re1,re2,res 85 character(len=500) :: msg 86 !arrays 87 integer,allocatable :: tmpflg(:,:,:,:) 88 real(dp),allocatable :: tmpval(:,:,:,:,:) 89 90 ! ************************************************************************* 91 92 !WARNING! Stupid patch in d2sym3 imposes these matrices to have size natom+2 93 natom = Cryst%natom 94 mpert = natom+2 95 96 ABI_ALLOCATE(tmpflg,(3,mpert,3,mpert)) 97 ABI_ALLOCATE(tmpval,(2,3,mpert,3,mpert)) 98 99 h1_mat_el_sq = zero 100 101 write (std_out,*) ' completeperts: shape(h1_mat_el_sq) = ', shape(h1_mat_el_sq) 102 103 do isppol=1,nsppol 104 write(std_out,*)'completeperts: isppol = ', isppol 105 ! 106 do ikpt_phon=1,nkpt 107 do ibb=1,nFSband**2 108 ! 109 tmpval = zero 110 tmpflg = 0 111 ! for a fixed k (q) band and sppol construct the gamma matrix for (3 natom)^2 perturbation pairs 112 do iatom1=1,natom 113 do idir1=1,3 114 ipert1 = (iatom1-1)*3+idir1 115 if (gkk_flag(ipert1,ipert1,ikpt_phon,isppol) < 0) cycle 116 re1 = h1_mat_el(1,ibb,ipert1,ikpt_phon,isppol) 117 im1 = h1_mat_el(2,ibb,ipert1,ikpt_phon,isppol) 118 119 do iatom2=1,natom 120 do idir2=1,3 121 ipert2 = (iatom2-1)*3+idir2 122 if (gkk_flag(ipert2,ipert2,ikpt_phon,isppol) < 0) cycle 123 tmpflg(idir1,iatom1,idir2,iatom2) = 1 124 re2 = h1_mat_el(1,ibb,ipert2,ikpt_phon,isppol) 125 im2 = h1_mat_el(2,ibb,ipert2,ikpt_phon,isppol) 126 ! 127 ! conjg(h1_mat_el_2) * h1_mat_el_1 128 res = re1*re2 + im1*im2 129 tmpval(1,idir1,iatom1,idir2,iatom2) = res 130 res = re1*im2 - im1*re2 131 tmpval(2,idir1,iatom1,idir2,iatom2) = res 132 133 end do !idir2 134 end do !iatom2 135 end do !idir1 136 end do !iatom1 137 138 ! matrix is symmetrized like a dynamical matrix. No change of band or k 139 ! in here. This should be checked (if we have to restrict further the symmetry operations) 140 call d2sym3(tmpflg,tmpval,Cryst%indsym,mpert,natom,Cryst%nsym,qpt,symq,Cryst%symrec,Cryst%symrel,qtimrev) 141 142 if (sum(tmpflg(:,1:natom,:,1:natom)) /= 3*natom*3*natom) then 143 write(msg,'(3a,4i0)')& 144 & 'A perturbation is missing after completion with d2sym3',ch10,& 145 & 'tmpflg, ikpt_phon, isppol: ',tmpflg,ikpt_phon,isppol 146 MSG_ERROR(msg) 147 end if 148 ! 149 ! Save values for calculation of |gkk|^2 150 do iatom1=1,natom 151 do idir1=1,3 152 ipert1 = (iatom1-1)*3+idir1 153 do iatom2=1,natom 154 do idir2=1,3 155 ! 156 ! mjv 29/10/2007 ipert2 now contains the composite index ip1*nperts+ip2 157 ipert2 = (iatom2-1)*3 + idir2 + (ipert1-1)*3*natom 158 h1_mat_el_sq(1,ibb,ipert2,ikpt_phon,isppol) = pi*tmpval(1,idir2,iatom2,idir1,iatom1) 159 h1_mat_el_sq(2,ibb,ipert2,ikpt_phon,isppol) = pi*tmpval(2,idir2,iatom2,idir1,iatom1) 160 end do 161 end do 162 end do 163 end do 164 ! 165 end do !end ibb band dos 166 ! 167 ! Set flags. 168 do ipert1=1,3*natom 169 do ipert2=1,3*natom 170 if (gkk_flag(ipert2,ipert1,ikpt_phon,isppol) < 0) gkk_flag(ipert2,ipert1,ikpt_phon,isppol) = 1 171 end do 172 end do 173 174 end do !end kpt_phon do 175 end do !end sppol do 176 177 ABI_DEALLOCATE(tmpflg) 178 ABI_DEALLOCATE(tmpval) 179 180 end subroutine completeperts