brysh1.f90 6.58 KB
Newer Older
1
2
3
4
5
6
7
8
9
MODULE m_brysh1
  USE m_juDFT
  !******************************************************
  !      shifts the charge density of the interstitial, m.t
  !      and vacuum part in one single vector
  !      in the spin polarized case the arrays consist of 
  !      spin up and spin down densities
  !******************************************************
CONTAINS
Daniel Wortmann's avatar
Daniel Wortmann committed
10
  SUBROUTINE brysh1(input,stars,atoms,sphhar,noco,vacuum,sym,oneD,&
11
                    intfac,vacfac,den,nmap,nmaph,mapmt,mapvac,mapvac2,sout) 
12
13
14
15
16
17
18
19
20
21
22

    USE m_types
    IMPLICIT NONE
    TYPE(t_oneD),INTENT(IN)    :: oneD
    TYPE(t_input),INTENT(IN)   :: input
    TYPE(t_vacuum),INTENT(IN)  :: vacuum
    TYPE(t_noco),INTENT(IN)    :: noco
    TYPE(t_sym),INTENT(IN)     :: sym
    TYPE(t_stars),INTENT(IN)   :: stars
    TYPE(t_sphhar),INTENT(IN)  :: sphhar
    TYPE(t_atoms),INTENT(IN)   :: atoms
23
24
25
26
27
28
29
30
31
32
    TYPE(t_potden),INTENT(IN)  :: den

    ! Scalar Arguments
    REAL,    INTENT (IN)  :: intfac,vacfac
    INTEGER, INTENT (OUT) :: mapmt,mapvac,mapvac2,nmap,nmaph

    ! Array Arguments
    REAL,    INTENT (OUT) :: sout(:)

    ! Local Scalars
33
    INTEGER i,iv,j,js,k,l,n,nall,na,nvaccoeff,nvaccoeff2,mapmtd
34

35
36
37
38
39
40
41
42
    !--->  put input into arrays sout 
    !      in the spin polarized case the arrays consist of 
    !      spin up and spin down densities

    j=0
    DO  js = 1,input%jspins
       DO i = 1,stars%ng3
          j = j + 1
43
          sout(j) = REAL(den%pw(i,js))
44
45
46
47
       END DO
       IF (.NOT.sym%invs) THEN
          DO i = 1,stars%ng3
             j = j + 1
48
             sout(j) = AIMAG(den%pw(i,js))
49
50
51
52
53
54
55
56
57
          END DO
       ENDIF
       mapmt=0
       na = 1
       DO n = 1,atoms%ntype
          DO l = 0,sphhar%nlh(atoms%ntypsy(na))
             DO i = 1,atoms%jri(n)
                mapmt = mapmt +1
                j = j + 1
58
                sout(j) = den%mt(i,l,n,js)
59
60
61
62
63
64
65
66
67
68
             END DO
          END DO
          na = na + atoms%neq(n)
       END DO
       mapvac=0
       IF (input%film) THEN
          DO iv = 1,vacuum%nvac
             DO k = 1,vacuum%nmz
                mapvac = mapvac + 1
                j = j + 1
69
                sout(j) = den%vacz(k,iv,js)
70
71
72
73
74
             END DO
             DO k = 1,oneD%odi%nq2-1
                DO i = 1,vacuum%nmzxy
                   mapvac = mapvac + 1
                   j = j + 1
75
                   sout(j) =  REAL(den%vacxy(i,k,iv,js))
76
77
78
79
80
81
82
                END DO
             END DO
             IF (.NOT.sym%invs2) THEN
                DO k = 1,oneD%odi%nq2-1
                   DO i = 1,vacuum%nmzxy
                      mapvac = mapvac + 1
                      j = j + 1
83
                      sout(j) =  AIMAG(den%vacxy(i,k,iv,js))
84
85
86
87
88
89
90
91
92
93
94
95
96
                   END DO
                END DO
             END IF
          END DO
       END IF
       IF (js .EQ. 1) nmaph = j
    ENDDO

    mapvac2=0
    IF (noco%l_noco) THEN
       !--->    off-diagonal part of the density matrix
       DO i = 1,stars%ng3
          j = j + 1
97
          sout(j) = REAL(den%cdom(i))
98
99
100
       END DO
       DO i = 1,stars%ng3
          j = j + 1
101
          sout(j) = AIMAG(den%cdom(i))
102
103
104
105
106
107
       END DO
       IF (input%film) THEN
          DO iv = 1,vacuum%nvac
             DO k = 1,vacuum%nmz
                mapvac2 = mapvac2 + 1
                j = j + 1
108
                sout(j) = REAL(den%cdomvz(k,iv))
109
110
111
112
113
             END DO
             DO k = 1,oneD%odi%nq2-1
                DO i = 1,vacuum%nmzxy
                   mapvac2 = mapvac2 + 1
                   j = j + 1
114
                   sout(j) =  REAL(den%cdomvxy(i,k,iv))
115
116
117
118
119
120
121
                END DO
             END DO
          END DO
          DO iv = 1,vacuum%nvac
             DO k = 1,vacuum%nmz
                mapvac2 = mapvac2 + 1
                j = j + 1
