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