u_mix.f90 5.99 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! Copyright (c) 2016 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.
!--------------------------------------------------------------------------------

7 8 9 10 11
MODULE m_umix
  USE m_juDFT
  !
  ! mix the old and new density matrix for the lda+U method
  !                                                 gb.2001
Daniel Wortmann's avatar
Daniel Wortmann committed
12 13
  ! --------------------------------------------------------
  ! Extension to multiple U per atom type by G.M. 2017
14
CONTAINS
15
  SUBROUTINE u_mix(input,atoms,n_mmp_in,n_mmp_out)
16

Daniel Wortmann's avatar
Daniel Wortmann committed
17
    USE m_types
18
    USE m_cdn_io
19
    USE m_nmat_rot
Daniel Wortmann's avatar
Daniel Wortmann committed
20 21
    USE m_xmlOutput

22
    ! ... Arguments
Daniel Wortmann's avatar
Daniel Wortmann committed
23

24
    IMPLICIT NONE
25
    TYPE(t_input),INTENT(IN)   :: input
26
    TYPE(t_atoms),INTENT(IN)   :: atoms
27 28
    COMPLEX, INTENT (INOUT)    :: n_mmp_out(-3:3,-3:3,atoms%n_u,input%jspins)
    COMPLEX, INTENT (INOUT)    :: n_mmp_in (-3:3,-3:3,atoms%n_u,input%jspins)
29 30
    !
    ! ... Locals ...
Daniel Wortmann's avatar
Daniel Wortmann committed
31 32 33
    INTEGER j,k,iofl,l,itype,ios,i_u,jsp,lty(atoms%n_u)
    REAL alpha,spinf,gam,del,sum1,sum2,mix_u, uParam, jParam
    REAL    theta(atoms%n_u),phi(atoms%n_u),zero(atoms%n_u)
34
    LOGICAL n_exist
Daniel Wortmann's avatar
Daniel Wortmann committed
35
    CHARACTER(LEN=20)   :: attributes(6)
36
    COMPLEX,ALLOCATABLE :: n_mmp(:,:,:,:)
37 38 39 40 41 42
    !
    ! check for possible rotation of n_mmp
    !
    INQUIRE (file='n_mmp_rot',exist=n_exist)
    IF (n_exist) THEN
       OPEN (68,file='n_mmp_rot',status='old',form='formatted')
Daniel Wortmann's avatar
Daniel Wortmann committed
43 44 45 46 47 48 49 50 51 52 53
       DO i_u = 1, atoms%n_u
          l = atoms%lda_u(i_u)%l
          READ(68,*,iostat=ios) theta(i_u),phi(i_u)
          IF (ios == 0) THEN
             lty(i_u) = l
          ELSE
             IF (i_u == 1)  CALL juDFT_error("ERROR reading n_mmp_rot", calledby ="u_mix")
             theta(i_u) = theta(i_u-1) ; phi(i_u) = phi(i_u-1)
             lty(i_u) = lty(i_u-1)
          END IF
       END DO
54 55
       CLOSE (68)
       zero = 0.0
56
       CALL nmat_rot(zero,-theta,-phi,3,atoms%n_u,input%jspins,lty,n_mmp_out)
Daniel Wortmann's avatar
Daniel Wortmann committed
57 58
    END IF

59
    ! Write out n_mmp_out to out.xml file
Daniel Wortmann's avatar
Daniel Wortmann committed
60 61

    CALL openXMLElementNoAttributes('ldaUDensityMatrix')
62
    DO jsp = 1, input%jspins
Daniel Wortmann's avatar
Daniel Wortmann committed
63 64 65 66 67 68 69 70 71 72 73 74 75 76
       DO i_u = 1, atoms%n_u
          l = atoms%lda_u(i_u)%l
          itype = atoms%lda_u(i_u)%atomType
          uParam = atoms%lda_u(i_u)%u
          jParam = atoms%lda_u(i_u)%j
          attributes = ''
          WRITE(attributes(1),'(i0)') jsp
          WRITE(attributes(2),'(i0)') itype
          WRITE(attributes(3),'(i0)') i_u
          WRITE(attributes(4),'(i0)') l
          WRITE(attributes(5),'(f15.8)') uParam
          WRITE(attributes(6),'(f15.8)') jParam
          CALL writeXMLElementMatrixPoly('densityMatrixFor',&
                                         (/'spin    ','atomType','uIndex  ','l       ','U       ','J       '/),&
77
                                         attributes,n_mmp_out(-l:l,-l:l,i_u,jsp))
Daniel Wortmann's avatar
Daniel Wortmann committed
78 79 80 81
       END DO
    END DO
    CALL closeXMLElement('ldaUDensityMatrix')

