types_xcpot_libxc.F90 12.3 KB
Newer Older
1 2 3 4 5 6 7 8
!--------------------------------------------------------------------------------
! 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.
!--------------------------------------------------------------------------------

!>This module contains the xcpot-type providing an interface to libxc
MODULE m_types_xcpot_libxc
9 10
#ifdef CPP_LIBXC
   USE xc_f03_lib_m
11
#endif
Matthias Redies's avatar
Matthias Redies committed
12 13 14
   USE m_types_xcpot
   USE m_judft
   IMPLICIT NONE
Matthias Redies's avatar
Matthias Redies committed
15 16 17


#ifdef CPP_LIBXC
Matthias Redies's avatar
Matthias Redies committed
18
   PRIVATE :: write_xc_info
Matthias Redies's avatar
Matthias Redies committed
19
#endif
20

Matthias Redies's avatar
Matthias Redies committed
21
   TYPE,EXTENDS(t_xcpot):: t_xcpot_libxc 
22
#ifdef CPP_LIBXC    
Matthias Redies's avatar
Matthias Redies committed
23 24
      TYPE(xc_f03_func_t)      :: vxc_func_x, vxc_func_c
      TYPE(xc_f03_func_t)      :: exc_func_x, exc_func_c
25
#endif
Matthias Redies's avatar
Matthias Redies committed
26 27 28
      INTEGER          :: func_vxc_id_c, func_vxc_id_x !> functionals to be used for potential & density convergence
      INTEGER          :: func_exc_id_c, func_exc_id_x !> functionals to be used in exc- & totale-calculations
      INTEGER          :: jspins
29
   CONTAINS
30 31
      PROCEDURE        :: vxc_is_gga          => xcpot_vxc_is_gga
      PROCEDURE        :: exc_is_gga          => xcpot_exc_is_gga
Matthias Redies's avatar
Matthias Redies committed
32
      PROCEDURE        :: exc_is_MetaGGA      => xcpot_exc_is_MetaGGA
Matthias Redies's avatar
Matthias Redies committed
33 34 35 36 37 38 39 40 41
      PROCEDURE        :: is_hybrid           => xcpot_is_hybrid 
      PROCEDURE        :: get_exchange_weight => xcpot_get_exchange_weight
      PROCEDURE        :: get_vxc             => xcpot_get_vxc
      PROCEDURE        :: get_exc             => xcpot_get_exc
      PROCEDURE,NOPASS :: alloc_gradients     => xcpot_alloc_gradients
      !Not             overloeaded...
      PROCEDURE        :: init                => xcpot_init
   END TYPE t_xcpot_libxc
   PUBLIC t_xcpot_libxc
42 43
CONTAINS

Matthias Redies's avatar
Matthias Redies committed
44 45 46 47 48 49
   SUBROUTINE xcpot_init(xcpot,jspins,vxc_id_x,vxc_id_c,exc_id_x,exc_id_c)
      USE m_judft
      IMPLICIT NONE
   CLASS(t_xcpot_libxc),INTENT(OUT)    :: xcpot
      INTEGER, INTENT(IN)                 :: jspins
      INTEGER, INTENT(IN)                 :: vxc_id_x, vxc_id_c ! potential functional
50
      INTEGER, INTENT(IN)                 :: exc_id_x, exc_id_c ! energy functionals
Matthias Redies's avatar
Matthias Redies committed
51 52
      LOGICAL                             :: same_functionals   ! are vxc and exc equal

53
#ifdef CPP_LIBXC
Matthias Redies's avatar
Matthias Redies committed
54 55 56 57 58
      INTEGER :: err
      xcpot%jspins = jspins

      xcpot%func_vxc_id_x = vxc_id_x
      xcpot%func_vxc_id_c = vxc_id_c
Matthias Redies's avatar
Matthias Redies committed
59
      
60 61
      xcpot%func_exc_id_x = exc_id_x
      xcpot%func_exc_id_c = exc_id_c
