TABLE OF CONTENTS
ABINIT/symgamma [ Functions ]
NAME
symgamma
FUNCTION
Symmetrize perturbations wrt atoms and reduced directions for a fixed qpoint.
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
elph_ds = datastructure for elph data (dimensions and eventually data) kphon_full2full = mapping btw kpoints under symops h1_mat_el = irreducible matrix elements to be completed indsym = mapping of atoms under symops natom = number of atoms nsym = number of syms symq = flags for symmetry elements conserving the present qpoint symrec = symmetry operations in reciprocal space
OUTPUT
h1_mat_el = changed on output
NOTES
scheduled for destruction: not called at present (2/2010) so should eventually dissappear
PARENTS
CHILDREN
SOURCE
40 #if defined HAVE_CONFIG_H 41 #include "config.h" 42 #endif 43 44 #include "abi_common.h" 45 46 47 subroutine symgamma(elph_ds,kphon_full2full,h1_mat_el,& 48 & indsym,natom,nsym,symq,symrec) 49 50 use defs_basis 51 use defs_elphon 52 use m_profiling_abi 53 use m_io_tools 54 55 !This section has been created automatically by the script Abilint (TD). 56 !Do not modify the following lines by hand. 57 #undef ABI_FUNC 58 #define ABI_FUNC 'symgamma' 59 !End of the abilint section 60 61 implicit none 62 63 !Arguments ------------------------------------ 64 !scalars 65 integer,intent(in) :: natom,nsym 66 type(elph_type),intent(in) :: elph_ds 67 !arrays 68 integer,intent(in) :: kphon_full2full(2,nsym,elph_ds%k_phon%nkpt),indsym(4,nsym,natom) 69 integer,intent(in) :: symq(4,2,nsym),symrec(3,3,nsym) 70 real(dp),intent(inout) :: h1_mat_el(2,elph_ds%nFSband,elph_ds%nFSband,elph_ds%nbranch,elph_ds%k_phon%nkpt) 71 72 !Local variables------------------------------- 73 !scalars 74 integer :: ikpt_phon,ipreatom,isym,isymatom,itim 75 integer :: nsymperts,symikpt_phon 76 real(dp) :: timsign 77 !arrays 78 real(dp) :: dsymrec(3,3) 79 real(dp) :: sym_mat_el(2,elph_ds%nFSband,elph_ds%nFSband,elph_ds%nbranch,elph_ds%k_phon%nkpt) 80 81 ! ************************************************************************* 82 83 !2. symmetrize the whole set of perturbations 84 sym_mat_el(:,:,:,:,:) = zero 85 86 nsymperts = 0 87 ! 88 !symrel(isym) sends ipreatom onto isymatom 89 !symrec(isym) sends ipredir onto isymdir 90 !symrec(isym) sends ikpt_phon onto symikpt_phon 91 ! 92 do isym=1,nsym 93 94 do itim=0,1 95 if (symq(4,itim+1,isym) == 0) cycle 96 timsign = one - two*itim 97 98 dsymrec(:,:) = timsign*dble(symrec(:,:,isym)) 99 100 nsymperts = nsymperts+1 101 102 ! loop over image perturbations 103 do isymatom=1,natom 104 ipreatom = indsym(4,isym,isymatom) 105 write(std_out,*) ' symgamma : ', isym,itim,isymatom,ipreatom,dsymrec 106 107 108 do symikpt_phon=1,elph_ds%k_phon%nkpt 109 ikpt_phon = kphon_full2full(itim+1,isym,symikpt_phon) 110 111 ! Real part 112 sym_mat_el (1,:,:,3*(isymatom-1)+1,symikpt_phon) = & 113 & sym_mat_el (1,:,:,3*(isymatom-1)+1,symikpt_phon) & 114 & + dsymrec(1,1)*h1_mat_el(1,:,:,3*(ipreatom-1)+1,ikpt_phon) & 115 & + dsymrec(1,2)*h1_mat_el(1,:,:,3*(ipreatom-1)+2,ikpt_phon) & 116 & + dsymrec(1,3)*h1_mat_el(1,:,:,3*(ipreatom-1)+3,ikpt_phon) 117 sym_mat_el (1,:,:,3*(isymatom-1)+2,symikpt_phon) = & 118 & sym_mat_el (1,:,:,3*(isymatom-1)+2,symikpt_phon) & 119 & + dsymrec(2,1)*h1_mat_el(1,:,:,3*(ipreatom-1)+1,ikpt_phon) & 120 & + dsymrec(2,2)*h1_mat_el(1,:,:,3*(ipreatom-1)+2,ikpt_phon) & 121 & + dsymrec(2,3)*h1_mat_el(1,:,:,3*(ipreatom-1)+3,ikpt_phon) 122 sym_mat_el (1,:,:,3*(isymatom-1)+3,symikpt_phon) = & 123 & sym_mat_el (1,:,:,3*(isymatom-1)+3,symikpt_phon) & 124 & + dsymrec(3,1)*h1_mat_el(1,:,:,3*(ipreatom-1)+1,ikpt_phon) & 125 & + dsymrec(3,2)*h1_mat_el(1,:,:,3*(ipreatom-1)+2,ikpt_phon) & 126 & + dsymrec(3,3)*h1_mat_el(1,:,:,3*(ipreatom-1)+3,ikpt_phon) 127 ! Imag part 128 sym_mat_el (2,:,:,3*(isymatom-1)+1,symikpt_phon) = & 129 & sym_mat_el (2,:,:,3*(isymatom-1)+1,symikpt_phon) & 130 & + (dsymrec(1,1)*h1_mat_el(2,:,:,3*(ipreatom-1)+1,ikpt_phon) & 131 & + dsymrec(1,2)*h1_mat_el(2,:,:,3*(ipreatom-1)+2,ikpt_phon) & 132 & + dsymrec(1,3)*h1_mat_el(2,:,:,3*(ipreatom-1)+3,ikpt_phon)) 133 sym_mat_el (2,:,:,3*(isymatom-1)+2,symikpt_phon) = & 134 & sym_mat_el (2,:,:,3*(isymatom-1)+2,symikpt_phon) & 135 & + (dsymrec(2,1)*h1_mat_el(2,:,:,3*(ipreatom-1)+1,ikpt_phon) & 136 & + dsymrec(2,2)*h1_mat_el(2,:,:,3*(ipreatom-1)+2,ikpt_phon) & 137 & + dsymrec(2,3)*h1_mat_el(2,:,:,3*(ipreatom-1)+3,ikpt_phon)) 138 sym_mat_el (2,:,:,3*(isymatom-1)+3,symikpt_phon) = & 139 & sym_mat_el (2,:,:,3*(isymatom-1)+3,symikpt_phon) & 140 & + (dsymrec(3,1)*h1_mat_el(2,:,:,3*(ipreatom-1)+1,ikpt_phon) & 141 & + dsymrec(3,2)*h1_mat_el(2,:,:,3*(ipreatom-1)+2,ikpt_phon) & 142 & + dsymrec(3,3)*h1_mat_el(2,:,:,3*(ipreatom-1)+3,ikpt_phon)) 143 end do 144 end do 145 end do 146 end do 147 !end isym and itim do 148 149 !commented to use un-symmetrized version of h1_mat_el 150 h1_mat_el(:,:,:,:,:) = sym_mat_el(:,:,:,:,:) / nsymperts 151 152 end subroutine symgamma