Line data Source code
1 : !--------------------------------------------------------------------------------
2 : ! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
3 : ! This file is part of FLEUR and available as free software under the conditions
4 : ! of the MIT license as expressed in the LICENSE file in more detail.
5 : !--------------------------------------------------------------------------------
6 :
7 : MODULE m_xcxal
8 : use m_juDFT
9 : !-----------------------------------------------------------------------
10 : ! Called in case of icorr=0 : nonspin-polarized exchange-correlation
11 : ! potential with the X-alpha method
12 :
13 : ! We calculate V_x = alpha V_{x,Slater} where alpha = 2/3
14 : ! i.e. the output corresponds to the Gaspar result for the
15 : ! local approximation to the Hartree-Fock method.
16 :
17 : ! krla=1: Relativistic correction of exchange energy and potential
18 : ! related to Dirac kinetic energy, according to:
19 : ! A.H. MacDonald and S.H. Vosko, J. Phys. C12, 2977 (1979)
20 :
21 : ! be careful: calculation in rydberg!
22 :
23 : ! vxcxal calculates the XC-potential and
24 : ! excxal calculates the XC-energy
25 :
26 : ! based on a subroutine by S. Bluegel; r.pentcheva 22.01.96
27 : !-----------------------------------------------------------------------
28 :
29 : USE m_constants, ONLY : pi_const
30 : USE m_relcor
31 : IMPLICIT NONE
32 :
33 : REAL, PARAMETER, PRIVATE :: one = 1.0 , three = 3.0 , four = 4.0
34 : REAL, PARAMETER, PRIVATE :: thrd = one/three , d_15 = 1.e-15
35 : REAL, PARAMETER, PRIVATE :: cvx = 1.221774115422 ! 2 * ( 3/(2*pi) )^(2/3)
36 :
37 : CONTAINS
38 : !************************************************************************
39 0 : SUBROUTINE vxcxal( &
40 : krla,jspins, &
41 0 : mgrid,ngrid,rh, &
42 0 : vx,vxc)
43 : !************************************************************************
44 :
45 : ! .. Scalar Arguments ..
46 : INTEGER, INTENT (IN) :: jspins
47 : INTEGER, INTENT (IN) :: krla ! switch for rel. corr.
48 : INTEGER, INTENT (IN) :: mgrid,ngrid ! mesh points
49 :
50 : ! .. Array Arguments ..
51 : REAL, INTENT (IN) :: rh(mgrid,jspins) ! charge density
52 : REAL, INTENT (OUT) :: vx (mgrid,jspins)
53 : REAL, INTENT (OUT) :: vxc(mgrid,jspins) ! xc potential
54 :
55 : ! .. Local Arrays ..
56 0 : REAL, ALLOCATABLE :: psi(:) ! relativistic exchange potential corr.
57 :
58 : REAL :: rs, rho, thfpi
59 : INTEGER :: ir
60 :
61 : !-----s Intrinsic Functions
62 : INTRINSIC alog,max
63 :
64 0 : thfpi = three / ( four * pi_const )
65 :
66 : !-----> evaluate relativistic corrections for exchange
67 :
68 0 : ALLOCATE ( psi(ngrid) )
69 : CALL relcor( &
70 : mgrid,ngrid,jspins,krla, .TRUE. ,rh, &
71 0 : psi)
72 :
73 0 : DO ir = 1,ngrid
74 0 : IF (jspins == 1) THEN
75 0 : rho = max(d_15,rh(ir,1))
76 : ELSE
77 0 : rho = max(d_15,rh(ir,1)) + max(d_15,rh(ir,jspins))
78 : ENDIF
79 : rs = (thfpi/rho)**thrd
80 : ! exc: exchange energy = -0.9163306/rs
81 0 : vxc(ir,1) = - psi(ir)*cvx/rs ! paramagnetic exchange potential = 4/3 exc
82 :
83 0 : vx (ir,1) = - psi(ir)*cvx/rs
84 :
85 : ENDDO
86 0 : IF ( jspins == 2 ) THEN ! spinpolarized calculation
87 0 : DO ir = 1,ngrid
88 0 : vxc(ir,jspins) = vxc(ir,1)
89 :
90 0 : vx(ir,jspins) = vx(ir,1)
91 : ENDDO
92 : ENDIF
93 :
94 0 : DEALLOCATE (psi)
95 0 : RETURN
96 :
97 : END SUBROUTINE vxcxal
98 : !***********************************************************************
99 0 : SUBROUTINE excxal( &
100 : iofile,krla,jspins, &
101 0 : mgrid,ngrid,rh, &
102 0 : exc)
103 : !***********************************************************************
104 :
105 : ! .. Scalar Arguments ..
106 : INTEGER, INTENT (IN) :: krla, jspins
107 : INTEGER, INTENT (IN) :: iofile ! file number for read and write
108 : INTEGER, INTENT (IN) :: mgrid,ngrid ! mesh points
109 :
110 : ! .. Array Arguments ..
111 : REAL, INTENT (IN) :: rh(mgrid,jspins) ! charge density
112 : REAL, INTENT (OUT) :: exc(mgrid) ! xc energy
113 :
114 : ! .. Local Scalars ..
115 : REAL :: cex
116 : REAL :: rs, rho, thfpi
117 : INTEGER :: ir
118 :
119 : ! .. Local Arrays ..
120 0 : REAL, ALLOCATABLE :: phi(:) ! relativistic exchange energy correct.
121 :
122 : !-----> Intrinsic Functions
123 : INTRINSIC alog,max
124 :
125 0 : cex = three * cvx / four
126 0 : thfpi = three / ( four * pi_const )
127 :
128 0 : ALLOCATE ( phi(ngrid) )
129 : CALL relcor( &
130 : mgrid,ngrid,jspins,krla, .FALSE. ,rh, &
131 0 : phi)
132 :
133 0 : DO ir = 1,ngrid
134 0 : IF (jspins == 1) THEN
135 0 : rho = max(d_15,rh(ir,1))
136 : ELSE
137 0 : rho = max(d_15,rh(ir,1)) + max(d_15,rh(ir,jspins))
138 : ENDIF
139 : rs = (thfpi/rho)**thrd
140 : ! exc: exchange energy = -0.9163306/rs
141 0 : exc(ir) = - phi(ir)*cex / rs
142 : ENDDO
143 0 : IF ( jspins == 2 ) WRITE (iofile,*) &
144 0 : 'WARNING: X-alpha XC & spin-polarized calculation'
145 0 : IF ( jspins > 2 ) THEN
146 0 : WRITE (iofile,'('' or error in jspins, jspins ='',i2)') jspins
147 0 : CALL juDFT_error("excxal",calledby="xcxal")
148 : END IF
149 :
150 0 : DEALLOCATE (phi)
151 0 : RETURN
152 :
153 : END SUBROUTINE excxal
154 : END MODULE m_xcxal
|