Matthias Redies's avatar
Matthias Redies committed
62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78

      if(xcpot%func_vxc_id_x == 0 .or. xcpot%func_exc_id_x == 0 ) then
         CALL judft_error("LibXC exchange- and correlation-function indicies need to be set"&
            ,hint='Try this: ' // ACHAR(10) //&
            '<xcFunctional name="libxc" relativisticCorrections="F">' // ACHAR(10) //& 
            '  <libXC  exchange="1" correlation="1" /> ' // ACHAR(10) //& 
            '</xcFunctional> ')
      endif 

      IF (jspins==1) THEN
         ! potential functionals      
         CALL xc_f03_func_init(xcpot%vxc_func_x, xcpot%func_vxc_id_x, XC_UNPOLARIZED)
         IF (xcpot%func_vxc_id_c>0) CALL xc_f03_func_init(xcpot%vxc_func_c, xcpot%func_vxc_id_c, XC_UNPOLARIZED)

         ! energy functionals
         CALL xc_f03_func_init(xcpot%exc_func_x, xcpot%func_exc_id_x, XC_UNPOLARIZED)
         IF (xcpot%func_exc_id_c>0) CALL xc_f03_func_init(xcpot%exc_func_c, xcpot%func_exc_id_c, XC_UNPOLARIZED)
Matthias Redies's avatar
Matthias Redies committed
79
      ELSE
Matthias Redies's avatar
Matthias Redies committed
80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96
         ! potential functionals  
         CALL xc_f03_func_init(xcpot%vxc_func_x, xcpot%func_vxc_id_x, XC_POLARIZED)
         IF (xcpot%func_vxc_id_c>0) CALL xc_f03_func_init(xcpot%vxc_func_c, xcpot%func_vxc_id_c, XC_POLARIZED)

         !energy functionals
         CALL xc_f03_func_init(xcpot%exc_func_x, xcpot%func_exc_id_x, XC_POLARIZED)
         IF (xcpot%func_exc_id_c>0) CALL xc_f03_func_init(xcpot%exc_func_c, xcpot%func_exc_id_c, XC_POLARIZED)
      END IF

      CALL write_xc_info(xcpot%vxc_func_x)

      IF (xcpot%func_vxc_id_c>0) THEN
         CALL write_xc_info(xcpot%vxc_func_c)
      ELSE
         WRITE(*,*) "No Correlation functional"
      END IF

97 98
      same_functionals =     (xcpot%func_vxc_id_x == xcpot%func_exc_id_x) &
                       .and. (xcpot%func_vxc_id_c == xcpot%func_exc_id_c)
Matthias Redies's avatar
Matthias Redies committed
99 100 101 102 103 104 105
      IF(.not. same_functionals) THEN
         CALL write_xc_info(xcpot%exc_func_x)
         IF (xcpot%func_exc_id_c>0) THEN
            CALL write_xc_info(xcpot%exc_func_c)
         ELSE
            WRITE(*,*) "No Correlation functional for TotalE"
         ENDIF
Matthias Redies's avatar
Matthias Redies committed
106 107
      ELSE
         write (*,*) "Using same functional for VXC and EXC"
Matthias Redies's avatar
Matthias Redies committed
108
      END IF
109
#else
110
      CALL judft_error("You specified a libxc-exchange correlation potential but FLEUR is not linked against libxc", &
Matthias Redies's avatar
Matthias Redies committed
111
         hint="Please recompile FLEUR with libxc support")
112
#endif
113
   END SUBROUTINE xcpot_init
114

115
   LOGICAL FUNCTION xcpot_vxc_is_gga(xcpot)
116
      IMPLICIT NONE
Matthias Redies's avatar
Matthias Redies committed
117
   CLASS(t_xcpot_libxc),INTENT(IN):: xcpot
118
#ifdef CPP_LIBXC
119 120
      TYPE(xc_f03_func_info_t)        :: xc_info

121
      xc_info = xc_f03_func_get_info(xcpot%vxc_func_x)
Matthias Redies's avatar
Matthias Redies committed
122
      xcpot_vxc_is_gga =  ANY([XC_FAMILY_GGA, XC_FAMILY_HYB_GGA]==xc_f03_func_info_get_family(xc_info))
123 124 125 126 127 128 129 130 131 132 133 134
#else
      xcpot_vxc_is_gga=.false.
#endif
   END FUNCTION xcpot_vxc_is_gga

   LOGICAL FUNCTION xcpot_exc_is_gga(xcpot)
      IMPLICIT NONE
   CLASS(t_xcpot_libxc),INTENT(IN):: xcpot
