1! 2! Copyright (C) 2005 PWSCF 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! 9!---------------------------------------------------------------------------- 10SUBROUTINE report_mag 11 !------------------------------------------------------------------------ 12 !! This subroutine prints out information about the local magnetization 13 !! and/or charge, integrated around the atomic positions at points which 14 !! are calculated in make_pointlists. 15 ! 16 USE kinds, ONLY: DP 17 USE ions_base, ONLY: nat, tau, ityp 18 USE io_global, ONLY: stdout 19 use constants, ONLY: pi 20 USE scf, ONLY: rho 21 USE noncollin_module, ONLY: noncolin, mcons, i_cons 22 USE lsda_mod, ONLY: nspin 23 ! 24 IMPLICIT NONE 25 ! 26 REAL(DP) :: theta, phi, norm, norm1 27 INTEGER :: ipol, iat 28 REAL(DP) :: r1_loc(nat), m1_loc(nspin-1,nat) 29 ! 30 ! get_local integrates on the previously determined points 31 ! 32 CALL get_locals( r1_loc, m1_loc, rho%of_r ) 33 ! 34 IF (nspin == 2) THEN 35 WRITE( stdout, * ) 36 WRITE( stdout, '(5X,"Magnetic moment per site:")' ) 37 DO iat = 1, nat 38 WRITE(stdout,1020) iat, r1_loc(iat), m1_loc(1,iat), mcons(1,ityp(iat)) 39 END DO 40 ! 41 ELSE IF (noncolin) THEN 42 DO iat = 1, nat 43 ! 44 ! norm is the length of the magnetic moment vector 45 ! 46 norm = DSQRT( m1_loc(1,iat)**2+m1_loc(2,iat)**2+m1_loc(3,iat)**2 ) 47 ! 48 ! norm1 is the length of the projection of the mm vector into 49 ! the xy plane 50 ! 51 norm1 = DSQRT( m1_loc(1,iat)**2+m1_loc(2,iat)**2 ) 52 ! 53 ! calculate the polar angles of the magnetic moment 54 ! 55 IF (norm > 1.d-10) THEN 56 theta = ACOS(m1_loc(3,iat)/norm) 57 IF (norm1 > 1.d-10) THEN 58 phi = ACOS(m1_loc(1,iat)/norm1) 59 IF (m1_loc(2,iat) < 0.d0) phi = - phi 60 ELSE 61 phi = 2.d0*pi 62 ENDIF 63 ELSE 64 theta = 2.d0*pi 65 phi = 2.d0*pi 66 ENDIF 67 ! 68 ! go to degrees 69 theta = theta*180.d0/pi 70 phi = phi*180.d0/pi 71 ! 72 WRITE( stdout,1010) 73 WRITE( stdout,1011) iat,(tau(ipol,iat),ipol=1,3) 74 WRITE( stdout,1014) r1_loc (iat) 75 WRITE( stdout,1012) (m1_loc(ipol,iat),ipol=1,3) 76 WRITE( stdout,1018) (m1_loc(ipol,iat)/r1_loc(iat),ipol=1,3) 77 WRITE( stdout,1013) norm,theta,phi 78 IF (i_cons==1) THEN 79 WRITE( stdout,1015) (mcons(ipol,ityp(iat)),ipol=1,3) 80 ELSEIF (i_cons==2) THEN 81 WRITE( stdout,1017) 180.d0 * ACOS(mcons(3,ityp(iat)))/pi 82 ENDIF 83 WRITE( stdout,1010) 84 ENDDO 85 ! 86 END IF 87 ! 88 ! 89 1010 FORMAT (/,1x,78('=')) 90 1011 FORMAT (5x,'atom number ',i4,' relative position : ',3f9.4) 91 1012 FORMAT (5x,'magnetization : ',3f12.6) 92 1013 FORMAT (5x,'polar coord.: r, theta, phi [deg] : ',3f12.6) 93 1014 FORMAT (5x,'charge : ',f12.6) 94 1018 FORMAT (5x,'magnetization/charge:',3f12.6) 95 1015 FORMAT (5x,'constrained moment : ',3f12.6) 96 1017 FORMAT (5x,'constrained theta [deg] : ',f12.6) 97 1020 FORMAT (5x,'atom: ',i4,4X,'charge: ',F9.4,4X,'magn: ',F9.4,4X,'constr: ',f9.4) 98 ! 99END SUBROUTINE report_mag 100