1! 2! Copyright (C) 2001-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!----------------------------------------------------------------------- 9subroutine symdynph_gq_new (xq, phi, s, invs, rtau, irt, nsymq, & 10 nat, irotmq, minus_q) 11 !----------------------------------------------------------------------- 12 ! 13 ! This routine receives as input an unsymmetrized dynamical 14 ! matrix expressed on the crystal axes and imposes the symmetry 15 ! of the small group of q. Furthermore it imposes also the symmetry 16 ! q -> -q+G if present. 17 ! February 2020: Update (A. Urru) to include the symmetry operations 18 ! that require the time reversal operator (meaning that TS is a 19 ! symmetry of the crystal). For more information please see: 20 ! Phys. Rev. B 100, 045115 (2019) 21 ! 22 ! 23 USE kinds, only : DP 24 USE constants, ONLY: tpi 25 USE symm_base, ONLY : t_rev 26 implicit none 27 ! 28 ! The dummy variables 29 ! 30 integer :: nat, s (3, 3, 48), irt (48, nat), invs (48), & 31 nsymq, irotmq 32 ! input: the number of atoms 33 ! input: the symmetry matrices 34 ! input: the rotated of each vector 35 ! input: the small group of q 36 ! input: the inverse of each matrix 37 ! input: the order of the small gro 38 ! input: the rotation sending q -> 39 real(DP) :: xq (3), rtau (3, 48, nat) 40 ! input: the q point 41 ! input: the R associated at each t 42 43 logical :: minus_q 44 ! input: true if a symmetry q->-q+G 45 complex(DP) :: phi (3, 3, nat, nat) 46 ! inp/out: the matrix to symmetrize 47 ! 48 ! local variables 49 ! 50 integer :: isymq, sna, snb, irot, na, nb, ipol, jpol, lpol, kpol, & 51 iflb (nat, nat) 52 ! counters, indices, work space 53 54 real(DP) :: arg 55 ! the argument of the phase 56 57 complex(DP) :: phip (3, 3, nat, nat), work (3, 3), fase, faseq (48) 58 ! work space, phase factors 59 ! 60 ! We start by imposing hermiticity 61 ! 62 do na = 1, nat 63 do nb = 1, nat 64 do ipol = 1, 3 65 do jpol = 1, 3 66 phi (ipol, jpol, na, nb) = 0.5d0 * (phi (ipol, jpol, na, nb) & 67 + CONJG(phi (jpol, ipol, nb, na) ) ) 68 phi (jpol, ipol, nb, na) = CONJG(phi (ipol, jpol, na, nb) ) 69 enddo 70 enddo 71 enddo 72 enddo 73 ! 74 ! If no other symmetry is present we quit here 75 ! 76 if ( (nsymq == 1) .and. (.not.minus_q) ) return 77 ! 78 ! Then we impose the symmetry q -> -q+G if present 79 ! 80 if (minus_q) then 81 do na = 1, nat 82 do nb = 1, nat 83 do ipol = 1, 3 84 do jpol = 1, 3 85 work(:,:) = (0.d0, 0.d0) 86 sna = irt (irotmq, na) 87 snb = irt (irotmq, nb) 88 arg = 0.d0 89 do kpol = 1, 3 90 arg = arg + (xq (kpol) * (rtau (kpol, irotmq, na) - & 91 rtau (kpol, irotmq, nb) ) ) 92 enddo 93 arg = arg * tpi 94 fase = CMPLX(cos (arg), sin (arg) ,kind=DP) 95 do kpol = 1, 3 96 do lpol = 1, 3 97 work (ipol, jpol) = work (ipol, jpol) + & 98 s (ipol, kpol, irotmq) * s (jpol, lpol, irotmq) & 99 * phi (kpol, lpol, sna, snb) * fase 100 enddo 101 enddo 102 phip (ipol, jpol, na, nb) = (phi (ipol, jpol, na, nb) + & 103 CONJG( work (ipol, jpol) ) ) * 0.5d0 104 enddo 105 enddo 106 enddo 107 enddo 108 phi = phip 109 endif 110 111 ! 112 ! Here we symmetrize with respect to the small group of q 113 ! 114 if (nsymq == 1) return 115 116 iflb (:, :) = 0 117 do na = 1, nat 118 do nb = 1, nat 119 if (iflb (na, nb) == 0) then 120 work(:,:) = (0.d0, 0.d0) 121 do isymq = 1, nsymq 122 irot = isymq 123 sna = irt (irot, na) 124 snb = irt (irot, nb) 125 arg = 0.d0 126 do ipol = 1, 3 127 arg = arg + (xq (ipol) * (rtau (ipol, irot, na) - & 128 rtau (ipol, irot, nb) ) ) 129 enddo 130 arg = arg * tpi 131 faseq (isymq) = CMPLX(cos (arg), sin (arg) ,kind=DP) 132 do ipol = 1, 3 133 do jpol = 1, 3 134 do kpol = 1, 3 135 do lpol = 1, 3 136 IF (t_rev(isymq)==1) THEN 137 work (ipol, jpol) = work (ipol, jpol) + & 138 s (ipol, kpol, irot) * s (jpol, lpol, irot) & 139 * CONJG(phi (kpol, lpol, sna, snb) * faseq (isymq)) 140 ELSE 141 work (ipol, jpol) = work (ipol, jpol) + & 142 s (ipol, kpol, irot) * s (jpol, lpol, irot) & 143 * phi (kpol, lpol, sna, snb) * faseq (isymq) 144 ENDIF 145 enddo 146 enddo 147 enddo 148 enddo 149 enddo 150 do isymq = 1, nsymq 151 irot = isymq 152 sna = irt (irot, na) 153 snb = irt (irot, nb) 154 do ipol = 1, 3 155 do jpol = 1, 3 156 phi (ipol, jpol, sna, snb) = (0.d0, 0.d0) 157 do kpol = 1, 3 158 do lpol = 1, 3 159 IF (t_rev(isymq)==1) THEN 160 phi(ipol,jpol,sna,snb)=phi(ipol,jpol,sna,snb) & 161 + s(ipol,kpol,invs(irot))*s(jpol,lpol,invs(irot))& 162 * CONJG(work (kpol, lpol)*faseq (isymq)) 163 ELSE 164 phi(ipol,jpol,sna,snb)=phi(ipol,jpol,sna,snb) & 165 + s(ipol,kpol,invs(irot))*s(jpol,lpol,invs(irot))& 166 * work (kpol, lpol) * CONJG(faseq (isymq) ) 167 ENDIF 168 enddo 169 enddo 170 enddo 171 enddo 172 iflb (sna, snb) = 1 173 enddo 174 endif 175 enddo 176 enddo 177 phi (:, :, :, :) = phi (:, :, :, :) / DBLE(nsymq) 178 return 179end subroutine symdynph_gq_new 180