#ifdef CPP_LIBXC
      TYPE(xc_f03_func_info_t)        :: xc_info

      xc_info = xc_f03_func_get_info(xcpot%vxc_func_x)
Matthias Redies's avatar
Matthias Redies committed
135
      xcpot_exc_is_gga =  ANY([XC_FAMILY_GGA, XC_FAMILY_HYB_GGA]==xc_f03_func_info_get_family(xc_info))
136
#else
137
      xcpot_exc_is_gga=.false.
138
#endif
139
   END FUNCTION xcpot_exc_is_gga
140

Matthias Redies's avatar
Matthias Redies committed
141
   LOGICAL FUNCTION xcpot_exc_is_MetaGGA(xcpot)
Matthias Redies's avatar
Matthias Redies committed
142 143
      IMPLICIT NONE
   CLASS(t_xcpot_libxc),INTENT(IN):: xcpot
Matthias Redies's avatar
Matthias Redies committed
144
#ifdef CPP_LIBXC    
Matthias Redies's avatar
Matthias Redies committed
145
      TYPE(xc_f03_func_info_t)        :: xc_info
Matthias Redies's avatar
Matthias Redies committed
146

Matthias Redies's avatar
Matthias Redies committed
147
      xc_info = xc_f03_func_get_info(xcpot%exc_func_x)
Matthias Redies's avatar
Matthias Redies committed
148
      xcpot_exc_is_MetaGGA=ANY([XC_FAMILY_MGGA, XC_FAMILY_HYB_MGGA]==xc_f03_func_info_get_family(xc_info))
149
#else
Matthias Redies's avatar
Matthias Redies committed
150
      xcpot_exc_is_MetaGGA =  .False.
Matthias Redies's avatar
Matthias Redies committed
151
#endif
Matthias Redies's avatar
Matthias Redies committed
152
   END FUNCTION xcpot_exc_is_MetaGGA
Matthias Redies's avatar
Matthias Redies committed
153

154 155
   LOGICAL FUNCTION xcpot_is_hybrid(xcpot)
      IMPLICIT NONE
Matthias Redies's avatar
Matthias Redies committed
156
   CLASS(t_xcpot_libxc),INTENT(IN):: xcpot
157
#ifdef CPP_LIBXC
Matthias Redies's avatar
Matthias Redies committed
158
      TYPE(xc_f03_func_info_t)        :: xc_info
Matthias Redies's avatar
Matthias Redies committed
159

Matthias Redies's avatar
Matthias Redies committed
160 161
      xc_info = xc_f03_func_get_info(xcpot%vxc_func_x)
      xcpot_is_hybrid=ANY([XC_FAMILY_HYB_MGGA, XC_FAMILY_HYB_GGA]==xc_f03_func_info_get_family(xc_info))
162
#else
Matthias Redies's avatar
Matthias Redies committed
163
      xcpot_is_hybrid = .False.
164
#endif
165
   END FUNCTION xcpot_is_hybrid
166

167 168 169
   FUNCTION xcpot_get_exchange_weight(xcpot) RESULT(a_ex)
      USE m_judft
      IMPLICIT NONE
Matthias Redies's avatar
Matthias Redies committed
170
   CLASS(t_xcpot_libxc),INTENT(IN):: xcpot
171

Matthias Redies's avatar
Matthias Redies committed
172
      REAL:: a_ex
173
#ifdef CPP_LIBXC    
Matthias Redies's avatar
Matthias Redies committed
174
      a_ex=xc_f03_hyb_exx_coef(xcpot%vxc_func_x)
175
#endif
Matthias Redies's avatar
Matthias Redies committed
176 177 178 179 180 181 182 183 184 185 186 187 188
   END FUNCTION xcpot_get_exchange_weight

   !***********************************************************************
   SUBROUTINE xcpot_get_vxc(xcpot,jspins,rh, vxc,vx, grad)
      !***********************************************************************
      IMPLICIT NONE
   CLASS(t_xcpot_libxc),INTENT(IN) :: xcpot
      INTEGER, INTENT (IN)     :: jspins
      REAL,INTENT (IN)         :: rh(:,:)   !Dimensions here
      REAL, INTENT (OUT)       :: vx (:,:)  !points,spin
      REAL, INTENT (OUT  )     :: vxc(:,:)  !
      ! optional arguments for GGA
      TYPE(t_gradients),OPTIONAL,INTENT(INOUT)::grad
