1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 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
|
!--------------------------------------------------------------------------------------------------!
! CP2K: A general program to perform molecular dynamics simulations !
! Copyright (C) 2000 - 2018 CP2K developers group !
!--------------------------------------------------------------------------------------------------!
! **************************************************************************************************
!> \brief Set disperson types for DFT calculations
!> \author JGH (04.2014)
! **************************************************************************************************
MODULE qs_dispersion_utils
USE atomic_kind_types, ONLY: atomic_kind_type,&
get_atomic_kind
USE cp_log_handling, ONLY: cp_get_default_logger,&
cp_logger_type
USE cp_output_handling, ONLY: cp_print_key_finished_output,&
cp_print_key_unit_nr
USE input_constants, ONLY: vdw_nl_DRSLL,&
vdw_nl_LMKLL,&
vdw_nl_RVV10,&
vdw_pairpot_dftd2,&
vdw_pairpot_dftd3,&
vdw_pairpot_dftd3bj,&
xc_vdw_fun_nonloc,&
xc_vdw_fun_pairpot
USE input_section_types, ONLY: section_vals_get_subs_vals,&
section_vals_type,&
section_vals_val_get
USE kinds, ONLY: dp
USE physcon, ONLY: bohr,&
kjmol
USE qs_dispersion_pairpot, ONLY: qs_scaling_dftd3,&
qs_scaling_dftd3bj,&
qs_scaling_init
USE qs_dispersion_types, ONLY: qs_atom_dispersion_type,&
qs_dispersion_type
USE qs_environment_types, ONLY: get_qs_env,&
qs_environment_type
USE qs_kind_types, ONLY: get_qs_kind,&
qs_kind_type
#include "./base/base_uses.f90"
IMPLICIT NONE
PRIVATE
CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_dispersion_utils'
PUBLIC :: qs_dispersion_env_set, qs_write_dispersion
! **************************************************************************************************
CONTAINS
! **************************************************************************************************
!> \brief ...
!> \param dispersion_env ...
!> \param xc_section ...
! **************************************************************************************************
SUBROUTINE qs_dispersion_env_set(dispersion_env, xc_section)
TYPE(qs_dispersion_type), POINTER :: dispersion_env
TYPE(section_vals_type), POINTER :: xc_section
CHARACTER(len=*), PARAMETER :: routineN = 'qs_dispersion_env_set', &
routineP = moduleN//':'//routineN
LOGICAL :: exfun, explicit
REAL(dp), POINTER :: params(:), scal(:)
TYPE(section_vals_type), POINTER :: nl_section, pp_section, vdw_section, &
xc_fun_section
CPASSERT(ASSOCIATED(dispersion_env))
! set general defaults
dispersion_env%doabc = .FALSE.
dispersion_env%c9cnst = .FALSE.
dispersion_env%lrc = .FALSE.
dispersion_env%srb = .FALSE.
dispersion_env%verbose = .FALSE.
NULLIFY (dispersion_env%c6ab, dispersion_env%maxci, dispersion_env%r0ab, dispersion_env%rcov, &
dispersion_env%r2r4, dispersion_env%cn, dispersion_env%cnkind, dispersion_env%cnlist)
NULLIFY (dispersion_env%q_mesh, dispersion_env%kernel, dispersion_env%d2phi_dk2, &
dispersion_env%d2y_dx2)
NULLIFY (dispersion_env%sab_vdw, dispersion_env%sab_cn)
NULLIFY (dispersion_env%dftd_section)
NULLIFY (vdw_section, xc_fun_section)
vdw_section => section_vals_get_subs_vals(xc_section, "vdw_potential")
xc_fun_section => section_vals_get_subs_vals(xc_section, "XC_FUNCTIONAL")
CALL section_vals_val_get(vdw_section, "POTENTIAL_TYPE", i_val=dispersion_env%type)
IF (dispersion_env%type == xc_vdw_fun_pairpot) THEN
NULLIFY (pp_section)
pp_section => section_vals_get_subs_vals(vdw_section, "PAIR_POTENTIAL")
CALL section_vals_val_get(pp_section, "VERBOSE_OUTPUT", l_val=dispersion_env%verbose)
CALL section_vals_val_get(pp_section, "TYPE", i_val=dispersion_env%pp_type)
IF (dispersion_env%pp_type == vdw_pairpot_dftd2) THEN
! functional parameters for Grimme D2 type
CALL section_vals_val_get(pp_section, "EXP_PRE", r_val=dispersion_env%exp_pre)
CALL section_vals_val_get(pp_section, "SCALING", explicit=explicit)
IF (.NOT. explicit) THEN
CALL section_vals_val_get(pp_section, "REFERENCE_FUNCTIONAL", explicit=exfun)
CPASSERT(exfun)
CALL qs_scaling_init(dispersion_env%scaling, vdw_section)
ELSE
CALL section_vals_val_get(pp_section, "SCALING", r_val=dispersion_env%scaling)
END IF
ELSE
dispersion_env%exp_pre = 0._dp
dispersion_env%scaling = 0._dp
END IF
IF (dispersion_env%pp_type == vdw_pairpot_dftd3 .OR. &
dispersion_env%pp_type == vdw_pairpot_dftd3bj) THEN
! functional parameters for Grimme DFT-D3 type
CALL section_vals_val_get(pp_section, "EPS_CN", r_val=dispersion_env%eps_cn)
CALL section_vals_val_get(pp_section, "CALCULATE_C9_TERM", l_val=dispersion_env%doabc)
CALL section_vals_val_get(pp_section, "REFERENCE_C9_TERM", l_val=dispersion_env%c9cnst)
CALL section_vals_val_get(pp_section, "LONG_RANGE_CORRECTION", l_val=dispersion_env%lrc)
CALL section_vals_val_get(pp_section, "SHORT_RANGE_CORRECTION", l_val=dispersion_env%srb)
! KG corrections
CALL section_vals_val_get(pp_section, "MOLECULE_CORRECTION", l_val=dispersion_env%domol)
CALL section_vals_val_get(pp_section, "MOLECULE_CORRECTION_C8", r_val=dispersion_env%kgc8)
IF (dispersion_env%pp_type == vdw_pairpot_dftd3) THEN
CALL section_vals_val_get(pp_section, "D3_SCALING", explicit=explicit)
ELSE IF (dispersion_env%pp_type == vdw_pairpot_dftd3bj) THEN
CALL section_vals_val_get(pp_section, "D3BJ_SCALING", explicit=explicit)
END IF
IF (.NOT. explicit) THEN
CALL section_vals_val_get(pp_section, "REFERENCE_FUNCTIONAL", explicit=exfun)
CPASSERT(exfun)
IF (dispersion_env%pp_type == vdw_pairpot_dftd3) THEN
CALL qs_scaling_dftd3(dispersion_env%s6, dispersion_env%sr6, dispersion_env%s8, vdw_section)
ELSE IF (dispersion_env%pp_type == vdw_pairpot_dftd3bj) THEN
CALL qs_scaling_dftd3bj(dispersion_env%s6, dispersion_env%a1, dispersion_env%s8, &
dispersion_env%a2, vdw_section)
END IF
ELSE
IF (dispersion_env%pp_type == vdw_pairpot_dftd3) THEN
! zero damping
CALL section_vals_val_get(pp_section, "D3_SCALING", r_vals=scal)
dispersion_env%s6 = scal(1)
dispersion_env%sr6 = scal(2)
dispersion_env%s8 = scal(3)
dispersion_env%a1 = 0.0_dp
dispersion_env%a2 = 0.0_dp
ELSE IF (dispersion_env%pp_type == vdw_pairpot_dftd3bj) THEN
! BJ damping
CALL section_vals_val_get(pp_section, "D3BJ_SCALING", r_vals=scal)
dispersion_env%s6 = scal(1)
dispersion_env%a1 = scal(2)
dispersion_env%s8 = scal(3)
dispersion_env%a2 = scal(4)
dispersion_env%sr6 = 0.0_dp
END IF
END IF
ELSE
dispersion_env%s6 = 0._dp
dispersion_env%sr6 = 0._dp
dispersion_env%s8 = 0._dp
dispersion_env%a1 = 0._dp
dispersion_env%a2 = 0._dp
dispersion_env%eps_cn = 0._dp
END IF
CALL section_vals_val_get(pp_section, "R_CUTOFF", r_val=dispersion_env%rc_disp)
CALL section_vals_val_get(pp_section, "PARAMETER_FILE_NAME", &
c_val=dispersion_env%parameter_file_name)
! set DFTD section for output handling
dispersion_env%dftd_section => pp_section
ELSE IF (dispersion_env%type == xc_vdw_fun_nonloc) THEN
NULLIFY (nl_section)
nl_section => section_vals_get_subs_vals(vdw_section, "NON_LOCAL")
CALL section_vals_val_get(nl_section, "VERBOSE_OUTPUT", l_val=dispersion_env%verbose)
CALL section_vals_val_get(nl_section, "KERNEL_FILE_NAME", &
c_val=dispersion_env%kernel_file_name)
CALL section_vals_val_get(nl_section, "TYPE", i_val=dispersion_env%nl_type)
CALL section_vals_val_get(nl_section, "CUTOFF", r_val=dispersion_env%pw_cutoff)
dispersion_env%pw_cutoff = 0.5_dp*dispersion_env%pw_cutoff
CALL section_vals_val_get(nl_section, "PARAMETERS", r_vals=params)
dispersion_env%b_value = params(1)
dispersion_env%c_value = params(2)
END IF
END SUBROUTINE qs_dispersion_env_set
! **************************************************************************************************
! **************************************************************************************************
!> \brief ...
!> \param qs_env ...
!> \param dispersion_env ...
!> \param ounit ...
! **************************************************************************************************
SUBROUTINE qs_write_dispersion(qs_env, dispersion_env, ounit)
TYPE(qs_environment_type), POINTER :: qs_env
TYPE(qs_dispersion_type), POINTER :: dispersion_env
INTEGER, INTENT(in), OPTIONAL :: ounit
CHARACTER(LEN=*), PARAMETER :: routineN = 'qs_write_dispersion', &
routineP = moduleN//':'//routineN
CHARACTER(LEN=2) :: symbol
INTEGER :: i, ikind, nkind, output_unit
TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
TYPE(cp_logger_type), POINTER :: logger
TYPE(qs_atom_dispersion_type), POINTER :: disp
TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
TYPE(section_vals_type), POINTER :: dft_section
IF (PRESENT(ounit)) THEN
output_unit = ounit
ELSE
NULLIFY (logger)
logger => cp_get_default_logger()
dft_section => section_vals_get_subs_vals(qs_env%input, "DFT")
output_unit = cp_print_key_unit_nr(logger, dft_section, &
"PRINT%DFT_CONTROL_PARAMETERS", extension=".Log")
END IF
IF (output_unit > 0) THEN
! vdW type specific output
IF (dispersion_env%type == xc_vdw_fun_pairpot) THEN
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T67,'Pair Potential')")
! Pair potentials
IF (dispersion_env%pp_type == vdw_pairpot_dftd2) THEN
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T35,'Potential Form: S. Grimme, JCC 27: 1787 (2006)')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T35,'Cutoff Radius [Bohr]:',T73,F8.2)") dispersion_env%rc_disp
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T35,'Scaling Factor:',T73,F8.4)") dispersion_env%scaling
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T35,'Exp Prefactor for Damping:',T73,F8.1)") dispersion_env%exp_pre
CALL get_qs_env(qs_env, atomic_kind_set=atomic_kind_set, qs_kind_set=qs_kind_set)
nkind = SIZE(atomic_kind_set)
DO ikind = 1, nkind
CALL get_atomic_kind(atomic_kind_set(ikind), element_symbol=symbol)
CALL get_qs_kind(qs_kind_set(ikind), dispersion=disp)
IF (disp%defined) THEN
WRITE (output_unit, fmt="(' vdW PARAMETER| ',T18,'Atom=',A2,"// &
"T28,'C6[J*nm^6*mol^-1]=',F8.4,T63,'r(vdW)[A]=',F8.4)") &
symbol, disp%c6/(1000._dp*bohr**6/kjmol), disp%vdw_radii/bohr
ELSE
WRITE (output_unit, fmt="(' vdW PARAMETER| ',T20,'Atom=',A2,T70,'not defined')")
END IF
END DO
ELSE IF (dispersion_env%pp_type == vdw_pairpot_dftd3) THEN
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'DFT-D3 (Version 3.1)')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Potential Form: S. Grimme et al, JCP 132: 154104 (2010)')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Zero Damping')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Cutoff Radius [Bohr]:',T73,F8.2)") dispersion_env%rc_disp
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'s6 Scaling Factor:',T73,F8.4)") dispersion_env%s6
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'sr6 Scaling Factor:',T73,F8.4)") dispersion_env%sr6
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'s8 Scaling Factor:',T73,F8.4)") dispersion_env%s8
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Cutoff for CN calculation:',T69,E12.4)") dispersion_env%eps_cn
ELSE IF (dispersion_env%pp_type == vdw_pairpot_dftd3bj) THEN
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'DFT-D3 (Version 3.1)')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Potential Form: S. Grimme et al, JCP 132: 154104 (2010)')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'BJ Damping: S. Grimme et al, JCC 32: 1456 (2011)')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Cutoff Radius [Bohr]:',T73,F8.2)") dispersion_env%rc_disp
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'s6 Scaling Factor:',T73,F8.4)") dispersion_env%s6
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'a1 Damping Factor:',T73,F8.4)") dispersion_env%a1
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'s8 Scaling Factor:',T73,F8.4)") dispersion_env%s8
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'a2 Damping Factor:',T73,F8.4)") dispersion_env%a2
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T26,'Cutoff for CN calculation:',T69,E12.4)") dispersion_env%eps_cn
END IF
ELSE IF (dispersion_env%type == xc_vdw_fun_nonloc) THEN
WRITE (output_unit, fmt="(' vdW POTENTIAL| ',T61,'Non-local Functional')")
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ','Implementation: G. Roman-Perez, J. Soler, PRL 103: 096102 (2009)')")
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ',T38,' T. Thonhauser et al, PRB 76: 125112 (2007)')")
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ',T22,' R. Sabatini et al, J.Phys:Condens Matter 24: 424209 (2012)')")
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ',T16,' Based on QE implementation by Brian Kolb, Timo Thonhauser (2009)')")
SELECT CASE (dispersion_env%nl_type)
CASE DEFAULT
! unknown functional
CPABORT("")
CASE (vdw_nl_DRSLL)
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ','DRSLL Functional: M. Dion et al, PRL 92: 246401 (2004)')")
CASE (vdw_nl_LMKLL)
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ','LMKLL Functional: K. Lee et al, PRB 82: 081101 (2010)')")
CASE (vdw_nl_RVV10)
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ','RVV10 Functional: R. Sabatini et al, PRB 87: 041108(R) (2013)')")
END SELECT
IF (dispersion_env%verbose) THEN
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ',' Carrying out vdW-DF run using the following parameters:')")
WRITE (output_unit, fmt="(' vdW POTENTIAL| ','Nqs =',I8,' Nr_points =',I8,' r_max =',F10.3)") &
dispersion_env%nqs, dispersion_env%nr_points, dispersion_env%r_max
WRITE (output_unit, fmt="(' vdW POTENTIAL| ','q_mesh =')")
WRITE (output_unit, fmt="(8X,4F18.8)") (dispersion_env%q_mesh(i), i=1, dispersion_env%nqs)
WRITE (output_unit, &
fmt="(' vdW POTENTIAL| ','Density cutoff for comvolution:',T68,F10.1,' Ry')") &
2.0_dp*dispersion_env%pw_cutoff
END IF
END IF
END IF
IF (.NOT. PRESENT(ounit)) THEN
CALL cp_print_key_finished_output(output_unit, logger, dft_section, &
"PRINT%DFT_CONTROL_PARAMETERS")
END IF
END SUBROUTINE qs_write_dispersion
! **************************************************************************************************
END MODULE qs_dispersion_utils
|