Line data Source code
1 : MODULE m_mkgz
2 : c.....------------------------------------------------------------------
3 : c by use of cartesian x,y,z components of charge density gradients,
4 : c make the quantities
5 : cc agrt,agru,agrd,g2rt,g2ru,g2rd,gggrt,gggru,gggrd,
6 : cc gzgr
7 : cc used to calculate gradient contribution to xc potential and
8 : cc energy.
9 : c.....------------------------------------------------------------------
10 : CONTAINS
11 5 : SUBROUTINE mkgz(
12 5 : > nmzdf,jspins,rh1,rh2,rhdz1,rhdz2,rhdzz1,rhdzz2,
13 : < grad)
14 : USE m_types
15 : IMPLICIT NONE
16 : INTEGER, INTENT (IN) :: nmzdf,jspins
17 : REAL, INTENT (IN) :: rh1(nmzdf),rhdz1(nmzdf),rhdzz1(nmzdf)
18 : REAL, INTENT (IN) :: rh2(nmzdf),rhdz2(nmzdf),rhdzz2(nmzdf)
19 : TYPE(t_gradients),INTENT(INOUT)::grad
20 :
21 : INTEGER i
22 : REAL vlt,dvzt,dvzzt,vlu,dvzu,dvzzu,vld,dvzd,dvzzd
23 : REAL dagrzt,dagrzu,dagrzd,dztadz,sml
24 :
25 5 : sml = 1.e-14
26 :
27 5 : IF (ALLOCATED(grad%sigma)) THEN
28 0 : IF(jspins==1) THEN
29 0 : DO i=1,nmzdf
30 0 : grad%sigma(1,i)=rhdz1(i)*rhdz1(i)
31 : ENDDO
32 : ELSE
33 0 : DO i=1,nmzdf
34 0 : grad%sigma(1,i)=rhdz1(i)*rhdz1(i)
35 0 : grad%sigma(2,i)=rhdz1(i)*rhdz2(i)
36 0 : grad%sigma(3,i)=rhdz2(i)*rhdz2(i)
37 : ENDDO
38 : ENDIF
39 : RETURN
40 : ENDIF
41 :
42 5 : IF (jspins == 1) THEN
43 :
44 602 : DO i = 1,nmzdf
45 :
46 300 : vlu = max(rh1(i)/2,sml)
47 300 : dvzu = rhdz1(i)/2
48 300 : dvzzu = rhdzz1(i)/2
49 300 : vld = vlu
50 300 : dvzd = dvzu
51 300 : dvzzd = dvzzu
52 :
53 300 : vlt = vlu+vld
54 300 : dvzt = dvzu+dvzd
55 300 : dvzzt = dvzzu+dvzzd
56 :
57 : c agrt: abs(grad(ro)), u,d for up and down.
58 :
59 300 : grad%agrt(i) = max(abs(dvzt),sml)
60 300 : grad%agru(i) = max(abs(dvzu),sml)
61 300 : grad%agrd(i) = max(abs(dvzd),sml)
62 :
63 300 : dagrzt= dvzt*dvzzt/grad%agrt(i)
64 300 : dagrzu= dvzu*dvzzu/grad%agru(i)
65 300 : dagrzd= dvzd*dvzzd/grad%agrd(i)
66 :
67 300 : grad%gggrt(i) = dvzt*dagrzt
68 300 : grad%gggru(i) = dvzu*dagrzu
69 300 : grad%gggrd(i) = dvzd*dagrzd
70 :
71 : c dztadz=d(zeta)/dz,..
72 :
73 300 : dztadz = (dvzu-dvzd)/vlt - (vlu-vld)*dvzt/vlt**2
74 :
75 : c gzgr=grad(zeta)*grad(ro).
76 :
77 300 : grad%gzgr(i) = dztadz*dvzt
78 :
79 : c g2rt: grad(grad(ro))
80 :
81 300 : grad%g2rt(i) = dvzzt
82 300 : grad%g2ru(i) = dvzzu
83 302 : grad%g2rd(i) = dvzzd
84 :
85 : ENDDO
86 :
87 : ELSE
88 :
89 903 : DO i = 1,nmzdf
90 :
91 450 : vlu = max(rh1(i),sml)
92 450 : dvzu = rhdz1(i)
93 450 : dvzzu = rhdzz1(i)
94 450 : vld = max(rh2(i),sml)
95 450 : dvzd = rhdz2(i)
96 450 : dvzzd = rhdzz2(i)
97 :
98 450 : vlt = vlu+vld
99 450 : dvzt = dvzu+dvzd
100 450 : dvzzt = dvzzu+dvzzd
101 :
102 : c agrt: abs(grad(ro)), u,d for up and down.
103 :
104 450 : grad%agrt(i) = max(abs(dvzt),sml)
105 450 : grad%agru(i) = max(abs(dvzu),sml)
106 450 : grad%agrd(i) = max(abs(dvzd),sml)
107 :
108 450 : dagrzt= dvzt*dvzzt/grad%agrt(i)
109 450 : dagrzu= dvzu*dvzzu/grad%agru(i)
110 450 : dagrzd= dvzd*dvzzd/grad%agrd(i)
111 :
112 450 : grad%gggrt(i) = dvzt*dagrzt
113 450 : grad%gggru(i) = dvzu*dagrzu
114 450 : grad%gggrd(i) = dvzd*dagrzd
115 :
116 : c dztadz=d(zeta)/dz,..
117 :
118 450 : dztadz = (dvzu-dvzd)/vlt - (vlu-vld)*dvzt/vlt**2
119 :
120 : c gzgr=grad(zeta)*grad(ro).
121 :
122 450 : grad%gzgr(i) = dztadz*dvzt
123 :
124 : c g2rt: grad(grad(ro))
125 :
126 450 : grad%g2rt(i) = dvzzt
127 450 : grad%g2ru(i) = dvzzu
128 453 : grad%g2rd(i) = dvzzd
129 :
130 :
131 : ENDDO
132 : ENDIF
133 :
134 : END SUBROUTINE mkgz
135 : END MODULE m_mkgz
|