189
#ifdef CPP_LIBXC    
Matthias Redies's avatar
Matthias Redies committed
190 191 192 193
      REAL,ALLOCATABLE::vxc_tmp(:,:),vx_tmp(:,:),vsigma(:,:)
      !libxc uses the spin as a first index, hence we have to transpose....
      ALLOCATE(vxc_tmp(SIZE(vxc,2),SIZE(vxc,1)));vxc_tmp=0.0
      ALLOCATE(vx_tmp(SIZE(vx,2),SIZE(vx,1)));vx_tmp=0.0
194
      IF (xcpot%vxc_is_gga()) THEN
Matthias Redies's avatar
Matthias Redies committed
195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215
         IF (.NOT.PRESENT(grad)) CALL judft_error("Bug: You called get_vxc for a GGA potential without providing derivatives")
         ALLOCATE(vsigma,mold=grad%vsigma)
         !where(abs(grad%sigma)<1E-9) grad%sigma=1E-9
         CALL xc_f03_gga_vxc(xcpot%vxc_func_x, SIZE(rh,1), TRANSPOSE(rh),grad%sigma,vx_tmp,vsigma)
         IF (xcpot%func_vxc_id_c>0) THEN
            CALL xc_f03_gga_vxc(xcpot%vxc_func_c, SIZE(rh,1), TRANSPOSE(rh),grad%sigma,vxc_tmp,grad%vsigma)
            grad%vsigma=grad%vsigma+vsigma
            vxc_tmp=vxc_tmp+vx_tmp
         ELSE     
            vxc_tmp=vx_tmp
         ENDIF
      ELSE  !LDA potentials
         CALL xc_f03_lda_vxc(xcpot%vxc_func_x, SIZE(rh,1), TRANSPOSE(rh), vx_tmp)
         IF (xcpot%func_vxc_id_c>0) THEN
            CALL xc_f03_lda_vxc(xcpot%vxc_func_c, SIZE(rh,1), TRANSPOSE(rh), vxc_tmp)
            vxc_tmp=vxc_tmp+vx_tmp
         ENDIF
      ENDIF
      vx=TRANSPOSE(vx_tmp) 
      vxc=TRANSPOSE(vxc_tmp) 

216
#endif
Matthias Redies's avatar
Matthias Redies committed
217 218 219 220 221 222 223 224 225
   END SUBROUTINE xcpot_get_vxc


   SUBROUTINE xcpot_get_exc(xcpot,jspins,rh,exc,grad)
      IMPLICIT NONE
   CLASS(t_xcpot_libxc),INTENT(IN) :: xcpot
      INTEGER, INTENT (IN)     :: jspins
      REAL,INTENT (IN) :: rh(:,:)  !points,spin
      REAL, INTENT (OUT) :: exc(:) !points
226
      ! optional arguments for GGA
Matthias Redies's avatar
Matthias Redies committed
227 228 229
      TYPE(t_gradients),OPTIONAL,INTENT(IN)::grad

      REAL  :: excc(SIZE(exc))
230
#ifdef CPP_LIBXC
231
      IF (xcpot%exc_is_gga()) THEN
Matthias Redies's avatar
Matthias Redies committed
232 233
         IF (.NOT.PRESENT(grad)) CALL judft_error("Bug: You called get_vxc for a GGA potential without providing derivatives")
         CALL xc_f03_gga_exc(xcpot%exc_func_x, SIZE(rh,1), TRANSPOSE(rh),grad%sigma,exc)
234
         IF (xcpot%func_exc_id_c>0) THEN
Matthias Redies's avatar
Matthias Redies committed
235 236 237 238 239
            CALL xc_f03_gga_exc(xcpot%exc_func_c, SIZE(rh,1), TRANSPOSE(rh),grad%sigma,excc)
            exc=exc+excc
         END IF
      ELSE  !LDA potentials
         CALL xc_f03_lda_exc(xcpot%exc_func_x, SIZE(rh,1), TRANSPOSE(rh), exc)
240
         IF (xcpot%func_exc_id_c>0) THEN
