1!
2! Copyright (C) 2001-2016 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 summarize_epsilon()
10  !-----------------------------------------------------------------------
11  !
12  !      write the dielectric tensor on output
13  !
14
15  USE kinds, only : DP
16  USE io_global,  ONLY : stdout
17  USE constants, ONLY: fpi, bohr_radius_angs
18  USE cell_base, ONLY: omega
19  USE noncollin_module, ONLY : npol
20  USE efield_mod, ONLY : epsilon
21  USE control_ph, ONLY : lgamma_gamma, lnoloc, done_epsil
22  USE control_lr, ONLY : lrpa
23
24  IMPLICIT NONE
25
26  INTEGER :: ipol, jpol
27  ! counter on polarizations
28  ! counter on records
29  ! counter on k points
30  REAL(DP) :: chi(3,3)
31  !
32  IF (.NOT. done_epsil) RETURN
33  !
34  IF (lnoloc) THEN
35      WRITE( stdout, '(/,10x,"Dielectric constant in cartesian axis (DV_Hxc=0)",/)')
36  ELSE IF (lrpa) THEN
37     WRITE( stdout, '(/,10x,"RPA dielectric constant in cartesian axis (DV_xc=0)",/)')
38  ELSE
39      WRITE( stdout, '(/,10x,"Dielectric constant in cartesian axis ",/)')
40  ENDIF
41
42  WRITE( stdout, '(10x,"(",3f18.9," )")') ((epsilon(ipol,jpol), ipol=1,3), jpol=1,3)
43
44  IF (lgamma_gamma) THEN
45!
46! The system is probably a molecule. Try to estimate the polarizability
47!
48     DO ipol=1,3
49        DO jpol=1,3
50           IF (ipol == jpol) THEN
51              chi(ipol,jpol) = (epsilon(ipol,jpol)-1.0_DP)*3.0_DP*omega/fpi &
52                               /(epsilon(ipol,jpol)+2.0_DP)
53           ELSE
54              chi(ipol,jpol) = epsilon(ipol,jpol)*omega/fpi
55           END IF
56        END DO
57     END DO
58
59     WRITE(stdout,'(/5x,"Polarizability (a.u.)^3",20x,"Polarizability (A^3)")')
60     WRITE(stdout,'(3f10.2,5x,3f14.4)') ( (chi(ipol,jpol), jpol=1,3), &
61                   (chi(ipol,jpol)*bohr_radius_angs**3, jpol=1,3), ipol=1,3)
62  ENDIF
63
64  RETURN
65END SUBROUTINE summarize_epsilon
66!
67!-----------------------------------------------------------------------
68SUBROUTINE summarize_zeu()
69  !-----------------------------------------------------------------------
70  !
71  !  write the zue effective charges on output
72  !
73  USE kinds,     ONLY : DP
74  USE ions_base, ONLY : nat, ityp, atm
75  USE io_global, ONLY : stdout
76
77  USE efield_mod,   ONLY : zstareu
78  USE control_ph,   ONLY : done_zeu
79
80  IMPLICIT NONE
81
82  INTEGER :: jpol, na
83  ! counters
84  !
85
86  IF (.NOT. done_zeu) RETURN
87
88  WRITE( stdout, '(/,10x,"Effective charges (d Force / dE) in cartesian axis",/)')
89  DO na = 1, nat
90     WRITE( stdout, '(10x," atom ",i6, a6)') na, atm(ityp(na))
91     WRITE( stdout, '(6x,"Ex  (",3f15.5," )")')  (zstareu (1, jpol, na), &
92            jpol = 1, 3)
93     WRITE( stdout, '(6x,"Ey  (",3f15.5," )")')  (zstareu (2, jpol, na), &
94            jpol = 1, 3)
95     WRITE( stdout, '(6x,"Ez  (",3f15.5," )")')  (zstareu (3, jpol, na), &
96            jpol = 1, 3)
97  ENDDO
98
99  RETURN
100END SUBROUTINE summarize_zeu
101
102!-----------------------------------------------------------------------
103SUBROUTINE summarize_zue
104!-----------------------------------------------------------------------
105!
106!  Write the zue effective charges on output
107!
108  USE kinds,      ONLY : DP
109  USE ions_base,  ONLY : nat, atm, ityp
110  USE io_global,  ONLY : stdout
111  USE efield_mod, ONLY : zstarue
112  USE control_ph, ONLY : done_zue
113
114  IMPLICIT NONE
115
116  INTEGER :: ipol, na
117  ! counter on polarization
118  ! counter on atoms
119  !
120  IF (.NOT. done_zue) RETURN
121
122  WRITE( stdout, '(/,10x,"Effective charges (d P / du) in cartesian axis ",/)')
123  !
124  DO na = 1, nat
125     WRITE( stdout, '(10x," atom ",i6,a6)') na, atm(ityp(na))
126     WRITE( stdout, '(6x,"Px  (",3f15.5," )")') (zstarue (ipol, na, 1), &
127            ipol = 1, 3)
128     WRITE( stdout, '(6x,"Py  (",3f15.5," )")') (zstarue (ipol, na, 2), &
129            ipol = 1, 3)
130     WRITE( stdout, '(6x,"Pz  (",3f15.5," )")') (zstarue (ipol, na, 3), &
131            ipol = 1, 3)
132  ENDDO
133  !
134  RETURN
135END SUBROUTINE summarize_zue
136!
137!-----------------------------------------------------------------------
138SUBROUTINE summarize_elopt()
139  !-----------------------------------------------------------------------
140  !
141  ! write the electro-optic tensor on output
142  !
143  USE io_global,  ONLY : stdout
144  USE ramanm,     ONLY : eloptns, done_elop
145IMPLICIT NONE
146
147INTEGER :: ipa, ipb, ipc
148
149IF (.NOT. done_elop) RETURN
150
151WRITE(stdout, '(/,10x,''    Electro-optic tensor is defined as '')' )
152WRITE(stdout, '(10x  ,''  the derivative of the dielectric tensor '')' )
153WRITE(stdout, '(10x  ,''    with respect to one electric field '')' )
154WRITE(stdout, '(10x  ,''       units are Rydberg a.u. '',/)' )
155WRITE(stdout, '(10x  ,''  to obtain the static chi^2 multiply by 1/2  '',/)' )
156WRITE(stdout, '(10x  ,''  to convert to pm/Volt multiply per 2.7502  '',/)' )
157WRITE(stdout, '(/,10x,''Electro-optic tensor in cartesian axis: '',/)' )
158
159DO ipc = 1, 3
160   DO ipb = 1, 3
161      WRITE(stdout,'(10x,''('',3f18.9,'' )'')')      &
162                   (eloptns (ipa, ipb, ipc), ipa = 1, 3)
163   ENDDO
164   WRITE(6,'(10x)')
165ENDDO
166
167RETURN
168END SUBROUTINE summarize_elopt
169
170SUBROUTINE summarize_fpol ()
171  !-----------------------------------------------------------------------
172
173  USE kinds, only : DP
174  USE freq_ph,    ONLY : nfs, done_iu
175
176  IMPLICIT NONE
177  ! input variables
178
179  INTEGER :: iu
180  ! counter on frequencies
181  !
182  DO iu=nfs, 1, -1
183     IF ( done_iu(iu) ) CALL write_polariz(iu)
184  END DO
185
186  RETURN
187END SUBROUTINE summarize_fpol
188!
189