magMoms.f90 2.88 KB
Newer Older
1 2 3 4 5 6 7 8 9 10
!--------------------------------------------------------------------------------
! Copyright (c) 2018 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

MODULE m_magMoms

CONTAINS

11
SUBROUTINE magMoms(dimension,input,atoms,noco,vTot,moments)
12 13 14

   USE m_types
   USE m_xmlOutput
15
   USE m_m_perp
16 17 18 19 20 21

   IMPLICIT NONE

   TYPE(t_dimension), INTENT(IN) :: dimension
   TYPE(t_input), INTENT(IN)     :: input
   TYPE(t_atoms), INTENT(IN)     :: atoms
22 23
   TYPE(t_noco), INTENT(INOUT)   :: noco
   TYPE(t_potden),INTENT(IN)     :: vTot
24
   TYPE(t_moments),INTENT(IN)    :: moments
25

26
   INTEGER                       :: iType, j, iRepAtom
27
   REAL                          :: sval,stot,scor,smom
28 29
   CHARACTER(LEN=20)             :: attributes(4)

30 31
   WRITE (6,FMT=8000)
   DO iType = 1,atoms%ntype
32 33
      sval = moments%svdn(iType,1) - moments%svdn(iType,input%jspins)
      stot = moments%stdn(iType,1) - moments%stdn(iType,input%jspins)
34
      scor = stot - sval
35
      WRITE (6,FMT=8010) iType,stot,sval,scor,moments%svdn(iType,1),moments%stdn(iType,1)
36 37 38
   END DO

   8000 FORMAT (/,/,10x,'spin density at the nucleus:',/,10x,'type',t25,&
39 40
                'total',t42,'valence',t65,'core',t90,&
                'majority valence and total density',/)
41 42
   8010 FORMAT (i13,2x,3e20.8,5x,2e20.8)

43
   WRITE (6,FMT=8020)
44
  
45
   CALL openXMLElement('magneticMomentsInMTSpheres',(/'units'/),(/'muBohr'/))
46
   iRepAtom = 1
47
   DO iType = 1, atoms%ntype
48 49
      smom = moments%chmom(iType,1) - moments%chmom(iType,input%jspins)
      WRITE (6,FMT=8030) iType,smom, (moments%chmom(iType,j),j=1,input%jspins)
50 51 52
      attributes = ''
      WRITE(attributes(1),'(i0)') iType
      WRITE(attributes(2),'(f15.10)') smom
53 54
      WRITE(attributes(3),'(f15.10)') moments%chmom(iType,1)
      WRITE(attributes(4),'(f15.10)') moments%chmom(iType,2)
55 56 57
      CALL writeXMLElementFormPoly('magneticMoment',(/'atomType      ','moment        ','spinUpCharge  ',&
                                                      'spinDownCharge'/),&
                                   attributes,reshape((/8,6,12,14,6,15,15,15/),(/4,2/)))
58 59 60 61 62

      IF (noco%l_mperp) THEN
         !calculate the perpendicular part of the local moment
         !and relax the angle of the local moment or calculate
         !the constraint B-field.
63
         CALL m_perp(atoms,iType,iRepAtom,noco,vTot%mt(:,0,:,:),moments%chmom,moments%qa21)
64 65
      END IF
      iRepAtom= iRepAtom + atoms%neq(iType)
66 67 68 69 70 71 72 73 74 75
   END DO
   CALL closeXMLElement('magneticMomentsInMTSpheres')

   8020 FORMAT (/,/,2x,'-->  magnetic moments in the spheres:',/,2x,&
                'mm -->   type',t22,'moment',t33,'spin-up',t43,'spin-down')
   8030 FORMAT (2x,'--> mm',i8,2x,3f12.5)

END SUBROUTINE magMoms

END MODULE m_magMoms