Matthias Redies's avatar
Matthias Redies committed
241 242 243 244 245
            CALL xc_f03_lda_exc(xcpot%exc_func_c, SIZE(rh,1), TRANSPOSE(rh), excc)
            exc=exc+excc
         END IF
      ENDIF

246 247
#endif
   END SUBROUTINE xcpot_get_exc
248

249 250 251 252 253 254 255 256 257
   SUBROUTINE xcpot_alloc_gradients(ngrid,jspins,grad)
      INTEGER, INTENT (IN)         :: jspins,ngrid
      TYPE(t_gradients),INTENT(INOUT):: grad
      !For libxc we only need the sigma array...
      IF (ALLOCATED(grad%sigma)) DEALLOCATE(grad%sigma,grad%gr,grad%laplace,grad%vsigma)
      ALLOCATE(grad%sigma(MERGE(1,3,jspins==1),ngrid))
      ALLOCATE(grad%gr(3,ngrid,jspins))
      ALLOCATE(grad%laplace(ngrid,jspins))
      ALLOCATE(grad%vsigma(MERGE(1,3,jspins==1),ngrid))
258

259
   END SUBROUTINE xcpot_alloc_gradients
260

Matthias Redies's avatar
Matthias Redies committed
261
#ifdef CPP_LIBXC  
Matthias Redies's avatar
Matthias Redies committed
262 263 264 265 266 267 268
   SUBROUTINE write_xc_info(xc_func, is_E_func)
      IMPLICIT NONE
      LOGICAL,INTENT(IN),OPTIONAL         :: is_E_func
      INTEGER                             :: i
      CHARACTER(len=120)                  :: kind, family
      LOGICAL                             :: is_energy_func

269 270 271
      TYPE(xc_f03_func_t),INTENT(IN)      :: xc_func
      TYPE(xc_f03_func_info_t)            :: xc_info

Matthias Redies's avatar
Matthias Redies committed
272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313
      xc_info = xc_f03_func_get_info(xc_func)
      is_energy_func = merge(is_E_func, .False., PRESENT(is_E_func))

      SELECT CASE(xc_f03_func_info_get_kind(xc_info))
      CASE (XC_EXCHANGE)
         WRITE(kind, '(a)') 'an exchange functional'
      CASE (XC_CORRELATION)
         WRITE(kind, '(a)') 'a correlation functional'
      CASE (XC_EXCHANGE_CORRELATION)
         WRITE(kind, '(a)') 'an exchange-correlation functional'
      CASE (XC_KINETIC)
         WRITE(kind, '(a)') 'a kinetic energy functional'
      CASE default
         WRITE(kind, '(a)') 'of unknown kind'
      END SELECT
      SELECT CASE (xc_f03_func_info_get_family(xc_info))
      CASE (XC_FAMILY_LDA);
         WRITE(family,'(a)') "LDA"
      CASE (XC_FAMILY_GGA);
         WRITE(family,'(a)') "GGA"
      CASE (XC_FAMILY_HYB_GGA);
         WRITE(family,'(a)') "Hybrid GGA"
      CASE (XC_FAMILY_MGGA);
         WRITE(family,'(a)') "MGGA"
      CASE (XC_FAMILY_HYB_MGGA);
         WRITE(family,'(a)') "Hybrid MGGA"
      CASE default;
         WRITE(family,'(a)') "unknown"
      END SELECT

      IF(.not. is_energy_func) THEN
         WRITE(*,'("The functional ''", a, "'' is ", a, ", it belongs to the ''", a, "'' family and is defined in the reference(s):")') &
            TRIM(xc_f03_func_info_get_name(xc_info)), TRIM(kind), TRIM(family)
      ELSE
         WRITE(*,'("The functional used for TotalE ''", a, "'' is ", a, ", it belongs to the ''", a, "'' family and is defined in the reference(s):")') &
            TRIM(xc_f03_func_info_get_name(xc_info)), TRIM(kind), TRIM(family)
      ENDIF

      i = 0
      DO WHILE(i >= 0)
         WRITE(*, '(a,i1,2a)') '[', i+1, '] ', TRIM(xc_f03_func_reference_get_ref(xc_f03_func_info_get_references(xc_info, i)))
      END DO
314
   END SUBROUTINE write_xc_info
Matthias Redies's avatar
Matthias Redies committed
315
#endif
316

317
END MODULE m_types_xcpot_libxc