1 ! 2 ! Copyright (C) 2010-2016 Samuel Ponce', Roxana Margine, Carla Verdi, Feliciano Giustino 3 ! Copyright (C) 2007-2009 Jesse Noffsinger, Brad Malone, Feliciano Giustino 4 ! 5 ! This file is distributed under the terms of the GNU General Public 6 ! License. See the file `LICENSE' in the root directory of the 7 ! present distribution, or http://www.gnu.org/copyleft.gpl.txt . 8 ! 9 ! Adapted from the code PH/phq_summary - Quantum-ESPRESSO group 10 ! 11 !----------------------------------------------------------------------- 12 SUBROUTINE epw_summary() 13 !----------------------------------------------------------------------- 14 !! 15 !! Output symmetry informations 16 !! 17 !! 27/03/2018 - Update with recent PHonon/PH/phq_summary.f90 routine - S. Ponce 18 !! 19 !! RM - Nov 2018: Updated based on QE 6.3 20 !! 21 USE kinds, ONLY : DP 22 USE ions_base, ONLY : nat, ityp, atm, tau, ntyp => nsp, amass 23 USE io_global, ONLY : stdout 24 USE cell_base, ONLY : at, bg, ibrav, alat, omega, celldm 25 USE klist, ONLY : lgauss, degauss, ngauss, nkstot, wk 26 USE klist_epw, ONLY : xk_all 27 USE gvect, ONLY : gcutm, ngm 28 USE gvecs, ONLY : dual, doublegrid, gcutms, ngms 29 USE gvecw, ONLY : ecutwfc 30 USE symm_base, ONLY : s, sname, ft, s_axis_to_cart, sr, t_rev 31 USE noncollin_module, ONLY : noncolin 32 USE spin_orb, ONLY : lspinorb, domag 33 USE funct, ONLY : write_dft_name 34 USE epwcom, ONLY : title 35 USE lr_symm_base, ONLY : irotmq, minus_q, nsymq 36 USE control_flags, ONLY : iverbosity 37 USE fft_base, ONLY : dfftp, dffts 38 USE constants, ONLY : amu_ry 39#if defined(__NAG) 40 USE f90_unix_io, ONLY : flush 41#endif 42 ! 43 IMPLICIT NONE 44 ! 45 INTEGER :: i 46 !! generic counter 47 INTEGER :: ipol 48 !! counter on polarizations 49 INTEGER :: apol 50 !! counter on polarizations 51 INTEGER :: na 52 !! counter on atoms 53 INTEGER :: isymq 54 !! counter on symmetries 55 INTEGER :: isym 56 !! counter on symmetries 57 INTEGER :: nsymtot 58 !! counter on symmetries 59 INTEGER :: ik 60 !! counter on k points 61 ! 62 REAL(KIND = DP) :: ft1, ft2, ft3 63 !! fractionary translation 64 REAL(KIND = DP) :: xkg(3) 65 !! k point in crystal coordinates 66 ! 67 WRITE(stdout, 100) title, ibrav, alat, omega, nat, ntyp, ecutwfc, ecutwfc * dual 68100 FORMAT(/,5x,a75,/,/,5x, & 69 & 'bravais-lattice index = ',i12,/,5x, & 70 & 'lattice parameter (a_0) = ',f12.4,' a.u.',/,5x, & 71 & 'unit-cell volume = ',f12.4,' (a.u.)^3',/,5x, & 72 & 'number of atoms/cell = ',i12,/,5x, & 73 & 'number of atomic types = ',i12,/,5x, & 74 & 'kinetic-energy cut-off = ',f12.4,' Ry',/,5x, & 75 & 'charge density cut-off = ',f12.4,' Ry') 76 77 CALL write_dft_name() 78 ! 79 ! Here add a message if this is a noncollinear or a spin_orbit calculation 80 ! 81 IF (noncolin) THEN 82 IF (lspinorb) THEN 83 IF (domag) THEN 84 WRITE(stdout, '(5x,"Noncollinear calculation with spin-orbit",/)') 85 ELSE 86 WRITE(stdout, '(5x,"Non magnetic calculation with spin-orbit",/)') 87 ENDIF 88 ELSE 89 WRITE(stdout, '(5x,"Noncollinear calculation without spin-orbit",/)') 90 ENDIF 91 ELSE 92 WRITE(stdout, '(/)') 93 ENDIF 94 ! 95 ! Description of the unit cell 96 ! 97 WRITE(stdout, '(2(3x,3(2x,"celldm(",i1,")=",f11.5),/))') & 98 (i, celldm(i), i = 1, 6) 99 WRITE(stdout, '(5x, & 100 & "crystal axes: (cart. coord. in units of a_0)",/, & 101 & 3(15x,"a(",i1,") = (",3f8.4," ) ",/ ) )') & 102 & (apol, (at(ipol, apol), ipol = 1, 3), apol = 1, 3) 103 WRITE(stdout, '(5x, & 104 & "reciprocal axes: (cart. coord. in units 2 pi/a_0)",/, & 105 & 3(15x,"b(",i1,") = (",3f8.4," ) ",/ ) )') & 106 & (apol, (bg(ipol, apol), ipol = 1, 3), apol = 1, 3) 107 ! 108 ! Description of the atoms inside the unit cell 109 ! 110 WRITE(stdout, '(/, 5x,"Atoms inside the unit cell: ")') 111 WRITE(stdout, '(/,3x,"Cartesian axes")') 112 WRITE(stdout, '(/,5x,"site n. atom mass ", & 113 & " positions (a_0 units)")') 114 115 WRITE(stdout, '(7x,i2,5x,a6,f8.4," tau(",i2, & 116 & ") = (",3f11.5," )")') & 117 & (na,atm(ityp(na)), amass(ityp(na))/amu_ry, na, & 118 & (tau(ipol,na), ipol = 1, 3), na = 1, nat) 119 ! 120 ! Description of symmetries 121 ! 122 WRITE(stdout, * ) 123 IF (nsymq <= 1 .AND. .NOT. minus_q) THEN 124 WRITE(stdout, '(5x,"No symmetry!")') 125 ELSE 126 IF (minus_q) THEN 127 WRITE(stdout, '(5x,i2," Sym.Ops. (with q -> -q+G )",/)') & 128 nsymq + 1 129 ELSE 130 WRITE(stdout, '(5x,i2," Sym.Ops. (no q -> -q+G )",/)') nsymq 131 ENDIF 132 ENDIF 133 134 IF (iverbosity == 1) THEN 135 WRITE(stdout, '(36x,"s",24x,"frac. trans.")') 136 IF (minus_q) THEN 137 nsymtot = nsymq + 1 138 ELSE 139 nsymtot = nsymq 140 ENDIF 141 DO isymq = 1, nsymtot 142 IF (isymq > nsymq) THEN 143 isym = irotmq 144 WRITE(stdout, '(/,5x,"This transformation sends q -> -q+G")') 145 ELSE 146 isym = isymq 147 ENDIF 148 WRITE(stdout, '(/6x,"isym = ",i2,5x,a45/)') isymq, sname (isym) 149 IF (noncolin .AND. domag) & 150 WRITE(stdout,'(1x, "Time Reversal",i3)') t_rev(isym) 151 ! 152 IF (ft(1, isym)**2 + ft(2, isym)**2 + ft(3, isym)**2 > 1.0d-8) THEN 153 ft1 = at(1, 1) * ft(1, isym) + at(1, 2) * ft(2, isym) + at(1, 3) * ft(3, isym) 154 ft2 = at(2, 1) * ft(1, isym) + at(2, 2) * ft(2, isym) + at(2, 3) * ft(3, isym) 155 ft3 = at(3, 1) * ft(1, isym) + at(3, 2) * ft(2, isym) + at(3, 3) * ft(3, isym) 156 WRITE(stdout, '(1x,"cryst.",3x,"s(",i2,") = (",3i6, & 157 & " ) f =( ",f10.7," )")') isymq, & 158 & (s(1, ipol, isym), ipol = 1, 3), ft(1,isym) 159 WRITE(stdout, '(17x," (",3(i6,5x), & 160 & " ) ( ",f10.7," )")') & 161 & (s(2, ipol, isym), ipol = 1, 3), ft(2,isym) 162 WRITE(stdout, '(17x," (",3(i6,5x), & 163 & " ) ( ",f10.7," )"/)') & 164 & (s(3,ipol,isym), ipol = 1, 3), ft(3,isym) 165 WRITE(stdout, '(1x,"cart.",3x,"s(",i2,") = (",3f11.7, & 166 & " ) f =( ",f10.7," )")') isymq, & 167 & (sr(1,ipol,isym), ipol = 1, 3), ft1 168 WRITE(stdout, '(17x," (",3f11.7, & 169 & " ) ( ",f10.7," )")') & 170 & (sr(2,ipol,isym), ipol = 1, 3), ft2 171 WRITE(stdout, '(17x," (",3f11.7, & 172 & " ) ( ",f10.7," )"/)') & 173 & (sr(3,ipol,isym), ipol = 1, 3), ft3 174 ELSE 175 WRITE(stdout, '(1x,"cryst.",3x,"s(",i2,") = (",3i6, & 176 & " )")') isymq, (s(1,ipol,isym), ipol = 1, 3) 177 WRITE(stdout, '(17x," (",3(i6,5x)," )")') (s(2,ipol,isym), ipol = 1, 3) 178 WRITE(stdout, '(17x," (",3(i6,5x)," )"/)') (s(3,ipol,isym), ipol = 1, 3) 179 WRITE(stdout, '(1x,"cart.",3x,"s(",i2,") = (",3f11.7, & 180 & " )")') isymq, (sr(1,ipol,isym), ipol = 1, 3) 181 WRITE(stdout, '(17x," (",3f11.7," )")') (sr(2,ipol,isym), ipol = 1, 3) 182 WRITE(stdout, '(17x," (",3f11.7," )"/)') (sr(3,ipol,isym), ipol = 1, 3) 183 ENDIF 184 ENDDO 185 ENDIF 186 ! 187 ! Description of the reciprocal lattice vectors 188 ! 189 WRITE(stdout, '(/5x,"G cutoff =",f10.4," (", & 190 & i7," G-vectors)"," FFT grid: (",i3, & 191 & ",",i3,",",i3,")")') gcutm, ngm, dfftp%nr1, dfftp%nr2, dfftp%nr3 192 193 IF (doublegrid) WRITE(stdout, '(5x,"G cutoff =",f10.4," (", & 194 & i7," G-vectors)"," smooth grid: (",i3, & 195 & ",",i3,",",i3,")")') gcutms, ngms, dffts%nr1, dffts%nr1, dffts%nr1 196 IF (.NOT.lgauss) THEN 197 WRITE(stdout, '(5x,"number of k points=",i5)') nkstot 198 ELSE 199 WRITE(stdout, '(5x,"number of k points=",i5, & 200 & " gaussian broad. (Ry)=",f8.4,5x, & 201 & "ngauss = ",i3)') nkstot, degauss, ngauss 202 ENDIF 203 IF (iverbosity == 1 .OR. nkstot < 10000) THEN 204 WRITE(stdout, '(23x,"cart. coord. in units 2pi/a_0")') 205 DO ik = 1, nkstot 206 WRITE(stdout, '(8x,"k(",i5,") = (",3f12.7,"), wk =",f12.7)') ik, & 207 (xk_all(ipol, ik) , ipol = 1, 3), wk(ik) 208 ENDDO 209 ENDIF 210 IF (iverbosity == 1) THEN 211 WRITE(stdout, '(/23x,"cryst. coord.")') 212 DO ik = 1, nkstot 213 DO ipol = 1, 3 214 xkg(ipol) = at(1, ipol) * xk_all(1, ik) + at(2, ipol) * xk_all(2, ik) & 215 + at(3, ipol) * xk_all(3, ik) 216 ! xkg are the components of xk in the reciprocal lattice basis 217 ENDDO 218 WRITE(stdout, '(8x,"k(",i5,") = (",3f12.7,"), wk =",f12.7)') & 219 ik, (xkg(ipol), ipol = 1, 3), wk(ik) 220 ENDDO 221 ENDIF 222 ! 223 CALL print_ps_info() 224 ! 225 FLUSH(stdout) 226 ! 227 !----------------------------------------------------------------------- 228 END SUBROUTINE epw_summary 229 !----------------------------------------------------------------------- 230