82 83 84 85 86
    ! exit subroutine if density matrix does nbot exist
    IF(.NOT.isDensityMatrixPresent()) THEN
       n_mmp_in = n_mmp_out
       RETURN
    END IF
87

88
    IF (input%ldauLinMix) THEN
89

90
       ! mix here straight with given mixing factors
91

92 93
       ALLOCATE (n_mmp(-3:3,-3:3,MAX(1,atoms%n_u),input%jspins))
       n_mmp = CMPLX(0.0,0.0)
94

95 96
       alpha = input%ldauMixParam
       spinf = input%ldauSpinf
97

98 99
       sum1 = 0.0
       IF (input%jspins.EQ.1) THEN
Daniel Wortmann's avatar
Daniel Wortmann committed
100
          DO i_u = 1, atoms%n_u
101 102
             DO j = -3,3
                DO k = -3,3
103
                   sum1 = sum1 + ABS(n_mmp_out(k,j,i_u,1) - n_mmp_in(k,j,i_u,1))
104
                   n_mmp(k,j,i_u,1) = alpha * n_mmp_out(k,j,i_u,1) + (1.0-alpha) * n_mmp_in(k,j,i_u,1)
Daniel Wortmann's avatar
Daniel Wortmann committed
105 106 107
                END DO
             END DO
          END DO
108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127
          WRITE (6,'(a16,f12.6)') 'n_mmp distance =',sum1
       ELSE
          sum2 = 0.0
          gam = 0.5 * alpha * (1.0 + spinf)
          del = 0.5 * alpha * (1.0 - spinf)
          DO i_u = 1,atoms%n_u
             DO j = -3,3
                DO k = -3,3
                   sum1 = sum1 + ABS(n_mmp_out(k,j,i_u,1) - n_mmp_in(k,j,i_u,1))
                   sum2 = sum2 + ABS(n_mmp_out(k,j,i_u,2) - n_mmp_in(k,j,i_u,2))

                   n_mmp(k,j,i_u,1) =       gam * n_mmp_out(k,j,i_u,1) + &
                                      (1.0-gam) * n_mmp_in (k,j,i_u,1) + &
                                            del * n_mmp_out(k,j,i_u,2) - &
                                            del * n_mmp_in (k,j,i_u,2)

                   n_mmp(k,j,i_u,2) =       gam * n_mmp_out(k,j,i_u,2) + &
                                      (1.0-gam) * n_mmp_in (k,j,i_u,2) + &
                                            del * n_mmp_out(k,j,i_u,1) - &
                                            del * n_mmp_in (k,j,i_u,1)
Daniel Wortmann's avatar
Daniel Wortmann committed
128 129
                END DO
             END DO
130 131 132 133 134 135 136 137 138 139 140 141 142 143 144
          END DO
          WRITE (6,'(a23,f12.6)') 'n_mmp distance spin 1 =',sum1
          WRITE (6,'(a23,f12.6)') 'n_mmp distance spin 2 =',sum2
       ENDIF
       n_mmp_in = n_mmp
       DEALLOCATE (n_mmp)
    ELSE ! input%ldauLinMix

       ! only calculate distance

       sum1 = 0.0
       DO i_u = 1, atoms%n_u
          DO j = -3,3
             DO k = -3,3
                sum1 = sum1 + ABS(n_mmp_out(k,j,i_u,1) - n_mmp_in(k,j,i_u,1))
Daniel Wortmann's avatar
Daniel Wortmann committed
145
             END DO
146 147 148 149 150 151 152 153 154 155 156 157
          END DO
       END DO
       IF (input%jspins.EQ.1) THEN
          WRITE (6,'(a16,f12.6)') 'n_mmp distance =',sum1
       ELSE
          sum2 = 0.0
          WRITE (6,'(a23,f12.6)') 'n_mmp distance spin 1 =',sum1
          DO i_u = 1, atoms%n_u
             DO j = -3,3
                DO k = -3,3
                   sum2 = sum2 + ABS(n_mmp_out(k,j,i_u,2) - n_mmp_in(k,j,i_u,2))
                END DO
Daniel Wortmann's avatar
Daniel Wortmann committed
158
             END DO
159 160 161 162 163 164 165 166 167 168
          END DO
          DO j=-3,3
             WRITE(6,'(14f12.6)') (n_mmp_in(k,j,1,2),k=-3,3)
          END DO
          WRITE (6,'(a23,f12.6)') 'n_mmp distance spin 2 =',sum2
          DO j=-3,3
             WRITE(6,'(14f12.6)') (n_mmp_out(k,j,1,2),k=-3,3)
          END DO
       END IF
    END IF ! input%ldauLinMix
169 170 171

  END SUBROUTINE u_mix
END MODULE m_umix