122
                sout(j) = AIMAG(den%cdomvz(k,iv))
123
124
125
126
127
             END DO
             DO k = 1,oneD%odi%nq2-1
                DO i = 1,vacuum%nmzxy
                   mapvac2 = mapvac2 + 1
                   j = j + 1
128
                   sout(j) =  AIMAG(den%cdomvxy(i,k,iv))
129
130
131
132
133
134
135
136
137
138
                END DO
             END DO
          END DO
          nvaccoeff2 = 2*vacuum%nmzxy*(oneD%odi%nq2-1)*vacuum%nvac + 2*vacuum%nmz*vacuum%nvac
          IF (mapvac2 .NE. nvaccoeff2) THEN
             WRITE (6,*)'The number of vaccum coefficients off the'
             WRITE (6,*)'off-diagonal part of the density matrix is'
             WRITE (6,*)'inconsitent:'
             WRITE (6,8000) mapvac2,nvaccoeff2
8000         FORMAT ('mapvac2= ',i12,'nvaccoeff2= ',i12)
Daniel Wortmann's avatar
Daniel Wortmann committed
139
             CALL juDFT_error("brysh1:# of vacuum coeff. inconsistent" ,calledby ="brysh1")
140
141
142
143
144
145
146
147
148
149
          ENDIF
       END IF
    ENDIF ! noco

    IF (atoms%n_u > 0 ) THEN     ! lda+U
       DO js = 1,input%jspins
          DO n = 1, atoms%n_u
             DO k = -3, 3
                DO i = -3, 3
                   j = j + 1 
150
                   sout(j) = REAL(den%mmpMat(i,k,n,js))
151
                   j = j + 1 
152
                   sout(j) = AIMAG(den%mmpMat(i,k,n,js))
153
154
155
156
157
158
159
160
161
162
163
164
165
                ENDDO
             ENDDO
          ENDDO
       ENDDO
    ENDIF

    IF (input%film) THEN
       nvaccoeff = vacfac*vacuum%nmzxy*(oneD%odi%nq2-1)*vacuum%nvac + vacuum%nmz*vacuum%nvac
       IF (mapvac .NE. nvaccoeff) THEN
          WRITE(6,*)'The number of vaccum coefficients is'
          WRITE(6,*)'inconsitent:'
          WRITE (6,8010) mapvac,nvaccoeff
8010      FORMAT ('mapvac= ',i12,'nvaccoeff= ',i12)
Daniel Wortmann's avatar
Daniel Wortmann committed
166
          CALL juDFT_error("brysh1: # of vacuum coeff. inconsistent" ,calledby ="brysh1")
167
168
169
       ENDIF
    ENDIF

Daniel Wortmann's avatar
Daniel Wortmann committed
170
    mapmtd = atoms%ntype*(sphhar%nlhd+1)*atoms%jmtd
171
172
173
174
175
    IF (mapmt .GT. mapmtd) THEN
       WRITE(6,*)'The number of mt coefficients is larger than the'
       WRITE(6,*)'dimensions:'
       WRITE (6,8040) mapmt,mapmtd
8040   FORMAT ('mapmt= ',i12,' > mapmtd= ',i12)
Daniel Wortmann's avatar
Daniel Wortmann committed
176
       CALL juDFT_error("brysh1: mapmt > mapmtd (dimensions)",calledby ="brysh1")
177
178
179
180
181
182
    ENDIF

    nmap = j
    nall = (intfac*stars%ng3 + mapmt + mapvac + 49*2*atoms%n_u )*input%jspins
    IF (noco%l_noco) nall = nall + 2*stars%ng3 + mapvac2
    IF (nall.NE.nmap) THEN
Daniel Wortmann's avatar
Daniel Wortmann committed
183
       WRITE(6,*)'The total number of charge density coefficients is'
184
185
186
187
       WRITE(6,*)'inconsitent:'
       WRITE (6,8020) nall,nmap
8020   FORMAT ('nall= ',i12,'not equal nmap= ',i12)
       WRITE (6,'(a,i5,a,i5)') 'nall = ',nall,' nmap = ',nmap
Daniel Wortmann's avatar
Daniel Wortmann committed
188
       CALL juDFT_error ("brysh1: input # of charge density coeff. inconsistent" ,calledby ="brysh1")
189
190
    ENDIF
    IF (nmap.GT.SIZE(sout)) THEN 
Daniel Wortmann's avatar
Daniel Wortmann committed
191
       WRITE(6,*)'The total number of charge density coefficients is'
192
193
194
       WRITE(6,*)'larger than the dimensions:'
       WRITE (6,8030) nmap,SIZE(sout)
8030   FORMAT ('nmap= ',i12,' > size(sout)= ',i12)
Daniel Wortmann's avatar
Daniel Wortmann committed
195
       CALL juDFT_error("brysh1: nmap > mmap (dimensions)",calledby ="brysh1")
196
197
198
199
    ENDIF

  END SUBROUTINE brysh1
END MODULE m_brysh1