Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief module that contains the definitions of the scf types
10 : !> \par History
11 : !> 02.2003 created [fawzi]
12 : !> \author fawzi
13 : ! **************************************************************************************************
14 : MODULE qs_density_mixing_types
15 : #if defined(__TBLITE)
16 : USE tblite_scf, ONLY: mixer_type
17 : #endif
18 : USE ao_util, ONLY: exp_radius
19 : USE input_constants, ONLY: broy_mix, &
20 : direct_p_mix, &
21 : gaussian, &
22 : kerker_mix, &
23 : modified_broy_mix, &
24 : multisec_mix, &
25 : new_pulay_mix, &
26 : no_mix, &
27 : pulay_mix
28 : USE input_keyword_types, ONLY: keyword_create, &
29 : keyword_release, &
30 : keyword_type
31 : USE input_section_types, ONLY: section_add_keyword, &
32 : section_create, &
33 : section_type, &
34 : section_vals_type, &
35 : section_vals_val_get
36 : USE input_val_types, ONLY: real_t
37 : USE kinds, ONLY: default_string_length, &
38 : dp
39 : USE qs_rho_atom_types, ONLY: rho_atom_coeff
40 : USE string_utilities, ONLY: s2a
41 : #include "./base/base_uses.f90"
42 :
43 : IMPLICIT NONE
44 : PRIVATE
45 :
46 : LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
47 :
48 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_density_mixing_types'
49 :
50 : INTEGER, PARAMETER, PUBLIC :: no_mixing_nr = 0, direct_mixing_nr = 1, &
51 : gspace_mixing_nr = 2, pulay_mixing_nr = 3, &
52 : broyden_mixing_nr = 4, &
53 : modified_broyden_mixing_nr = 5, &
54 : multisecant_mixing_nr = 6, &
55 : new_pulay_mixing_nr = 7
56 : PUBLIC :: cp_1d_z_p_type, mixing_storage_create, mixing_storage_type, mixing_storage_release, create_mixing_section
57 :
58 : TYPE cp_1d_z_p_type
59 : COMPLEX(dp), DIMENSION(:), POINTER :: cc => NULL()
60 : END TYPE cp_1d_z_p_type
61 :
62 : TYPE mixing_storage_type
63 : ! Replicated density coordinates are appended on world rank zero only.
64 : INTEGER :: naux = 0
65 : COMPLEX(KIND=dp), ALLOCATABLE :: aux_buffer(:, :)
66 : INTEGER :: ig_max = -1, ncall = -1, ncall_p(2) = -1, nbuffer = -1, n_simple_mix = -1, &
67 : nskip_mixing = -1, p_metric_method = -1
68 : INTEGER, POINTER, DIMENSION(:) :: ig_global_index => NULL()
69 : LOGICAL :: gmix_p = .FALSE.
70 : LOGICAL, POINTER, DIMENSION(:) :: paw => NULL()
71 : CHARACTER(len=15) :: iter_method = ""
72 : REAL(KIND=dp) :: alpha = -1.0_dp, bconst = -1.0_dp, beta = -1.0_dp, broy_w0 = -1.0_dp, &
73 : max_g2 = -1.0_dp, max_gvec_exp = -1.0_dp, pulay_alpha = -1.0_dp, &
74 : pulay_beta = -1.0_dp, r_step = -1.0_dp, reg_par = -1.0_dp, &
75 : sigma_max = -1.0_dp, wc = -1.0_dp, wmax = -1.0_dp, qm = -1.0_dp, &
76 : qk = -1.0_dp, qkappa = -1.0_dp
77 : ! Spin-channel-specific mixing parameters (for nspin==2, ispin=2 is the magnetization channel)
78 : REAL(KIND=dp) :: alpha_mag = -1.0_dp, beta_mag = -1.0_dp
79 : REAL(KIND=dp), DIMENSION(:), POINTER :: p_metric => NULL()
80 : REAL(KIND=dp), DIMENSION(:), POINTER :: kerker_factor => NULL()
81 : REAL(KIND=dp), DIMENSION(:), POINTER :: kerker_factor_mag => NULL()
82 : REAL(KIND=dp), DIMENSION(:), POINTER :: special_metric => NULL()
83 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: weight => NULL()
84 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: norm_res_buffer => NULL()
85 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: fmat => NULL(), gmat => NULL(), pulay_matrix => NULL(), smat => NULL()
86 : !
87 : INTEGER :: nat_local = -1, max_shell = -1
88 : INTEGER :: tb_scc_mixer_memory = 0, tb_scc_mixer_natom = 0, &
89 : tb_scc_mixer_ns = 0, tb_scc_mixer_step = 0
90 : REAL(KIND=dp) :: tb_scc_mixer_error = 0.0_dp
91 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: acharge => NULL()
92 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: dacharge => NULL()
93 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: dfbroy => NULL()
94 : REAL(KIND=dp), DIMENSION(:, :, :), POINTER :: ubroy => NULL()
95 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: abroy => NULL()
96 : REAL(KIND=dp), DIMENSION(:), POINTER :: wbroy => NULL()
97 : INTEGER, DIMENSION(:), POINTER :: atlist => NULL()
98 : !
99 : TYPE(cp_1d_z_p_type), DIMENSION(:), POINTER :: last_res => NULL(), rhoin => NULL(), rhoin_old => NULL()
100 : TYPE(cp_1d_z_p_type), DIMENSION(:, :), POINTER :: delta_res => NULL(), u_vec => NULL(), z_vec => NULL()
101 : TYPE(cp_1d_z_p_type), DIMENSION(:, :), POINTER :: drho_buffer => NULL(), rhoin_buffer => NULL(), res_buffer => NULL()
102 : !
103 : TYPE(rho_atom_coeff), DIMENSION(:, :), POINTER :: cpc_h_lastres => NULL(), cpc_s_lastres => NULL()
104 : TYPE(rho_atom_coeff), DIMENSION(:, :), POINTER :: cpc_h_in => NULL(), cpc_s_in => NULL()
105 : TYPE(rho_atom_coeff), DIMENSION(:, :), POINTER :: cpc_h_old => NULL(), cpc_s_old => NULL()
106 : TYPE(rho_atom_coeff), DIMENSION(:, :, :), POINTER :: cpc_h_in_buffer => NULL(), cpc_s_in_buffer => NULL()
107 : TYPE(rho_atom_coeff), DIMENSION(:, :, :), POINTER :: cpc_h_res_buffer => NULL(), cpc_s_res_buffer => NULL()
108 : TYPE(rho_atom_coeff), DIMENSION(:, :, :), POINTER :: dcpc_h_in => NULL(), dcpc_s_in => NULL()
109 : #if defined(__TBLITE)
110 : CLASS(mixer_type), ALLOCATABLE :: tb_scc_mixer
111 : #endif
112 : END TYPE mixing_storage_type
113 :
114 : CONTAINS
115 :
116 : ! **************************************************************************************************
117 : !> \brief creates a mixing_storage
118 : !> \param mixing_store ...
119 : !> \param mixing_section ...
120 : !> \param mixing_method ...
121 : !> \param ecut ...
122 : !> \par History
123 : !> 05.2009 created [MI]
124 : !> \author [MI]
125 : ! **************************************************************************************************
126 34320 : SUBROUTINE mixing_storage_create(mixing_store, mixing_section, mixing_method, ecut)
127 : TYPE(mixing_storage_type), INTENT(OUT) :: mixing_store
128 : TYPE(section_vals_type), POINTER :: mixing_section
129 : INTEGER, INTENT(IN) :: mixing_method
130 : REAL(dp), INTENT(IN) :: ecut
131 :
132 : REAL(dp) :: alpha, eps, gcut
133 :
134 8580 : mixing_store%nbuffer = 0
135 8580 : mixing_store%n_simple_mix = 0
136 8580 : mixing_store%ncall = 0
137 25740 : mixing_store%ncall_p = 0
138 8580 : mixing_store%alpha = 1.0_dp
139 8580 : mixing_store%pulay_beta = 1.0_dp
140 8580 : mixing_store%beta = 1.0_dp
141 8580 : mixing_store%alpha_mag = -1.0_dp
142 8580 : mixing_store%beta_mag = -1.0_dp
143 8580 : mixing_store%iter_method = "NoMix"
144 8580 : mixing_store%max_g2 = 2._dp*ecut
145 8580 : mixing_store%gmix_p = .FALSE.
146 8580 : mixing_store%tb_scc_mixer_error = 0.0_dp
147 : mixing_store%tb_scc_mixer_memory = 0
148 : mixing_store%tb_scc_mixer_natom = 0
149 : mixing_store%tb_scc_mixer_ns = 0
150 : mixing_store%tb_scc_mixer_step = 0
151 :
152 8580 : NULLIFY (mixing_store%p_metric)
153 8580 : NULLIFY (mixing_store%kerker_factor)
154 8580 : NULLIFY (mixing_store%kerker_factor_mag)
155 8580 : NULLIFY (mixing_store%special_metric)
156 8580 : NULLIFY (mixing_store%pulay_matrix)
157 8580 : NULLIFY (mixing_store%weight)
158 8580 : NULLIFY (mixing_store%fmat)
159 8580 : NULLIFY (mixing_store%gmat)
160 8580 : NULLIFY (mixing_store%smat)
161 8580 : NULLIFY (mixing_store%acharge)
162 8580 : NULLIFY (mixing_store%dacharge)
163 8580 : NULLIFY (mixing_store%dfbroy)
164 8580 : NULLIFY (mixing_store%ubroy)
165 8580 : NULLIFY (mixing_store%abroy)
166 8580 : NULLIFY (mixing_store%wbroy)
167 8580 : NULLIFY (mixing_store%atlist)
168 8580 : NULLIFY (mixing_store%last_res)
169 8580 : NULLIFY (mixing_store%rhoin)
170 8580 : NULLIFY (mixing_store%rhoin_old)
171 8580 : NULLIFY (mixing_store%delta_res)
172 8580 : NULLIFY (mixing_store%u_vec)
173 8580 : NULLIFY (mixing_store%z_vec)
174 8580 : NULLIFY (mixing_store%drho_buffer)
175 8580 : NULLIFY (mixing_store%rhoin_buffer)
176 8580 : NULLIFY (mixing_store%res_buffer)
177 8580 : NULLIFY (mixing_store%norm_res_buffer)
178 8580 : NULLIFY (mixing_store%ig_global_index)
179 8580 : NULLIFY (mixing_store%paw)
180 8580 : NULLIFY (mixing_store%cpc_h_in)
181 8580 : NULLIFY (mixing_store%cpc_s_in)
182 8580 : NULLIFY (mixing_store%cpc_h_old)
183 8580 : NULLIFY (mixing_store%cpc_s_old)
184 8580 : NULLIFY (mixing_store%dcpc_h_in)
185 8580 : NULLIFY (mixing_store%dcpc_s_in)
186 8580 : NULLIFY (mixing_store%cpc_h_lastres)
187 8580 : NULLIFY (mixing_store%cpc_s_lastres)
188 8580 : NULLIFY (mixing_store%cpc_h_in_buffer)
189 8580 : NULLIFY (mixing_store%cpc_s_in_buffer)
190 8580 : NULLIFY (mixing_store%cpc_h_res_buffer)
191 8580 : NULLIFY (mixing_store%cpc_s_res_buffer)
192 :
193 8580 : CALL section_vals_val_get(mixing_section, "ALPHA", r_val=mixing_store%alpha)
194 8580 : CALL section_vals_val_get(mixing_section, "BETA", r_val=mixing_store%beta)
195 8580 : CALL section_vals_val_get(mixing_section, "ALPHA_MAG", r_val=mixing_store%alpha_mag)
196 8580 : CALL section_vals_val_get(mixing_section, "BETA_MAG", r_val=mixing_store%beta_mag)
197 : ! Fall back to charge-channel values if magnetization parameters are not set
198 8580 : IF (mixing_store%alpha_mag < 0.0_dp) mixing_store%alpha_mag = mixing_store%alpha
199 8580 : IF (mixing_store%beta_mag < 0.0_dp) mixing_store%beta_mag = mixing_store%beta
200 8580 : CALL section_vals_val_get(mixing_section, "N_SIMPLE_MIX", i_val=mixing_store%n_simple_mix)
201 8580 : CALL section_vals_val_get(mixing_section, "NBUFFER", i_val=mixing_store%nbuffer)
202 8580 : CALL section_vals_val_get(mixing_section, "NSKIP", i_val=mixing_store%nskip_mixing)
203 8580 : CALL section_vals_val_get(mixing_section, "MAX_GVEC_EXP", r_val=mixing_store%max_gvec_exp)
204 8580 : CALL section_vals_val_get(mixing_section, "GMIX_P", l_val=mixing_store%gmix_p)
205 :
206 8580 : IF (mixing_store%max_gvec_exp > 0._dp) THEN
207 0 : alpha = 0.25_dp/mixing_store%max_gvec_exp
208 0 : eps = 1.e-4_dp
209 0 : gcut = exp_radius(3, alpha, eps, 1.0_dp)
210 0 : mixing_store%max_g2 = gcut*gcut
211 : END IF
212 :
213 8592 : SELECT CASE (mixing_method)
214 : CASE (gspace_mixing_nr)
215 12 : mixing_store%nbuffer = 1
216 : CASE (pulay_mixing_nr)
217 36 : CALL section_vals_val_get(mixing_section, "PULAY_ALPHA", r_val=mixing_store%pulay_alpha)
218 36 : CALL section_vals_val_get(mixing_section, "PULAY_BETA", r_val=mixing_store%pulay_beta)
219 : CASE (broyden_mixing_nr)
220 942 : CALL section_vals_val_get(mixing_section, "BROY_W0", r_val=mixing_store%broy_w0)
221 942 : mixing_store%bconst = 20.0_dp
222 : CASE (modified_broyden_mixing_nr)
223 8 : CALL section_vals_val_get(mixing_section, "BROY_W0", r_val=mixing_store%broy_w0)
224 8 : CALL section_vals_val_get(mixing_section, "BROY_WREF", r_val=mixing_store%wc)
225 8 : CALL section_vals_val_get(mixing_section, "BROY_WMAX", r_val=mixing_store%wmax)
226 8 : mixing_store%bconst = 20.0_dp
227 : CASE (multisecant_mixing_nr)
228 0 : CALL section_vals_val_get(mixing_section, "REGULARIZATION", r_val=mixing_store%reg_par)
229 0 : CALL section_vals_val_get(mixing_section, "MAX_STEP", r_val=mixing_store%sigma_max)
230 0 : CALL section_vals_val_get(mixing_section, "R_FACTOR", r_val=mixing_store%r_step)
231 : CASE (new_pulay_mixing_nr)
232 10 : CALL section_vals_val_get(mixing_section, "PULAY_ALPHA", r_val=mixing_store%pulay_alpha)
233 10 : CALL section_vals_val_get(mixing_section, "PULAY_BETA", r_val=mixing_store%pulay_beta)
234 10 : CALL section_vals_val_get(mixing_section, "QM", r_val=mixing_store%qm)
235 10 : CALL section_vals_val_get(mixing_section, "QK", r_val=mixing_store%qk)
236 8590 : CALL section_vals_val_get(mixing_section, "QKAPPA", r_val=mixing_store%qkappa)
237 : END SELECT
238 :
239 8580 : END SUBROUTINE mixing_storage_create
240 :
241 : ! **************************************************************************************************
242 : !> \brief releases a mixing_storage
243 : !> \param mixing_store ...
244 : !> \par History
245 : !> 05.2009 created [MI]
246 : !> \author [MI]
247 : ! **************************************************************************************************
248 8581 : SUBROUTINE mixing_storage_release(mixing_store)
249 : TYPE(mixing_storage_type), INTENT(INOUT) :: mixing_store
250 :
251 : INTEGER :: i, j, k
252 :
253 8581 : IF (ALLOCATED(mixing_store%aux_buffer)) DEALLOCATE (mixing_store%aux_buffer)
254 8581 : mixing_store%naux = 0
255 :
256 8581 : IF (ASSOCIATED(mixing_store%kerker_factor)) THEN
257 612 : DEALLOCATE (mixing_store%kerker_factor)
258 : END IF
259 :
260 8581 : IF (ASSOCIATED(mixing_store%kerker_factor_mag)) THEN
261 612 : DEALLOCATE (mixing_store%kerker_factor_mag)
262 : END IF
263 :
264 8581 : IF (ASSOCIATED(mixing_store%special_metric)) THEN
265 612 : DEALLOCATE (mixing_store%special_metric)
266 : END IF
267 :
268 8581 : IF (ASSOCIATED(mixing_store%pulay_matrix)) THEN
269 44 : DEALLOCATE (mixing_store%pulay_matrix)
270 : END IF
271 :
272 8581 : IF (ASSOCIATED(mixing_store%rhoin_buffer)) THEN
273 98 : DO i = 1, SIZE(mixing_store%rhoin_buffer, 2)
274 416 : DO j = 1, SIZE(mixing_store%rhoin_buffer, 1)
275 372 : DEALLOCATE (mixing_store%rhoin_buffer(j, i)%cc)
276 : END DO
277 : END DO
278 44 : DEALLOCATE (mixing_store%rhoin_buffer)
279 : END IF
280 :
281 8581 : IF (ASSOCIATED(mixing_store%paw)) THEN
282 16 : DEALLOCATE (mixing_store%paw)
283 : END IF
284 8581 : IF (ASSOCIATED(mixing_store%cpc_h_in)) THEN
285 38 : DO j = 1, SIZE(mixing_store%cpc_h_in, 2)
286 214 : DO k = 1, SIZE(mixing_store%cpc_h_in, 1)
287 198 : IF (ASSOCIATED(mixing_store%cpc_h_in(k, j)%r_coef)) THEN
288 106 : DEALLOCATE (mixing_store%cpc_h_in(k, j)%r_coef)
289 106 : DEALLOCATE (mixing_store%cpc_s_in(k, j)%r_coef)
290 : END IF
291 : END DO
292 : END DO
293 16 : DEALLOCATE (mixing_store%cpc_h_in)
294 16 : DEALLOCATE (mixing_store%cpc_s_in)
295 : END IF
296 8581 : IF (ASSOCIATED(mixing_store%cpc_h_old)) THEN
297 30 : DO j = 1, SIZE(mixing_store%cpc_h_old, 2)
298 174 : DO k = 1, SIZE(mixing_store%cpc_h_old, 1)
299 162 : IF (ASSOCIATED(mixing_store%cpc_h_old(k, j)%r_coef)) THEN
300 102 : DEALLOCATE (mixing_store%cpc_h_old(k, j)%r_coef)
301 102 : DEALLOCATE (mixing_store%cpc_s_old(k, j)%r_coef)
302 : END IF
303 : END DO
304 : END DO
305 12 : DEALLOCATE (mixing_store%cpc_h_old)
306 12 : DEALLOCATE (mixing_store%cpc_s_old)
307 : END IF
308 8581 : IF (ASSOCIATED(mixing_store%cpc_h_in_buffer)) THEN
309 4 : DO i = 1, SIZE(mixing_store%cpc_h_in_buffer, 3)
310 20 : DO j = 1, SIZE(mixing_store%cpc_h_in_buffer, 2)
311 98 : DO k = 1, SIZE(mixing_store%cpc_h_in_buffer, 1)
312 96 : IF (ASSOCIATED(mixing_store%cpc_h_in_buffer(k, j, i)%r_coef)) THEN
313 10 : DEALLOCATE (mixing_store%cpc_h_in_buffer(k, j, i)%r_coef)
314 10 : DEALLOCATE (mixing_store%cpc_s_in_buffer(k, j, i)%r_coef)
315 : END IF
316 : END DO
317 : END DO
318 : END DO
319 2 : DEALLOCATE (mixing_store%cpc_h_in_buffer)
320 2 : DEALLOCATE (mixing_store%cpc_s_in_buffer)
321 : END IF
322 8581 : IF (ASSOCIATED(mixing_store%cpc_h_res_buffer)) THEN
323 4 : DO i = 1, SIZE(mixing_store%cpc_h_res_buffer, 3)
324 20 : DO j = 1, SIZE(mixing_store%cpc_h_res_buffer, 2)
325 98 : DO k = 1, SIZE(mixing_store%cpc_h_res_buffer, 1)
326 96 : IF (ASSOCIATED(mixing_store%cpc_h_res_buffer(k, j, i)%r_coef)) THEN
327 10 : DEALLOCATE (mixing_store%cpc_h_res_buffer(k, j, i)%r_coef)
328 10 : DEALLOCATE (mixing_store%cpc_s_res_buffer(k, j, i)%r_coef)
329 : END IF
330 : END DO
331 : END DO
332 : END DO
333 2 : DEALLOCATE (mixing_store%cpc_h_res_buffer)
334 2 : DEALLOCATE (mixing_store%cpc_s_res_buffer)
335 : END IF
336 :
337 8581 : IF (ASSOCIATED(mixing_store%dcpc_h_in)) THEN
338 30 : DO i = 1, SIZE(mixing_store%dcpc_h_in, 3)
339 174 : DO j = 1, SIZE(mixing_store%dcpc_h_in, 2)
340 1266 : DO k = 1, SIZE(mixing_store%dcpc_h_in, 1)
341 1248 : IF (ASSOCIATED(mixing_store%dcpc_h_in(k, j, i)%r_coef)) THEN
342 810 : DEALLOCATE (mixing_store%dcpc_h_in(k, j, i)%r_coef)
343 810 : DEALLOCATE (mixing_store%dcpc_s_in(k, j, i)%r_coef)
344 : END IF
345 : END DO
346 : END DO
347 : END DO
348 12 : DEALLOCATE (mixing_store%dcpc_h_in)
349 12 : DEALLOCATE (mixing_store%dcpc_s_in)
350 : END IF
351 8581 : IF (ASSOCIATED(mixing_store%cpc_h_lastres)) THEN
352 30 : DO j = 1, SIZE(mixing_store%cpc_h_lastres, 2)
353 174 : DO k = 1, SIZE(mixing_store%cpc_h_lastres, 1)
354 162 : IF (ASSOCIATED(mixing_store%cpc_h_lastres(k, j)%r_coef)) THEN
355 102 : DEALLOCATE (mixing_store%cpc_h_lastres(k, j)%r_coef)
356 102 : DEALLOCATE (mixing_store%cpc_s_lastres(k, j)%r_coef)
357 : END IF
358 : END DO
359 : END DO
360 12 : DEALLOCATE (mixing_store%cpc_h_lastres)
361 12 : DEALLOCATE (mixing_store%cpc_s_lastres)
362 : END IF
363 :
364 8581 : IF (ASSOCIATED(mixing_store%res_buffer)) THEN
365 1300 : DO i = 1, SIZE(mixing_store%res_buffer, 2)
366 5918 : DO j = 1, SIZE(mixing_store%res_buffer, 1)
367 5318 : DEALLOCATE (mixing_store%res_buffer(j, i)%cc)
368 : END DO
369 : END DO
370 600 : DEALLOCATE (mixing_store%res_buffer)
371 : END IF
372 :
373 8581 : IF (ASSOCIATED(mixing_store%norm_res_buffer)) THEN
374 0 : DEALLOCATE (mixing_store%norm_res_buffer)
375 : END IF
376 :
377 8581 : IF (ASSOCIATED(mixing_store%ig_global_index)) THEN
378 0 : DEALLOCATE (mixing_store%ig_global_index)
379 : END IF
380 :
381 8581 : IF (ASSOCIATED(mixing_store%drho_buffer)) THEN
382 1202 : DO i = 1, SIZE(mixing_store%drho_buffer, 2)
383 5502 : DO j = 1, SIZE(mixing_store%drho_buffer, 1)
384 4946 : DEALLOCATE (mixing_store%drho_buffer(j, i)%cc)
385 : END DO
386 : END DO
387 556 : DEALLOCATE (mixing_store%drho_buffer)
388 : END IF
389 :
390 8581 : IF (ASSOCIATED(mixing_store%last_res)) THEN
391 1202 : DO i = 1, SIZE(mixing_store%last_res)
392 1202 : DEALLOCATE (mixing_store%last_res(i)%cc)
393 : END DO
394 556 : DEALLOCATE (mixing_store%last_res)
395 : END IF
396 :
397 8581 : IF (ASSOCIATED(mixing_store%rhoin)) THEN
398 1324 : DO i = 1, SIZE(mixing_store%rhoin)
399 1324 : DEALLOCATE (mixing_store%rhoin(i)%cc)
400 : END DO
401 612 : DEALLOCATE (mixing_store%rhoin)
402 : END IF
403 :
404 8581 : IF (ASSOCIATED(mixing_store%rhoin_old)) THEN
405 1202 : DO i = 1, SIZE(mixing_store%rhoin_old)
406 1202 : DEALLOCATE (mixing_store%rhoin_old(i)%cc)
407 : END DO
408 556 : DEALLOCATE (mixing_store%rhoin_old)
409 : END IF
410 :
411 8581 : IF (ASSOCIATED(mixing_store%p_metric)) THEN
412 556 : DEALLOCATE (mixing_store%p_metric)
413 : END IF
414 :
415 8581 : IF (ASSOCIATED(mixing_store%weight)) THEN
416 2 : DEALLOCATE (mixing_store%weight)
417 : END IF
418 :
419 8581 : IF (ASSOCIATED(mixing_store%fmat)) THEN
420 0 : DEALLOCATE (mixing_store%fmat)
421 : END IF
422 :
423 8581 : IF (ASSOCIATED(mixing_store%acharge)) THEN
424 48 : DEALLOCATE (mixing_store%acharge)
425 : END IF
426 8581 : IF (ASSOCIATED(mixing_store%dacharge)) THEN
427 48 : DEALLOCATE (mixing_store%dacharge)
428 : END IF
429 8581 : IF (ASSOCIATED(mixing_store%dfbroy)) THEN
430 48 : DEALLOCATE (mixing_store%dfbroy)
431 : END IF
432 8581 : IF (ASSOCIATED(mixing_store%ubroy)) THEN
433 48 : DEALLOCATE (mixing_store%ubroy)
434 : END IF
435 8581 : IF (ASSOCIATED(mixing_store%abroy)) THEN
436 48 : DEALLOCATE (mixing_store%abroy)
437 : END IF
438 8581 : IF (ASSOCIATED(mixing_store%wbroy)) THEN
439 48 : DEALLOCATE (mixing_store%wbroy)
440 : END IF
441 8581 : IF (ASSOCIATED(mixing_store%atlist)) THEN
442 48 : DEALLOCATE (mixing_store%atlist)
443 : END IF
444 : #if defined(__TBLITE)
445 8581 : IF (ALLOCATED(mixing_store%tb_scc_mixer)) THEN
446 4 : DEALLOCATE (mixing_store%tb_scc_mixer)
447 : END IF
448 : #endif
449 :
450 8581 : IF (ASSOCIATED(mixing_store%delta_res)) THEN
451 0 : DO i = 1, SIZE(mixing_store%delta_res, 2)
452 0 : DO j = 1, SIZE(mixing_store%delta_res, 1)
453 0 : DEALLOCATE (mixing_store%delta_res(j, i)%cc)
454 : END DO
455 : END DO
456 0 : DEALLOCATE (mixing_store%delta_res)
457 : END IF
458 :
459 8581 : IF (ASSOCIATED(mixing_store%u_vec)) THEN
460 0 : DO i = 1, SIZE(mixing_store%u_vec, 2)
461 0 : DO j = 1, SIZE(mixing_store%u_vec, 1)
462 0 : DEALLOCATE (mixing_store%u_vec(j, i)%cc)
463 : END DO
464 : END DO
465 0 : DEALLOCATE (mixing_store%u_vec)
466 : END IF
467 :
468 8581 : IF (ASSOCIATED(mixing_store%z_vec)) THEN
469 0 : DO i = 1, SIZE(mixing_store%z_vec, 2)
470 0 : DO j = 1, SIZE(mixing_store%z_vec, 1)
471 0 : DEALLOCATE (mixing_store%z_vec(j, i)%cc)
472 : END DO
473 : END DO
474 0 : DEALLOCATE (mixing_store%z_vec)
475 : END IF
476 :
477 8581 : END SUBROUTINE mixing_storage_release
478 :
479 : ! **************************************************************************************************
480 : !> \brief Create CP2K input section for SCF or active-space density mixing
481 : !> \param section ...
482 : !> \param ls_scf ...
483 : !> \param active_space Customize descriptions and deprecations for active-space mixing
484 : !> \date 20.02.2009
485 : !> \par History
486 : !> 02.2015 moved here from input_cp2k_dft.F, modified for use in LS SCF
487 : !> [Patrick Seewald]
488 : !> \author MI
489 : !> \version 1.0
490 : ! **************************************************************************************************
491 18463 : SUBROUTINE create_mixing_section(section, ls_scf, active_space)
492 :
493 : TYPE(section_type), POINTER :: section
494 : LOGICAL, INTENT(IN), OPTIONAL :: ls_scf, active_space
495 :
496 18463 : CHARACTER(LEN=:), ALLOCATABLE :: description
497 : CHARACTER(LEN=default_string_length) :: section_name
498 : INTEGER :: default_mix
499 : LOGICAL :: as, ls
500 : TYPE(keyword_type), POINTER :: keyword
501 :
502 18463 : CPASSERT(.NOT. ASSOCIATED(section))
503 :
504 18463 : IF (PRESENT(ls_scf)) THEN
505 2929 : IF (ls_scf) THEN
506 : ls = .TRUE.
507 : ELSE
508 1455 : ls = .FALSE.
509 : END IF
510 : ELSE
511 : ls = .FALSE.
512 : END IF
513 :
514 18463 : as = .FALSE.
515 18463 : IF (PRESENT(active_space)) as = active_space
516 1474 : CPASSERT(.NOT. (ls .AND. as))
517 :
518 18463 : IF (ls) THEN
519 1474 : section_name = "RHO_MIXING"
520 : ELSE
521 16989 : section_name = "MIXING"
522 : END IF
523 :
524 18463 : IF (as) THEN
525 : description = "Define type and parameters for mixing the active-space density matrix "// &
526 : "during self-consistent embedding. This mixing is independent of the SCF "// &
527 : "method used for the reference calculation, including OT. NBUFFER controls "// &
528 1474 : "the mixing history; N_SIMPLE_MIX controls the initial direct-mixing steps."
529 : ELSE
530 : description = "Define type and parameters for mixing "// &
531 : "procedures to be applied to the density matrix. Normally, "// &
532 : "only one type of mixing method should be accepted. The mixing "// &
533 : "procedures activated by this section are only active for diagonalization "// &
534 : "methods and linear scaling SCF, i.e. not with minimization methods based "// &
535 16989 : "on OT."
536 : END IF
537 :
538 : CALL section_create(section, __LOCATION__, &
539 : name=section_name, &
540 : description=description, &
541 : n_keywords=16, &
542 : n_subsections=0, &
543 18463 : repeats=.FALSE.)
544 :
545 18463 : NULLIFY (keyword)
546 :
547 : CALL keyword_create(keyword, __LOCATION__, &
548 : name="_SECTION_PARAMETERS_", &
549 : description="Controls the activation of the mixing procedure", &
550 : usage="&MIXING ON", &
551 : default_l_val=.TRUE., &
552 18463 : lone_keyword_l_val=.TRUE.)
553 18463 : CALL section_add_keyword(section, keyword)
554 18463 : CALL keyword_release(keyword)
555 :
556 18463 : IF (.NOT. ls) THEN
557 16989 : default_mix = direct_p_mix
558 : ELSE
559 1474 : default_mix = broy_mix
560 : END IF
561 :
562 : CALL keyword_create(keyword, __LOCATION__, &
563 : name="METHOD", &
564 : description="Mixing method to be applied", &
565 : repeats=.FALSE., &
566 : usage="METHOD KERKER_MIXING", &
567 : default_i_val=default_mix, &
568 : enum_c_vals=s2a("NONE", &
569 : "DIRECT_P_MIXING", &
570 : "KERKER_MIXING", &
571 : "PULAY_MIXING", &
572 : "BROYDEN_MIXING", &
573 : "MODIFIED_BROYDEN_MIXING", &
574 : "MULTISECANT_MIXING", &
575 : "NEW_PULAY_MIXING"), &
576 : enum_i_vals=[no_mix, direct_p_mix, kerker_mix, pulay_mix, broy_mix, &
577 : modified_broy_mix, &
578 : multisec_mix, new_pulay_mix], &
579 : enum_desc=s2a("No mixing is applied", &
580 : "Direct mixing of new and old density matrices", &
581 : "Mixing of the potential in reciprocal space using the Kerker damping", &
582 : "Pulay mixing", &
583 : "Original CP2K Broyden mixing with a constant BROY_W0 regularization", &
584 : "Modified Broyden mixing with dynamic residual weights controlled by BROY_W0, "// &
585 : "BROY_WREF, and BROY_WMAX", &
586 : "Multisecant scheme for mixing", "New Pulay mixing using Sundararaman et al.'s "// &
587 : "metric and preconditioner, with improved convergence behavior and suitable for "// &
588 18463 : "grand canonical SCF"))
589 :
590 18463 : CALL section_add_keyword(section, keyword)
591 18463 : CALL keyword_release(keyword)
592 :
593 : CALL keyword_create(keyword, __LOCATION__, &
594 : name="ALPHA", &
595 : description="Fraction of new density to be included", &
596 : repeats=.FALSE., &
597 : n_var=1, &
598 : type_of_var=real_t, &
599 : default_r_val=0.4_dp, &
600 18463 : usage="ALPHA 0.2")
601 18463 : CALL section_add_keyword(section, keyword)
602 18463 : CALL keyword_release(keyword)
603 :
604 : CALL keyword_create(keyword, __LOCATION__, &
605 : name="ALPHA_MAG", &
606 : description="Fraction of new magnetization density to be included "// &
607 : "(for spin-polarized calculations, ispin=2 channel after rho_total/m transform). "// &
608 : "A negative value (default) means: use the same value as ALPHA. "// &
609 : "For magnetic transition-metal systems, a larger value (e.g. 0.8-1.6) "// &
610 : "than ALPHA often improves convergence.", &
611 : repeats=.FALSE., &
612 : n_var=1, &
613 : type_of_var=real_t, &
614 : default_r_val=-1.0_dp, &
615 18463 : usage="ALPHA_MAG 0.8")
616 18463 : CALL section_add_keyword(section, keyword)
617 18463 : CALL keyword_release(keyword)
618 :
619 : CALL keyword_create(keyword, __LOCATION__, &
620 : name="BETA", &
621 : description="Denominator parameter in Kerker damping "// &
622 : "introduced to suppress charge sloshing: "// &
623 : "rho_mix(g) = rho_in(g) + alpha*g^2/(g^2 + beta^2)*(rho_out(g)-rho_in(g))", &
624 : repeats=.FALSE., &
625 : n_var=1, &
626 : type_of_var=real_t, &
627 : default_r_val=0.5_dp, &
628 : unit_str="bohr^-1", &
629 18463 : usage="BETA 1.5")
630 18463 : CALL section_add_keyword(section, keyword)
631 18463 : CALL keyword_release(keyword)
632 :
633 : CALL keyword_create(keyword, __LOCATION__, &
634 : name="BETA_MAG", &
635 : description="Kerker damping parameter for the magnetization channel "// &
636 : "(for spin-polarized calculations). A negative value (default) means: "// &
637 : "use the same value as BETA. Set to 0.0 to disable Kerker screening "// &
638 : "on the magnetization density, which avoids suppression of long-range "// &
639 : "magnetic order formation in transition-metal systems.", &
640 : repeats=.FALSE., &
641 : n_var=1, &
642 : type_of_var=real_t, &
643 : default_r_val=-1.0_dp, &
644 : unit_str="bohr^-1", &
645 18463 : usage="BETA_MAG 0.0")
646 18463 : CALL section_add_keyword(section, keyword)
647 18463 : CALL keyword_release(keyword)
648 :
649 : CALL keyword_create(keyword, __LOCATION__, &
650 : name="PULAY_ALPHA", &
651 : description="Fraction of new density to be added to the Pulay expansion", &
652 : repeats=.FALSE., &
653 : n_var=1, &
654 : type_of_var=real_t, &
655 : default_r_val=0.0_dp, &
656 18463 : usage="PULAY_ALPHA 0.2")
657 18463 : CALL section_add_keyword(section, keyword)
658 18463 : CALL keyword_release(keyword)
659 :
660 : CALL keyword_create(keyword, __LOCATION__, &
661 : name="PULAY_BETA", &
662 : description="Fraction of residual contribution to be added to Pulay expansion", &
663 : repeats=.FALSE., &
664 : n_var=1, &
665 : type_of_var=real_t, &
666 : default_r_val=1.0_dp, &
667 18463 : usage="PULAY_BETA 0.2")
668 18463 : CALL section_add_keyword(section, keyword)
669 18463 : CALL keyword_release(keyword)
670 :
671 : CALL keyword_create(keyword, __LOCATION__, &
672 : name="QM", &
673 : description="The control parameter in the numerator of the Pulay metric "// &
674 : "used in the new Pulay mixing, introduced to suppress charge sloshing: "// &
675 : "Pulay metric: M(g) = (g^2 + qm^2 + qkapa^2)/(g^2 + qkapa^2)", &
676 : repeats=.FALSE., &
677 : n_var=1, &
678 : type_of_var=real_t, &
679 : default_r_val=0.75_dp, &
680 : unit_str="bohr^-1", &
681 18463 : usage="QM 0.75")
682 18463 : CALL section_add_keyword(section, keyword)
683 18463 : CALL keyword_release(keyword)
684 :
685 : CALL keyword_create(keyword, __LOCATION__, &
686 : name="QK", &
687 : description="The control parameter in the denominator of the Kerker preconditioner "// &
688 : "used in the new Pulay mixing, introduced to suppress charge sloshing: "// &
689 : "Kerker preconditioner: K(g) = alpha * (g^2 + qkapa^2)/(g^2 + qk^2 + qkapa^2)", &
690 : repeats=.FALSE., &
691 : n_var=1, &
692 : type_of_var=real_t, &
693 : default_r_val=3.0_dp, &
694 : unit_str="bohr^-1", &
695 18463 : usage="QK 3.0")
696 18463 : CALL section_add_keyword(section, keyword)
697 18463 : CALL keyword_release(keyword)
698 :
699 : CALL keyword_create(keyword, __LOCATION__, &
700 : name="QKAPPA", &
701 : description="The control parameter in the numerator and the denominator of "// &
702 : "the Pulay metric and the Kerker preconditioner "// &
703 : "used in the new Pulay mixing, introduced to ensure a finite and well-defined "// &
704 : "Pulay metric at g=0 and a non-zero and well-defined Kerker preconditioner at g=0", &
705 : repeats=.FALSE., &
706 : n_var=1, &
707 : type_of_var=real_t, &
708 : default_r_val=0.25_dp, &
709 : unit_str="bohr^-1", &
710 18463 : usage="QKAPPA 0.25")
711 18463 : CALL section_add_keyword(section, keyword)
712 18463 : CALL keyword_release(keyword)
713 :
714 18463 : IF (as) THEN
715 : CALL keyword_create(keyword, __LOCATION__, name="NMIXING", &
716 : description="Legacy SCF DIIS keyword with no effect on active-space mixing", &
717 : deprecation_notice="NMIXING has no effect in ACTIVE_SPACE%MIXING. "// &
718 : "Use NBUFFER for the history length or N_SIMPLE_MIX for initial "// &
719 : "direct-mixing steps; these have different meanings.", &
720 1474 : usage="NMIXING 1", default_i_val=2)
721 : ELSE
722 : CALL keyword_create(keyword, __LOCATION__, name="NMIXING", &
723 : description="Minimal number of density mixing (should be greater than 0), "// &
724 : "before starting DIIS", &
725 16989 : usage="NMIXING 1", default_i_val=2)
726 : END IF
727 18463 : CALL section_add_keyword(section, keyword)
728 18463 : CALL keyword_release(keyword)
729 :
730 : CALL keyword_create(keyword, __LOCATION__, name="NBUFFER", &
731 : variants=s2a("NPULAY", "NBROYDEN", "NMULTISECANT"), &
732 : description="Number of previous steps stored for the actual mixing scheme", &
733 18463 : usage="NBUFFER 2", default_i_val=4)
734 18463 : CALL section_add_keyword(section, keyword)
735 18463 : CALL keyword_release(keyword)
736 :
737 : CALL keyword_create(keyword, __LOCATION__, &
738 : name="BROY_W0", &
739 : description="Regularization weight used in Broyden mixing. For the original "// &
740 : "BROYDEN_MIXING method this is the constant diagonal regularization of the "// &
741 : "small Broyden system. For MODIFIED_BROYDEN_MIXING it is the corresponding "// &
742 : "diagonal regularization of the dynamically weighted Broyden system. "// &
743 : "The default follows tblite.", &
744 : repeats=.FALSE., &
745 : n_var=1, &
746 : type_of_var=real_t, &
747 : default_r_val=0.01_dp, &
748 18463 : usage="BROY_W0 0.03")
749 18463 : CALL section_add_keyword(section, keyword)
750 18463 : CALL keyword_release(keyword)
751 :
752 : CALL keyword_create(keyword, __LOCATION__, &
753 : name="BROY_WREF", &
754 : description="Reference factor for the dynamic residual weight. This keyword "// &
755 : "is only used by MODIFIED_BROYDEN_MIXING; the original BROYDEN_MIXING path "// &
756 : "is unchanged. The effective history weight is proportional to BROY_WREF "// &
757 : "divided by the residual norm, clipped to the interval [1, BROY_WMAX]. "// &
758 : "The default follows tblite.", &
759 : repeats=.FALSE., &
760 : n_var=1, &
761 : type_of_var=real_t, &
762 : default_r_val=0.01_dp, &
763 18463 : usage="BROY_WREF 0.01")
764 18463 : CALL section_add_keyword(section, keyword)
765 18463 : CALL keyword_release(keyword)
766 :
767 : CALL keyword_create(keyword, __LOCATION__, &
768 : name="BROY_WMAX", &
769 : description="Upper bound for the dynamic residual weight. This keyword is "// &
770 : "only used by MODIFIED_BROYDEN_MIXING; the original BROYDEN_MIXING path is "// &
771 : "unchanged. The lower bound is fixed to 1.0. The default follows tblite.", &
772 : repeats=.FALSE., &
773 : n_var=1, &
774 : type_of_var=real_t, &
775 : default_r_val=100000.0_dp, &
776 18463 : usage="BROY_WMAX 100000.0")
777 18463 : CALL section_add_keyword(section, keyword)
778 18463 : CALL keyword_release(keyword)
779 :
780 : CALL keyword_create(keyword, __LOCATION__, &
781 : name="REGULARIZATION", &
782 : description="Regularization parameter to stabilize "// &
783 : "the inversion of the residual matrix {Yn^t Yn} in the "// &
784 : "multisecant mixing scheme (noise)", &
785 : repeats=.FALSE., &
786 : n_var=1, &
787 : type_of_var=real_t, &
788 : default_r_val=0.00001_dp, &
789 18463 : usage="REGULARIZATION 0.000001")
790 18463 : CALL section_add_keyword(section, keyword)
791 18463 : CALL keyword_release(keyword)
792 :
793 : CALL keyword_create(keyword, __LOCATION__, &
794 : name="MAX_STEP", &
795 : description="Upper bound for the magnitude of the "// &
796 : "unpredicted step size in the update by the "// &
797 : "multisecant mixing scheme", &
798 : repeats=.FALSE., &
799 : n_var=1, &
800 : type_of_var=real_t, &
801 : default_r_val=0.1_dp, &
802 18463 : usage="MAX_STEP .2")
803 18463 : CALL section_add_keyword(section, keyword)
804 18463 : CALL keyword_release(keyword)
805 :
806 : CALL keyword_create(keyword, __LOCATION__, &
807 : name="R_FACTOR", &
808 : description="Control factor for the magnitude of the "// &
809 : "unpredicted step size in the update by the "// &
810 : "multisecant mixing scheme", &
811 : repeats=.FALSE., &
812 : n_var=1, &
813 : type_of_var=real_t, &
814 : default_r_val=0.05_dp, &
815 18463 : usage="R_FACTOR .12")
816 18463 : CALL section_add_keyword(section, keyword)
817 18463 : CALL keyword_release(keyword)
818 :
819 : CALL keyword_create(keyword, __LOCATION__, name="NSKIP", &
820 : variants=["NSKIP_MIXING"], &
821 : description="Number of initial iteration for which the mixing is skipped", &
822 36926 : usage="NSKIP 10", default_i_val=0)
823 18463 : CALL section_add_keyword(section, keyword)
824 18463 : CALL keyword_release(keyword)
825 :
826 18463 : IF (as) THEN
827 : description = "Number of initial direct-mixing iterations after NSKIP before starting "// &
828 1474 : "the selected mixing procedure"
829 : ELSE
830 16989 : description = "Number of kerker damping iterations before starting other mixing procedures"
831 : END IF
832 : CALL keyword_create(keyword, __LOCATION__, name="N_SIMPLE_MIX", &
833 : variants=["NSIMPLEMIX"], &
834 : description=description, &
835 36926 : usage="NSIMPLEMIX", default_i_val=0)
836 18463 : CALL section_add_keyword(section, keyword)
837 18463 : CALL keyword_release(keyword)
838 :
839 : CALL keyword_create(keyword, __LOCATION__, name="MAX_GVEC_EXP", &
840 : description="Restricts the G-space mixing to lower part of G-vector spectrum,"// &
841 : " up to a G0, by assigning the exponent of the Gaussian that can be "// &
842 : "represented by vectors smaller than G0 within a certain accuracy. ", &
843 : repeats=.FALSE., &
844 : n_var=1, &
845 : type_of_var=real_t, &
846 : default_r_val=-1._dp, &
847 18463 : usage="MAX_GVEC_EXP 3.")
848 18463 : CALL section_add_keyword(section, keyword)
849 18463 : CALL keyword_release(keyword)
850 :
851 : CALL keyword_create(keyword, __LOCATION__, name="GMIX_P", &
852 : description="Activate the mixing of the density matrix, using the same"// &
853 : " mixing coefficient applied for the g-space mixing.", &
854 : repeats=.FALSE., &
855 : lone_keyword_l_val=.TRUE., &
856 : default_l_val=.FALSE., &
857 18463 : usage="GMIX_P")
858 18463 : CALL section_add_keyword(section, keyword)
859 18463 : CALL keyword_release(keyword)
860 :
861 18463 : END SUBROUTINE create_mixing_section
862 :
863 0 : END MODULE qs_density_mixing_types
|
