Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief orbital transformations
10 : !> \par History
11 : !> Added Taylor expansion based computation of the matrix functions (01.2004)
12 : !> added additional rotation variables for non-equivalent occupied orbs (08.2004)
13 : !> \author Joost VandeVondele (06.2002)
14 : ! **************************************************************************************************
15 : MODULE qs_ot_types
16 : USE bibliography, ONLY: VandeVondele2003,&
17 : Weber2008,&
18 : cite_reference
19 : USE cp_blacs_env, ONLY: cp_blacs_env_release,&
20 : cp_blacs_env_type
21 : USE cp_dbcsr_api, ONLY: dbcsr_init_p,&
22 : dbcsr_p_type,&
23 : dbcsr_release_p,&
24 : dbcsr_set,&
25 : dbcsr_type,&
26 : dbcsr_type_no_symmetry
27 : USE cp_dbcsr_operations, ONLY: cp_dbcsr_m_by_n_from_row_template,&
28 : cp_dbcsr_m_by_n_from_template,&
29 : dbcsr_allocate_matrix_set,&
30 : dbcsr_deallocate_matrix_set
31 : USE cp_fm_struct, ONLY: cp_fm_struct_get,&
32 : cp_fm_struct_type
33 : USE input_constants, ONLY: &
34 : cholesky_reduce, ls_2pnt, ls_3pnt, ls_adapt, ls_gold, ls_none, ot_algo_irac, &
35 : ot_algo_taylor_or_diag, ot_chol_irac, ot_lwdn_irac, ot_mini_broyden, ot_mini_cg, &
36 : ot_mini_diis, ot_mini_sd, ot_poly_irac, ot_precond_full_all, ot_precond_full_kinetic, &
37 : ot_precond_full_single, ot_precond_full_single_inverse, ot_precond_none, &
38 : ot_precond_s_inverse, ot_precond_solver_default, ot_precond_solver_direct, &
39 : ot_precond_solver_inv_chol, ot_precond_solver_update
40 : USE input_section_types, ONLY: section_vals_type,&
41 : section_vals_val_get
42 : USE kinds, ONLY: dp
43 : USE message_passing, ONLY: mp_para_env_release,&
44 : mp_para_env_type
45 : USE preconditioner_types, ONLY: preconditioner_type
46 : #include "./base/base_uses.f90"
47 :
48 : IMPLICIT NONE
49 :
50 : PRIVATE
51 :
52 : PUBLIC :: qs_ot_type
53 : PUBLIC :: qs_ot_settings_type
54 : PUBLIC :: qs_ot_destroy
55 : PUBLIC :: qs_ot_allocate
56 : PUBLIC :: qs_ot_init
57 : PUBLIC :: qs_ot_settings_init
58 : PUBLIC :: ot_readwrite_input
59 :
60 : ! **************************************************************************************************
61 : !> \brief notice, this variable needs to be copyable, needed for spins as e.g. in qs_ot_scf
62 : ! **************************************************************************************************
63 : TYPE qs_ot_settings_type
64 : LOGICAL :: do_rotation = .FALSE., do_ener = .FALSE.
65 : LOGICAL :: ks = .FALSE.
66 : CHARACTER(LEN=4) :: ot_method = ""
67 : CHARACTER(LEN=3) :: ot_algorithm = ""
68 : CHARACTER(LEN=4) :: line_search_method = ""
69 : CHARACTER(LEN=20) :: preconditioner_name = ""
70 : INTEGER :: preconditioner_type = -1
71 : INTEGER :: cholesky_type = -1
72 : INTEGER :: ot_state = 0
73 : CHARACTER(LEN=20) :: precond_solver_name = ""
74 : INTEGER :: precond_solver_type = -1
75 : LOGICAL :: safer_diis = .FALSE.
76 : REAL(KIND=dp) :: ds_min = -1.0_dp
77 : REAL(KIND=dp) :: energy_gap = -1.0_dp
78 : INTEGER :: diis_m = -1
79 : INTEGER :: max_scf_diis = 0
80 : REAL(KIND=dp) :: gold_target = -1.0_dp
81 : REAL(KIND=dp) :: eps_taylor = -1.0_dp ! minimum accuracy of Taylor expansion
82 : INTEGER :: max_taylor = -1 ! maximum order of Taylor expansion before switching to diagonalization
83 : INTEGER :: irac_degree = -1 ! used to control the refinement polynomial degree
84 : INTEGER :: max_irac = -1 ! maximum number of iteration for refinement
85 : REAL(KIND=dp) :: eps_irac = -1.0_dp ! target accuracy for refinement
86 : REAL(KIND=dp) :: eps_irac_quick_exit = -1.0_dp
87 : REAL(KIND=dp) :: eps_irac_filter_matrix = -1.0_dp
88 : REAL(KIND=dp) :: eps_irac_switch = -1.0_dp
89 : LOGICAL :: on_the_fly_loc = .FALSE.
90 : CHARACTER(LEN=4) :: ortho_irac = ""
91 : LOGICAL :: occupation_preconditioner = .FALSE., add_nondiag_energy = .FALSE.
92 : REAL(KIND=dp) :: nondiag_energy_strength = -1.0_dp
93 : REAL(KIND=dp) :: broyden_beta = -1.0_dp, broyden_gamma = -1.0_dp, broyden_sigma = -1.0_dp
94 : REAL(KIND=dp) :: broyden_eta = -1.0_dp, broyden_omega = -1.0_dp, broyden_sigma_decrease = -1.0_dp
95 : REAL(KIND=dp) :: broyden_sigma_min = -1.0_dp
96 : LOGICAL :: broyden_forget_history = .FALSE., broyden_adaptive_sigma = .FALSE.
97 : LOGICAL :: broyden_enable_flip = .FALSE.
98 : END TYPE qs_ot_settings_type
99 :
100 : ! **************************************************************************************************
101 : TYPE qs_ot_type
102 : ! this sets the method to be used
103 : TYPE(qs_ot_settings_type) :: settings = qs_ot_settings_type()
104 : LOGICAL :: restricted = .FALSE.
105 :
106 : ! first part of the variables, for occupied subspace invariant optimisation
107 :
108 : ! add a preconditioner matrix. should be symmetric and positive definite
109 : ! the type of this matrix might change in the future
110 : TYPE(preconditioner_type), POINTER :: preconditioner => NULL()
111 :
112 : ! these will/might change during iterations
113 :
114 : ! OT / TOD
115 : TYPE(dbcsr_type), POINTER :: matrix_p => NULL()
116 : TYPE(dbcsr_type), POINTER :: matrix_r => NULL()
117 : TYPE(dbcsr_type), POINTER :: matrix_sinp => NULL()
118 : TYPE(dbcsr_type), POINTER :: matrix_cosp => NULL()
119 : TYPE(dbcsr_type), POINTER :: matrix_sinp_b => NULL()
120 : TYPE(dbcsr_type), POINTER :: matrix_cosp_b => NULL()
121 : TYPE(dbcsr_type), POINTER :: matrix_buf1 => NULL()
122 : TYPE(dbcsr_type), POINTER :: matrix_buf2 => NULL()
123 : TYPE(dbcsr_type), POINTER :: matrix_buf3 => NULL()
124 : TYPE(dbcsr_type), POINTER :: matrix_buf4 => NULL()
125 : TYPE(dbcsr_type), POINTER :: matrix_os => NULL()
126 : TYPE(dbcsr_type), POINTER :: matrix_buf1_ortho => NULL()
127 : TYPE(dbcsr_type), POINTER :: matrix_buf2_ortho => NULL()
128 :
129 : REAL(KIND=dp), DIMENSION(:), POINTER :: evals => NULL()
130 : REAL(KIND=dp), DIMENSION(:), POINTER :: dum => NULL()
131 :
132 : ! matrix os valid
133 : LOGICAL :: os_valid = .FALSE.
134 :
135 : ! for efficient/parallel writing to the blacs_matrix
136 : TYPE(mp_para_env_type), POINTER :: para_env => NULL()
137 : TYPE(cp_blacs_env_type), POINTER :: blacs_env => NULL()
138 :
139 : ! mo-like vectors
140 : TYPE(dbcsr_type), POINTER :: matrix_c0 => NULL(), matrix_sc0 => NULL(), matrix_psc0 => NULL()
141 :
142 : ! OT / IR
143 : TYPE(dbcsr_type), POINTER :: buf1_k_k_nosym => NULL(), buf2_k_k_nosym => NULL(), &
144 : buf3_k_k_nosym => NULL(), buf4_k_k_nosym => NULL(), &
145 : buf1_k_k_sym => NULL(), buf2_k_k_sym => NULL(), &
146 : buf3_k_k_sym => NULL(), buf4_k_k_sym => NULL(), &
147 : p_k_k_sym => NULL(), buf1_n_k => NULL(), buf1_n_k_dp => NULL()
148 :
149 : ! only here for the ease of programming. These will have to be supplied
150 : ! explicitly at all times
151 : TYPE(dbcsr_type), POINTER :: matrix_x => NULL(), matrix_sx => NULL(), matrix_gx => NULL()
152 : TYPE(dbcsr_type), POINTER :: matrix_dx => NULL(), matrix_gx_old => NULL()
153 :
154 : LOGICAL :: use_gx_old = .FALSE., use_dx = .FALSE.
155 :
156 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_h_e => NULL(), matrix_h_x => NULL()
157 :
158 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: ls_diis => NULL()
159 : REAL(KIND=dp), DIMENSION(:, :), POINTER :: lss_diis => NULL()
160 : REAL(KIND=dp), DIMENSION(:), POINTER :: c_diis => NULL()
161 : REAL(KIND=dp), DIMENSION(:), POINTER :: c_broy => NULL()
162 : REAL(KIND=dp), DIMENSION(:), POINTER :: energy_h => NULL()
163 : INTEGER, DIMENSION(:), POINTER :: ipivot => NULL()
164 :
165 : REAL(KIND=dp) :: ot_pos(53) = -1.0_dp, ot_energy(53) = -1.0_dp, ot_grad(53) = -1.0_dp ! HARD LIMIT FOR THE LS
166 : INTEGER :: line_search_left = -1, line_search_right = -1, line_search_mid = -1
167 : INTEGER :: line_search_count = -1
168 : LOGICAL :: line_search_might_be_done = .FALSE.
169 : REAL(KIND=dp) :: delta = -1.0_dp, gnorm = -1.0_dp, gnorm_old = -1.0_dp, etotal = -1.0_dp, gradient = -1.0_dp
170 : LOGICAL :: energy_only = .FALSE.
171 : INTEGER :: diis_iter = -1
172 : CHARACTER(LEN=8) :: OT_METHOD_FULL = ""
173 : INTEGER :: OT_count = -1
174 : REAL(KIND=dp) :: ds_min = -1.0_dp
175 : REAL(KIND=dp) :: broyden_adaptive_sigma = -1.0_dp
176 :
177 : LOGICAL :: do_taylor = .FALSE.
178 : INTEGER :: taylor_order = -1
179 : REAL(KIND=dp) :: largest_eval_upper_bound = -1.0_dp
180 :
181 : ! second part of the variables, if an explicit rotation is required as well
182 : TYPE(dbcsr_type), POINTER :: rot_mat_u => NULL() ! rotation matrix
183 : TYPE(dbcsr_type), POINTER :: rot_mat_x => NULL() ! antisymmetric matrix that parametrises rot_matrix_u
184 : TYPE(dbcsr_type), POINTER :: rot_mat_dedu => NULL() ! derivative of the total energy wrt to u
185 : TYPE(dbcsr_type), POINTER :: rot_mat_chc => NULL() ! for convencience, the matrix c^T H c
186 :
187 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rot_mat_h_e => NULL()
188 : TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rot_mat_h_x => NULL()
189 : TYPE(dbcsr_type), POINTER :: rot_mat_gx => NULL()
190 : TYPE(dbcsr_type), POINTER :: rot_mat_gx_old => NULL()
191 : TYPE(dbcsr_type), POINTER :: rot_mat_dx => NULL()
192 :
193 : REAL(KIND=dp), DIMENSION(:), POINTER :: rot_mat_evals => NULL()
194 : TYPE(dbcsr_type), POINTER :: rot_mat_evec_re => NULL()
195 : TYPE(dbcsr_type), POINTER :: rot_mat_evec_im => NULL()
196 :
197 : ! third part of the variables, if we need to optimize orbital energies
198 : REAL(KIND=dp), POINTER, DIMENSION(:) :: ener_x => NULL()
199 : REAL(KIND=dp), POINTER, DIMENSION(:) :: ener_dx => NULL()
200 : REAL(KIND=dp), POINTER, DIMENSION(:) :: ener_gx => NULL()
201 : REAL(KIND=dp), POINTER, DIMENSION(:) :: ener_gx_old => NULL()
202 : REAL(KIND=dp), POINTER, DIMENSION(:, :) :: ener_h_e => NULL()
203 : REAL(KIND=dp), POINTER, DIMENSION(:, :) :: ener_h_x => NULL()
204 : END TYPE qs_ot_type
205 :
206 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_ot_types'
207 :
208 : CONTAINS
209 :
210 : ! **************************************************************************************************
211 : !> \brief sets default values for the settings type
212 : !> \param settings ...
213 : !> \par History
214 : !> 10.2004 created [Joost VandeVondele]
215 : ! **************************************************************************************************
216 14459 : SUBROUTINE qs_ot_settings_init(settings)
217 : TYPE(qs_ot_settings_type) :: settings
218 :
219 14459 : settings%ot_method = "CG"
220 14459 : settings%ot_algorithm = "TOD"
221 14459 : settings%diis_m = 7
222 14459 : settings%preconditioner_name = "FULL_KINETIC"
223 14459 : settings%preconditioner_type = ot_precond_full_kinetic
224 14459 : settings%cholesky_type = cholesky_reduce
225 14459 : settings%precond_solver_name = "CHOLESKY_INVERSE"
226 14459 : settings%precond_solver_type = ot_precond_solver_inv_chol
227 14459 : settings%line_search_method = "2PNT"
228 14459 : settings%ds_min = 0.15_dp
229 14459 : settings%safer_diis = .TRUE.
230 14459 : settings%energy_gap = 0.2_dp
231 14459 : settings%eps_taylor = 1.0E-16_dp
232 14459 : settings%max_taylor = 4
233 14459 : settings%gold_target = 0.01_dp
234 14459 : settings%do_rotation = .FALSE.
235 14459 : settings%do_ener = .FALSE.
236 14459 : settings%irac_degree = 4
237 14459 : settings%max_irac = 50
238 14459 : settings%max_scf_diis = 0
239 14459 : settings%eps_irac = 1.0E-10_dp
240 14459 : settings%eps_irac_quick_exit = 1.0E-5_dp
241 14459 : settings%eps_irac_switch = 1.0E-2
242 14459 : settings%eps_irac_filter_matrix = 0.0_dp
243 14459 : settings%on_the_fly_loc = .FALSE.
244 14459 : settings%ortho_irac = "CHOL"
245 14459 : settings%ks = .TRUE.
246 14459 : settings%occupation_preconditioner = .FALSE.
247 14459 : settings%add_nondiag_energy = .FALSE.
248 14459 : settings%nondiag_energy_strength = 0.0_dp
249 :
250 14459 : END SUBROUTINE qs_ot_settings_init
251 :
252 : ! **************************************************************************************************
253 : !> \brief init matrices, needs c0 and sc0 so that c0*sc0=1
254 : !> \param qs_ot_env ...
255 : ! **************************************************************************************************
256 8368 : SUBROUTINE qs_ot_init(qs_ot_env)
257 : TYPE(qs_ot_type) :: qs_ot_env
258 :
259 451872 : qs_ot_env%OT_energy(:) = 0.0_dp
260 451872 : qs_ot_env%OT_pos(:) = 0.0_dp
261 451872 : qs_ot_env%OT_grad(:) = 0.0_dp
262 8368 : qs_ot_env%line_search_count = 0
263 :
264 8368 : qs_ot_env%energy_only = .FALSE.
265 8368 : qs_ot_env%gnorm_old = 1.0_dp
266 8368 : qs_ot_env%diis_iter = 0
267 8368 : qs_ot_env%ds_min = qs_ot_env%settings%ds_min
268 8368 : qs_ot_env%os_valid = .FALSE.
269 :
270 8368 : CALL dbcsr_set(qs_ot_env%matrix_gx, 0.0_dp)
271 :
272 8368 : IF (qs_ot_env%use_dx) &
273 3799 : CALL dbcsr_set(qs_ot_env%matrix_dx, 0.0_dp)
274 :
275 8368 : IF (qs_ot_env%use_gx_old) &
276 3799 : CALL dbcsr_set(qs_ot_env%matrix_gx_old, 0.0_dp)
277 :
278 8368 : IF (qs_ot_env%settings%do_rotation) THEN
279 212 : CALL dbcsr_set(qs_ot_env%rot_mat_gx, 0.0_dp)
280 212 : IF (qs_ot_env%use_dx) &
281 182 : CALL dbcsr_set(qs_ot_env%rot_mat_dx, 0.0_dp)
282 212 : IF (qs_ot_env%use_gx_old) &
283 182 : CALL dbcsr_set(qs_ot_env%rot_mat_gx_old, 0.0_dp)
284 : END IF
285 8368 : IF (qs_ot_env%settings%do_ener) THEN
286 0 : qs_ot_env%ener_gx(:) = 0.0_dp
287 0 : IF (qs_ot_env%use_dx) &
288 0 : qs_ot_env%ener_dx(:) = 0.0_dp
289 0 : IF (qs_ot_env%use_gx_old) &
290 0 : qs_ot_env%ener_gx_old(:) = 0.0_dp
291 : END IF
292 :
293 8368 : END SUBROUTINE qs_ot_init
294 :
295 : ! **************************************************************************************************
296 : !> \brief allocates the data in qs_ot_env, for a calculation with fm_struct_ref
297 : !> ortho_k allows for specifying an additional orthogonal subspace (i.e. c will
298 : !> be kept orthogonal provided c0 was, used in qs_ot_eigensolver)
299 : !> \param qs_ot_env ...
300 : !> \param matrix_s ...
301 : !> \param fm_struct_ref ...
302 : !> \param ortho_k ...
303 : ! **************************************************************************************************
304 8368 : SUBROUTINE qs_ot_allocate(qs_ot_env, matrix_s, fm_struct_ref, ortho_k)
305 : TYPE(qs_ot_type) :: qs_ot_env
306 : TYPE(dbcsr_type), POINTER :: matrix_s
307 : TYPE(cp_fm_struct_type), POINTER :: fm_struct_ref
308 : INTEGER, OPTIONAL :: ortho_k
309 :
310 : INTEGER :: i, k, m_diis, my_ortho_k, n, ncoef
311 : TYPE(cp_blacs_env_type), POINTER :: context
312 : TYPE(mp_para_env_type), POINTER :: para_env
313 :
314 8368 : CALL cite_reference(VandeVondele2003)
315 :
316 8368 : NULLIFY (qs_ot_env%preconditioner)
317 8368 : NULLIFY (qs_ot_env%matrix_psc0)
318 8368 : NULLIFY (qs_ot_env%para_env)
319 8368 : NULLIFY (qs_ot_env%blacs_env)
320 :
321 : CALL cp_fm_struct_get(fm_struct_ref, nrow_global=n, ncol_global=k, &
322 8368 : para_env=para_env, context=context)
323 :
324 8368 : qs_ot_env%para_env => para_env
325 8368 : qs_ot_env%blacs_env => context
326 8368 : CALL para_env%retain()
327 8368 : CALL context%retain()
328 :
329 8368 : IF (PRESENT(ortho_k)) THEN
330 534 : my_ortho_k = ortho_k
331 : ELSE
332 7834 : my_ortho_k = k
333 : END IF
334 :
335 8368 : m_diis = qs_ot_env%settings%diis_m
336 :
337 8368 : qs_ot_env%use_gx_old = .FALSE.
338 8368 : qs_ot_env%use_dx = .FALSE.
339 :
340 3799 : SELECT CASE (qs_ot_env%settings%ot_method)
341 : CASE ("SD")
342 : ! nothing
343 : CASE ("CG")
344 3799 : qs_ot_env%use_gx_old = .TRUE.
345 3799 : qs_ot_env%use_dx = .TRUE.
346 : CASE ("DIIS", "BROY")
347 4559 : IF (m_diis .LT. 1) CPABORT("m_diis less than one")
348 : CASE DEFAULT
349 8368 : CPABORT("Unknown option")
350 : END SELECT
351 :
352 8368 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS" .OR. &
353 : qs_ot_env%settings%ot_method .EQ. "BROY") THEN
354 18236 : ALLOCATE (qs_ot_env%ls_diis(m_diis + 1, m_diis + 1))
355 344017 : qs_ot_env%ls_diis = 0.0_dp
356 13677 : ALLOCATE (qs_ot_env%lss_diis(m_diis + 1, m_diis + 1))
357 13677 : ALLOCATE (qs_ot_env%c_diis(m_diis + 1))
358 13677 : ALLOCATE (qs_ot_env%c_broy(m_diis))
359 9118 : ALLOCATE (qs_ot_env%energy_h(m_diis))
360 13677 : ALLOCATE (qs_ot_env%ipivot(m_diis + 1))
361 : END IF
362 :
363 24898 : ALLOCATE (qs_ot_env%evals(k))
364 16530 : ALLOCATE (qs_ot_env%dum(k))
365 :
366 8368 : NULLIFY (qs_ot_env%matrix_os)
367 8368 : NULLIFY (qs_ot_env%matrix_buf1_ortho)
368 8368 : NULLIFY (qs_ot_env%matrix_buf2_ortho)
369 8368 : NULLIFY (qs_ot_env%matrix_p)
370 8368 : NULLIFY (qs_ot_env%matrix_r)
371 8368 : NULLIFY (qs_ot_env%matrix_sinp)
372 8368 : NULLIFY (qs_ot_env%matrix_cosp)
373 8368 : NULLIFY (qs_ot_env%matrix_sinp_b)
374 8368 : NULLIFY (qs_ot_env%matrix_cosp_b)
375 8368 : NULLIFY (qs_ot_env%matrix_buf1)
376 8368 : NULLIFY (qs_ot_env%matrix_buf2)
377 8368 : NULLIFY (qs_ot_env%matrix_buf3)
378 8368 : NULLIFY (qs_ot_env%matrix_buf4)
379 8368 : NULLIFY (qs_ot_env%matrix_c0)
380 8368 : NULLIFY (qs_ot_env%matrix_sc0)
381 8368 : NULLIFY (qs_ot_env%matrix_x)
382 8368 : NULLIFY (qs_ot_env%matrix_sx)
383 8368 : NULLIFY (qs_ot_env%matrix_gx)
384 8368 : NULLIFY (qs_ot_env%matrix_gx_old)
385 8368 : NULLIFY (qs_ot_env%matrix_dx)
386 8368 : NULLIFY (qs_ot_env%buf1_k_k_nosym)
387 8368 : NULLIFY (qs_ot_env%buf2_k_k_nosym)
388 8368 : NULLIFY (qs_ot_env%buf3_k_k_nosym)
389 8368 : NULLIFY (qs_ot_env%buf4_k_k_nosym)
390 8368 : NULLIFY (qs_ot_env%buf1_k_k_sym)
391 8368 : NULLIFY (qs_ot_env%buf2_k_k_sym)
392 8368 : NULLIFY (qs_ot_env%buf3_k_k_sym)
393 8368 : NULLIFY (qs_ot_env%buf4_k_k_sym)
394 8368 : NULLIFY (qs_ot_env%buf1_n_k)
395 8368 : NULLIFY (qs_ot_env%buf1_n_k_dp)
396 8368 : NULLIFY (qs_ot_env%p_k_k_sym)
397 :
398 : ! COMMON MATRICES
399 8368 : CALL dbcsr_init_p(qs_ot_env%matrix_c0)
400 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_c0, template=matrix_s, n=k, &
401 8368 : sym=dbcsr_type_no_symmetry)
402 :
403 8368 : CALL dbcsr_init_p(qs_ot_env%matrix_sc0)
404 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_sc0, template=matrix_s, n=my_ortho_k, &
405 8368 : sym=dbcsr_type_no_symmetry)
406 :
407 8368 : CALL dbcsr_init_p(qs_ot_env%matrix_x)
408 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_x, template=matrix_s, n=k, &
409 8368 : sym=dbcsr_type_no_symmetry)
410 :
411 8368 : CALL dbcsr_init_p(qs_ot_env%matrix_sx)
412 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_sx, template=matrix_s, n=k, &
413 8368 : sym=dbcsr_type_no_symmetry)
414 :
415 8368 : CALL dbcsr_init_p(qs_ot_env%matrix_gx)
416 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_gx, template=matrix_s, n=k, &
417 8368 : sym=dbcsr_type_no_symmetry)
418 :
419 8368 : IF (qs_ot_env%use_dx) THEN
420 3799 : CALL dbcsr_init_p(qs_ot_env%matrix_dx)
421 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_dx, template=matrix_s, n=k, &
422 3799 : sym=dbcsr_type_no_symmetry)
423 : END IF
424 :
425 8368 : IF (qs_ot_env%use_gx_old) THEN
426 3799 : CALL dbcsr_init_p(qs_ot_env%matrix_gx_old)
427 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_gx_old, template=matrix_s, n=k, &
428 3799 : sym=dbcsr_type_no_symmetry)
429 : END IF
430 :
431 7292 : SELECT CASE (qs_ot_env%settings%ot_algorithm)
432 : CASE ("TOD")
433 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_p)
434 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_p, template=matrix_s, m=k, n=k, &
435 7292 : sym=dbcsr_type_no_symmetry)
436 :
437 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_r)
438 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_r, template=matrix_s, m=k, n=k, &
439 7292 : sym=dbcsr_type_no_symmetry)
440 :
441 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_sinp)
442 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_sinp, template=matrix_s, m=k, n=k, &
443 7292 : sym=dbcsr_type_no_symmetry)
444 :
445 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_cosp)
446 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_cosp, template=matrix_s, m=k, n=k, &
447 7292 : sym=dbcsr_type_no_symmetry)
448 :
449 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_sinp_b)
450 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_sinp_b, template=matrix_s, m=k, n=k, &
451 7292 : sym=dbcsr_type_no_symmetry)
452 :
453 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_cosp_b)
454 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_cosp_b, template=matrix_s, m=k, n=k, &
455 7292 : sym=dbcsr_type_no_symmetry)
456 :
457 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_buf1)
458 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf1, template=matrix_s, m=k, n=k, &
459 7292 : sym=dbcsr_type_no_symmetry)
460 :
461 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_buf2)
462 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf2, template=matrix_s, m=k, n=k, &
463 7292 : sym=dbcsr_type_no_symmetry)
464 :
465 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_buf3)
466 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf3, template=matrix_s, m=k, n=k, &
467 7292 : sym=dbcsr_type_no_symmetry)
468 :
469 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_buf4)
470 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf4, template=matrix_s, m=k, n=k, &
471 7292 : sym=dbcsr_type_no_symmetry)
472 :
473 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_os)
474 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_os, template=matrix_s, m=my_ortho_k, n=my_ortho_k, &
475 7292 : sym=dbcsr_type_no_symmetry)
476 :
477 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_buf1_ortho)
478 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf1_ortho, template=matrix_s, m=my_ortho_k, n=k, &
479 7292 : sym=dbcsr_type_no_symmetry)
480 :
481 7292 : CALL dbcsr_init_p(qs_ot_env%matrix_buf2_ortho)
482 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf2_ortho, template=matrix_s, m=my_ortho_k, n=k, &
483 7292 : sym=dbcsr_type_no_symmetry)
484 :
485 : CASE ("REF")
486 1076 : CALL dbcsr_init_p(qs_ot_env%buf1_k_k_nosym)
487 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf1_k_k_nosym, template=matrix_s, m=k, n=k, &
488 1076 : sym=dbcsr_type_no_symmetry)
489 :
490 1076 : CALL dbcsr_init_p(qs_ot_env%buf2_k_k_nosym)
491 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf2_k_k_nosym, template=matrix_s, m=k, n=k, &
492 1076 : sym=dbcsr_type_no_symmetry)
493 :
494 1076 : CALL dbcsr_init_p(qs_ot_env%buf3_k_k_nosym)
495 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf3_k_k_nosym, template=matrix_s, m=k, n=k, &
496 1076 : sym=dbcsr_type_no_symmetry)
497 :
498 1076 : CALL dbcsr_init_p(qs_ot_env%buf4_k_k_nosym)
499 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf4_k_k_nosym, template=matrix_s, m=k, n=k, &
500 1076 : sym=dbcsr_type_no_symmetry)
501 :
502 : ! It claims to be symmetric but to avoid dbcsr confusion nonsymmetric is kept
503 1076 : CALL dbcsr_init_p(qs_ot_env%buf1_k_k_sym)
504 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf1_k_k_sym, template=matrix_s, m=k, n=k, &
505 1076 : sym=dbcsr_type_no_symmetry)
506 :
507 1076 : CALL dbcsr_init_p(qs_ot_env%buf2_k_k_sym)
508 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf2_k_k_sym, template=matrix_s, m=k, n=k, &
509 1076 : sym=dbcsr_type_no_symmetry)
510 :
511 1076 : CALL dbcsr_init_p(qs_ot_env%buf3_k_k_sym)
512 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf3_k_k_sym, template=matrix_s, m=k, n=k, &
513 1076 : sym=dbcsr_type_no_symmetry)
514 : !
515 1076 : CALL dbcsr_init_p(qs_ot_env%buf4_k_k_sym)
516 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%buf4_k_k_sym, template=matrix_s, m=k, n=k, &
517 1076 : sym=dbcsr_type_no_symmetry)
518 : !
519 1076 : CALL dbcsr_init_p(qs_ot_env%p_k_k_sym)
520 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%p_k_k_sym, template=matrix_s, m=k, n=k, &
521 1076 : sym=dbcsr_type_no_symmetry)
522 : !
523 1076 : CALL dbcsr_init_p(qs_ot_env%buf1_n_k)
524 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%buf1_n_k, template=matrix_s, n=k, &
525 1076 : sym=dbcsr_type_no_symmetry)
526 : !
527 1076 : CALL dbcsr_init_p(qs_ot_env%matrix_buf1)
528 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%matrix_buf1, template=matrix_s, m=k, n=k, &
529 9444 : sym=dbcsr_type_no_symmetry)
530 :
531 : END SELECT
532 :
533 8368 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS" .OR. &
534 : qs_ot_env%settings%ot_method .EQ. "BROY") THEN
535 4559 : NULLIFY (qs_ot_env%matrix_h_e)
536 4559 : NULLIFY (qs_ot_env%matrix_h_x)
537 4559 : CALL dbcsr_allocate_matrix_set(qs_ot_env%matrix_h_e, m_diis)
538 4559 : CALL dbcsr_allocate_matrix_set(qs_ot_env%matrix_h_x, m_diis)
539 37062 : DO i = 1, m_diis
540 32503 : CALL dbcsr_init_p(qs_ot_env%matrix_h_x(i)%matrix)
541 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_h_x(i)%matrix, template=matrix_s, n=k, &
542 32503 : sym=dbcsr_type_no_symmetry)
543 :
544 32503 : CALL dbcsr_init_p(qs_ot_env%matrix_h_e(i)%matrix)
545 : CALL cp_dbcsr_m_by_n_from_row_template(qs_ot_env%matrix_h_e(i)%matrix, template=matrix_s, n=k, &
546 40871 : sym=dbcsr_type_no_symmetry)
547 : END DO
548 : END IF
549 :
550 8368 : NULLIFY (qs_ot_env%rot_mat_u, qs_ot_env%rot_mat_x, qs_ot_env%rot_mat_h_e, qs_ot_env%rot_mat_h_x, &
551 8368 : qs_ot_env%rot_mat_gx, qs_ot_env%rot_mat_gx_old, qs_ot_env%rot_mat_dx, &
552 8368 : qs_ot_env%rot_mat_evals, qs_ot_env%rot_mat_dedu, qs_ot_env%rot_mat_chc, &
553 8368 : qs_ot_env%rot_mat_evec_re, qs_ot_env%rot_mat_evec_im)
554 :
555 8368 : IF (qs_ot_env%settings%do_rotation) THEN
556 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_u)
557 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_u, template=matrix_s, m=k, n=k, &
558 212 : sym=dbcsr_type_no_symmetry)
559 :
560 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_x)
561 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_x, template=matrix_s, m=k, n=k, &
562 212 : sym=dbcsr_type_no_symmetry)
563 :
564 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_dedu)
565 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_dedu, template=matrix_s, m=k, n=k, &
566 212 : sym=dbcsr_type_no_symmetry)
567 :
568 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_chc)
569 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_chc, template=matrix_s, m=k, n=k, &
570 212 : sym=dbcsr_type_no_symmetry)
571 :
572 212 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS") THEN
573 30 : CALL dbcsr_allocate_matrix_set(qs_ot_env%rot_mat_h_e, m_diis)
574 30 : CALL dbcsr_allocate_matrix_set(qs_ot_env%rot_mat_h_x, m_diis)
575 240 : DO i = 1, m_diis
576 210 : CALL dbcsr_init_p(qs_ot_env%rot_mat_h_e(i)%matrix)
577 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_h_e(i)%matrix, template=matrix_s, m=k, n=k, &
578 210 : sym=dbcsr_type_no_symmetry)
579 :
580 210 : CALL dbcsr_init_p(qs_ot_env%rot_mat_h_x(i)%matrix)
581 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_h_x(i)%matrix, template=matrix_s, m=k, n=k, &
582 240 : sym=dbcsr_type_no_symmetry)
583 : END DO
584 : END IF
585 :
586 632 : ALLOCATE (qs_ot_env%rot_mat_evals(k))
587 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_evec_re)
588 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_evec_re, template=matrix_s, m=k, n=k, &
589 212 : sym=dbcsr_type_no_symmetry)
590 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_evec_im)
591 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_evec_im, template=matrix_s, m=k, n=k, &
592 212 : sym=dbcsr_type_no_symmetry)
593 :
594 212 : CALL dbcsr_init_p(qs_ot_env%rot_mat_gx)
595 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_gx, template=matrix_s, m=k, n=k, &
596 212 : sym=dbcsr_type_no_symmetry)
597 :
598 212 : IF (qs_ot_env%use_gx_old) THEN
599 182 : CALL dbcsr_init_p(qs_ot_env%rot_mat_gx_old)
600 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_gx_old, template=matrix_s, m=k, n=k, &
601 182 : sym=dbcsr_type_no_symmetry)
602 : END IF
603 :
604 212 : IF (qs_ot_env%use_dx) THEN
605 182 : CALL dbcsr_init_p(qs_ot_env%rot_mat_dx)
606 : CALL cp_dbcsr_m_by_n_from_template(qs_ot_env%rot_mat_dx, template=matrix_s, m=k, n=k, &
607 182 : sym=dbcsr_type_no_symmetry)
608 : END IF
609 :
610 : END IF
611 :
612 8368 : IF (qs_ot_env%settings%do_ener) THEN
613 0 : ncoef = k
614 0 : ALLOCATE (qs_ot_env%ener_x(ncoef))
615 :
616 0 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS") THEN
617 0 : ALLOCATE (qs_ot_env%ener_h_e(m_diis, ncoef))
618 0 : ALLOCATE (qs_ot_env%ener_h_x(m_diis, ncoef))
619 : END IF
620 :
621 0 : ALLOCATE (qs_ot_env%ener_gx(ncoef))
622 :
623 0 : IF (qs_ot_env%use_gx_old) THEN
624 0 : ALLOCATE (qs_ot_env%ener_gx_old(ncoef))
625 : END IF
626 :
627 0 : IF (qs_ot_env%use_dx) THEN
628 0 : ALLOCATE (qs_ot_env%ener_dx(ncoef))
629 0 : qs_ot_env%ener_dx = 0.0_dp
630 : END IF
631 : END IF
632 :
633 8368 : END SUBROUTINE qs_ot_allocate
634 :
635 : ! **************************************************************************************************
636 : !> \brief deallocates data
637 : !> \param qs_ot_env ...
638 : ! **************************************************************************************************
639 8368 : SUBROUTINE qs_ot_destroy(qs_ot_env)
640 : TYPE(qs_ot_type) :: qs_ot_env
641 :
642 8368 : CALL mp_para_env_release(qs_ot_env%para_env)
643 8368 : CALL cp_blacs_env_release(qs_ot_env%blacs_env)
644 :
645 8368 : DEALLOCATE (qs_ot_env%evals)
646 8368 : DEALLOCATE (qs_ot_env%dum)
647 :
648 8368 : IF (ASSOCIATED(qs_ot_env%matrix_os)) CALL dbcsr_release_p(qs_ot_env%matrix_os)
649 8368 : IF (ASSOCIATED(qs_ot_env%matrix_p)) CALL dbcsr_release_p(qs_ot_env%matrix_p)
650 8368 : IF (ASSOCIATED(qs_ot_env%matrix_cosp)) CALL dbcsr_release_p(qs_ot_env%matrix_cosp)
651 8368 : IF (ASSOCIATED(qs_ot_env%matrix_sinp)) CALL dbcsr_release_p(qs_ot_env%matrix_sinp)
652 8368 : IF (ASSOCIATED(qs_ot_env%matrix_r)) CALL dbcsr_release_p(qs_ot_env%matrix_r)
653 8368 : IF (ASSOCIATED(qs_ot_env%matrix_cosp_b)) CALL dbcsr_release_p(qs_ot_env%matrix_cosp_b)
654 8368 : IF (ASSOCIATED(qs_ot_env%matrix_sinp_b)) CALL dbcsr_release_p(qs_ot_env%matrix_sinp_b)
655 8368 : IF (ASSOCIATED(qs_ot_env%matrix_buf1)) CALL dbcsr_release_p(qs_ot_env%matrix_buf1)
656 8368 : IF (ASSOCIATED(qs_ot_env%matrix_buf2)) CALL dbcsr_release_p(qs_ot_env%matrix_buf2)
657 8368 : IF (ASSOCIATED(qs_ot_env%matrix_buf3)) CALL dbcsr_release_p(qs_ot_env%matrix_buf3)
658 8368 : IF (ASSOCIATED(qs_ot_env%matrix_buf4)) CALL dbcsr_release_p(qs_ot_env%matrix_buf4)
659 8368 : IF (ASSOCIATED(qs_ot_env%matrix_buf1_ortho)) CALL dbcsr_release_p(qs_ot_env%matrix_buf1_ortho)
660 8368 : IF (ASSOCIATED(qs_ot_env%matrix_buf2_ortho)) CALL dbcsr_release_p(qs_ot_env%matrix_buf2_ortho)
661 8368 : IF (ASSOCIATED(qs_ot_env%matrix_c0)) CALL dbcsr_release_p(qs_ot_env%matrix_c0)
662 8368 : IF (ASSOCIATED(qs_ot_env%matrix_sc0)) CALL dbcsr_release_p(qs_ot_env%matrix_sc0)
663 8368 : IF (ASSOCIATED(qs_ot_env%matrix_psc0)) CALL dbcsr_release_p(qs_ot_env%matrix_psc0)
664 8368 : IF (ASSOCIATED(qs_ot_env%matrix_x)) CALL dbcsr_release_p(qs_ot_env%matrix_x)
665 8368 : IF (ASSOCIATED(qs_ot_env%matrix_sx)) CALL dbcsr_release_p(qs_ot_env%matrix_sx)
666 8368 : IF (ASSOCIATED(qs_ot_env%matrix_gx)) CALL dbcsr_release_p(qs_ot_env%matrix_gx)
667 8368 : IF (ASSOCIATED(qs_ot_env%matrix_dx)) CALL dbcsr_release_p(qs_ot_env%matrix_dx)
668 8368 : IF (ASSOCIATED(qs_ot_env%matrix_gx_old)) CALL dbcsr_release_p(qs_ot_env%matrix_gx_old)
669 8368 : IF (ASSOCIATED(qs_ot_env%buf1_k_k_nosym)) CALL dbcsr_release_p(qs_ot_env%buf1_k_k_nosym)
670 8368 : IF (ASSOCIATED(qs_ot_env%buf2_k_k_nosym)) CALL dbcsr_release_p(qs_ot_env%buf2_k_k_nosym)
671 8368 : IF (ASSOCIATED(qs_ot_env%buf3_k_k_nosym)) CALL dbcsr_release_p(qs_ot_env%buf3_k_k_nosym)
672 8368 : IF (ASSOCIATED(qs_ot_env%buf4_k_k_nosym)) CALL dbcsr_release_p(qs_ot_env%buf4_k_k_nosym)
673 8368 : IF (ASSOCIATED(qs_ot_env%p_k_k_sym)) CALL dbcsr_release_p(qs_ot_env%p_k_k_sym)
674 8368 : IF (ASSOCIATED(qs_ot_env%buf1_k_k_sym)) CALL dbcsr_release_p(qs_ot_env%buf1_k_k_sym)
675 8368 : IF (ASSOCIATED(qs_ot_env%buf2_k_k_sym)) CALL dbcsr_release_p(qs_ot_env%buf2_k_k_sym)
676 8368 : IF (ASSOCIATED(qs_ot_env%buf3_k_k_sym)) CALL dbcsr_release_p(qs_ot_env%buf3_k_k_sym)
677 8368 : IF (ASSOCIATED(qs_ot_env%buf4_k_k_sym)) CALL dbcsr_release_p(qs_ot_env%buf4_k_k_sym)
678 8368 : IF (ASSOCIATED(qs_ot_env%buf1_n_k)) CALL dbcsr_release_p(qs_ot_env%buf1_n_k)
679 8368 : IF (ASSOCIATED(qs_ot_env%buf1_n_k_dp)) CALL dbcsr_release_p(qs_ot_env%buf1_n_k_dp)
680 :
681 8368 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS" .OR. &
682 : qs_ot_env%settings%ot_method .EQ. "BROY") THEN
683 4559 : CALL dbcsr_deallocate_matrix_set(qs_ot_env%matrix_h_x)
684 4559 : CALL dbcsr_deallocate_matrix_set(qs_ot_env%matrix_h_e)
685 4559 : DEALLOCATE (qs_ot_env%ls_diis)
686 4559 : DEALLOCATE (qs_ot_env%lss_diis)
687 4559 : DEALLOCATE (qs_ot_env%c_diis)
688 4559 : DEALLOCATE (qs_ot_env%c_broy)
689 4559 : DEALLOCATE (qs_ot_env%energy_h)
690 4559 : DEALLOCATE (qs_ot_env%ipivot)
691 : END IF
692 :
693 8368 : IF (qs_ot_env%settings%do_rotation) THEN
694 :
695 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_u)) CALL dbcsr_release_p(qs_ot_env%rot_mat_u)
696 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_x)) CALL dbcsr_release_p(qs_ot_env%rot_mat_x)
697 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_dedu)) CALL dbcsr_release_p(qs_ot_env%rot_mat_dedu)
698 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_chc)) CALL dbcsr_release_p(qs_ot_env%rot_mat_chc)
699 :
700 212 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS") THEN
701 30 : CALL dbcsr_deallocate_matrix_set(qs_ot_env%rot_mat_h_x)
702 30 : CALL dbcsr_deallocate_matrix_set(qs_ot_env%rot_mat_h_e)
703 : END IF
704 :
705 212 : DEALLOCATE (qs_ot_env%rot_mat_evals)
706 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_evec_re)) CALL dbcsr_release_p(qs_ot_env%rot_mat_evec_re)
707 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_evec_im)) CALL dbcsr_release_p(qs_ot_env%rot_mat_evec_im)
708 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_gx)) CALL dbcsr_release_p(qs_ot_env%rot_mat_gx)
709 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_gx_old)) CALL dbcsr_release_p(qs_ot_env%rot_mat_gx_old)
710 212 : IF (ASSOCIATED(qs_ot_env%rot_mat_dx)) CALL dbcsr_release_p(qs_ot_env%rot_mat_dx)
711 : END IF
712 :
713 8368 : IF (qs_ot_env%settings%do_ener) THEN
714 0 : DEALLOCATE (qs_ot_env%ener_x)
715 0 : DEALLOCATE (qs_ot_env%ener_gx)
716 0 : IF (qs_ot_env%settings%ot_method .EQ. "DIIS") THEN
717 0 : DEALLOCATE (qs_ot_env%ener_h_x)
718 0 : DEALLOCATE (qs_ot_env%ener_h_e)
719 : END IF
720 0 : IF (qs_ot_env%use_dx) THEN
721 0 : DEALLOCATE (qs_ot_env%ener_dx)
722 : END IF
723 0 : IF (qs_ot_env%use_gx_old) THEN
724 0 : DEALLOCATE (qs_ot_env%ener_gx_old)
725 : END IF
726 : END IF
727 :
728 8368 : END SUBROUTINE qs_ot_destroy
729 :
730 : ! **************************************************************************************************
731 : !> \brief ...
732 : !> \param settings ...
733 : !> \param ot_section ...
734 : !> \param output_unit ...
735 : ! **************************************************************************************************
736 33045 : SUBROUTINE ot_readwrite_input(settings, ot_section, output_unit)
737 : TYPE(qs_ot_settings_type) :: settings
738 : TYPE(section_vals_type), POINTER :: ot_section
739 : INTEGER, INTENT(IN) :: output_unit
740 :
741 : CHARACTER(len=*), PARAMETER :: routineN = 'ot_readwrite_input'
742 :
743 : INTEGER :: handle, ls_method, ot_algorithm, &
744 : ot_method, ot_ortho_irac
745 :
746 6609 : CALL timeset(routineN, handle)
747 :
748 : ! choose algorithm
749 6609 : CALL section_vals_val_get(ot_section, "ALGORITHM", i_val=ot_algorithm)
750 6039 : SELECT CASE (ot_algorithm)
751 : CASE (ot_algo_taylor_or_diag)
752 6039 : settings%ot_algorithm = "TOD"
753 : CASE (ot_algo_irac)
754 570 : CALL cite_reference(Weber2008)
755 570 : settings%ot_algorithm = "REF"
756 : CASE DEFAULT
757 6609 : CPABORT("Value unknown")
758 : END SELECT
759 :
760 : ! irac input
761 6609 : CALL section_vals_val_get(ot_section, "IRAC_DEGREE", i_val=settings%irac_degree)
762 6609 : IF (settings%irac_degree < 2 .OR. settings%irac_degree > 4) THEN
763 0 : CPABORT("READ OT IRAC_DEGREE: Value unknown")
764 : END IF
765 6609 : CALL section_vals_val_get(ot_section, "MAX_IRAC", i_val=settings%max_irac)
766 6609 : IF (settings%max_irac < 1) THEN
767 0 : CPABORT("READ OT MAX_IRAC: VALUE MUST BE GREATER THAN ZERO")
768 : END IF
769 6609 : CALL section_vals_val_get(ot_section, "EPS_IRAC_FILTER_MATRIX", r_val=settings%eps_irac_filter_matrix)
770 6609 : CALL section_vals_val_get(ot_section, "EPS_IRAC", r_val=settings%eps_irac)
771 6609 : IF (settings%eps_irac < 0.0_dp) THEN
772 0 : CPABORT("READ OT EPS_IRAC: VALUE MUST BE GREATER THAN ZERO")
773 : END IF
774 6609 : CALL section_vals_val_get(ot_section, "EPS_IRAC_QUICK_EXIT", r_val=settings%eps_irac_quick_exit)
775 6609 : IF (settings%eps_irac_quick_exit < 0.0_dp) THEN
776 0 : CPABORT("READ OT EPS_IRAC_QUICK_EXIT: VALUE MUST BE GREATER THAN ZERO")
777 : END IF
778 :
779 6609 : CALL section_vals_val_get(ot_section, "EPS_IRAC_SWITCH", r_val=settings%eps_irac_switch)
780 6609 : IF (settings%eps_irac_switch < 0.0_dp) THEN
781 0 : CPABORT("READ OT EPS_IRAC_SWITCH: VALUE MUST BE GREATER THAN ZERO")
782 : END IF
783 :
784 6609 : CALL section_vals_val_get(ot_section, "ORTHO_IRAC", i_val=ot_ortho_irac)
785 6525 : SELECT CASE (ot_ortho_irac)
786 : CASE (ot_chol_irac)
787 6525 : settings%ortho_irac = "CHOL"
788 : CASE (ot_poly_irac)
789 34 : settings%ortho_irac = "POLY"
790 : CASE (ot_lwdn_irac)
791 50 : settings%ortho_irac = "LWDN"
792 : CASE DEFAULT
793 6609 : CPABORT("READ OT ORTHO_IRAC: Value unknown")
794 : END SELECT
795 :
796 6609 : CALL section_vals_val_get(ot_section, "ON_THE_FLY_LOC", l_val=settings%on_the_fly_loc)
797 :
798 6609 : CALL section_vals_val_get(ot_section, "MINIMIZER", i_val=ot_method)
799 : ! overwrite input if ot_state is set
800 6609 : IF (settings%ot_state == 1) THEN
801 4 : ot_method = ot_mini_diis
802 : END IF
803 : ! compatibility
804 10 : SELECT CASE (ot_method)
805 : CASE (ot_mini_sd)
806 10 : settings%ot_method = "SD"
807 : CASE (ot_mini_cg)
808 2318 : settings%ot_method = "CG"
809 : CASE (ot_mini_diis)
810 4267 : settings%ot_method = "DIIS"
811 4267 : CALL section_vals_val_get(ot_section, "N_HISTORY_VEC", i_val=settings%diis_m)
812 : CASE (ot_mini_broyden)
813 14 : CALL section_vals_val_get(ot_section, "N_HISTORY_VEC", i_val=settings%diis_m)
814 14 : CALL section_vals_val_get(ot_section, "BROYDEN_BETA", r_val=settings%broyden_beta)
815 14 : CALL section_vals_val_get(ot_section, "BROYDEN_GAMMA", r_val=settings%broyden_gamma)
816 14 : CALL section_vals_val_get(ot_section, "BROYDEN_SIGMA", r_val=settings%broyden_sigma)
817 14 : CALL section_vals_val_get(ot_section, "BROYDEN_ETA", r_val=settings%broyden_eta)
818 14 : CALL section_vals_val_get(ot_section, "BROYDEN_OMEGA", r_val=settings%broyden_omega)
819 14 : CALL section_vals_val_get(ot_section, "BROYDEN_SIGMA_DECREASE", r_val=settings%broyden_sigma_decrease)
820 14 : CALL section_vals_val_get(ot_section, "BROYDEN_SIGMA_MIN", r_val=settings%broyden_sigma_min)
821 14 : CALL section_vals_val_get(ot_section, "BROYDEN_FORGET_HISTORY", l_val=settings%broyden_forget_history)
822 14 : CALL section_vals_val_get(ot_section, "BROYDEN_ADAPTIVE_SIGMA", l_val=settings%broyden_adaptive_sigma)
823 14 : CALL section_vals_val_get(ot_section, "BROYDEN_ENABLE_FLIP", l_val=settings%broyden_enable_flip)
824 14 : settings%ot_method = "BROY"
825 : CASE DEFAULT
826 6609 : CPABORT("READ OTSCF MINIMIZER: Value unknown")
827 : END SELECT
828 6609 : CALL section_vals_val_get(ot_section, "SAFER_DIIS", l_val=settings%safer_diis)
829 6609 : CALL section_vals_val_get(ot_section, "LINESEARCH", i_val=ls_method)
830 2 : SELECT CASE (ls_method)
831 : CASE (ls_none)
832 2 : settings%line_search_method = "NONE"
833 : CASE (ls_2pnt)
834 6541 : settings%line_search_method = "2PNT"
835 : CASE (ls_3pnt)
836 42 : settings%line_search_method = "3PNT"
837 : CASE (ls_adapt)
838 12 : settings%line_search_method = "ADPT"
839 : CASE (ls_gold)
840 12 : settings%line_search_method = "GOLD"
841 12 : CALL section_vals_val_get(ot_section, "GOLD_TARGET", r_val=settings%gold_target)
842 : CASE DEFAULT
843 6609 : CPABORT("READ OTSCF LS: Value unknown")
844 : END SELECT
845 :
846 6609 : CALL section_vals_val_get(ot_section, "PRECOND_SOLVER", i_val=settings%precond_solver_type)
847 13196 : SELECT CASE (settings%precond_solver_type)
848 : CASE (ot_precond_solver_default)
849 6587 : settings%precond_solver_name = "DEFAULT"
850 : CASE (ot_precond_solver_inv_chol)
851 16 : settings%precond_solver_name = "INVERSE_CHOLESKY"
852 : CASE (ot_precond_solver_direct)
853 0 : settings%precond_solver_name = "DIRECT"
854 : CASE (ot_precond_solver_update)
855 6 : settings%precond_solver_name = "INVERSE_UPDATE"
856 : CASE DEFAULT
857 6609 : CPABORT("READ OTSCF SOLVER: Value unknown")
858 : END SELECT
859 :
860 6609 : CALL section_vals_val_get(ot_section, "MAX_SCF_DIIS", i_val=settings%max_scf_diis)
861 :
862 : !If these values are negative we will set them "optimal" for a given precondtioner below
863 6609 : CALL section_vals_val_get(ot_section, "STEPSIZE", r_val=settings%ds_min)
864 6609 : CALL section_vals_val_get(ot_section, "ENERGY_GAP", r_val=settings%energy_gap)
865 :
866 6609 : CALL section_vals_val_get(ot_section, "PRECONDITIONER", i_val=settings%preconditioner_type)
867 7311 : SELECT CASE (settings%preconditioner_type)
868 : CASE (ot_precond_none)
869 702 : settings%preconditioner_name = "NONE"
870 702 : IF (settings%ds_min < 0.0_dp) settings%ds_min = 0.15_dp
871 702 : IF (settings%energy_gap < 0.0_dp) settings%energy_gap = 0.2_dp
872 : CASE (ot_precond_full_single)
873 30 : settings%preconditioner_name = "FULL_SINGLE"
874 30 : IF (settings%ds_min < 0.0_dp) settings%ds_min = 0.15_dp
875 30 : IF (settings%energy_gap < 0.0_dp) settings%energy_gap = 0.2_dp
876 : CASE (ot_precond_full_single_inverse)
877 2756 : settings%preconditioner_name = "FULL_SINGLE_INVERSE"
878 2756 : IF (settings%ds_min < 0.0_dp) settings%ds_min = 0.08_dp
879 2756 : IF (settings%energy_gap < 0.0_dp) settings%energy_gap = 0.08_dp
880 : CASE (ot_precond_full_all)
881 1864 : settings%preconditioner_name = "FULL_ALL"
882 1864 : IF (settings%ds_min < 0.0_dp) settings%ds_min = 0.15_dp
883 1864 : IF (settings%energy_gap < 0.0_dp) settings%energy_gap = 0.08_dp
884 : CASE (ot_precond_full_kinetic)
885 1221 : settings%preconditioner_name = "FULL_KINETIC"
886 1221 : IF (settings%ds_min < 0.0_dp) settings%ds_min = 0.15_dp
887 1221 : IF (settings%energy_gap < 0.0_dp) settings%energy_gap = 0.2_dp
888 : CASE (ot_precond_s_inverse)
889 36 : settings%preconditioner_name = "FULL_S_INVERSE"
890 36 : IF (settings%ds_min < 0.0_dp) settings%ds_min = 0.15_dp
891 36 : IF (settings%energy_gap < 0.0_dp) settings%energy_gap = 0.2_dp
892 : CASE DEFAULT
893 6609 : CPABORT("READ OTSCF PRECONDITIONER: Value unknown")
894 : END SELECT
895 6609 : CALL section_vals_val_get(ot_section, "CHOLESKY", i_val=settings%cholesky_type)
896 6609 : CALL section_vals_val_get(ot_section, "EPS_TAYLOR", r_val=settings%eps_taylor)
897 6609 : CALL section_vals_val_get(ot_section, "MAX_TAYLOR", i_val=settings%max_taylor)
898 6609 : CALL section_vals_val_get(ot_section, "ROTATION", l_val=settings%do_rotation)
899 6609 : CALL section_vals_val_get(ot_section, "ENERGIES", l_val=settings%do_ener)
900 : CALL section_vals_val_get(ot_section, "OCCUPATION_PRECONDITIONER", &
901 6609 : l_val=settings%occupation_preconditioner)
902 6609 : CALL section_vals_val_get(ot_section, "NONDIAG_ENERGY", l_val=settings%add_nondiag_energy)
903 : CALL section_vals_val_get(ot_section, "NONDIAG_ENERGY_STRENGTH", &
904 6609 : r_val=settings%nondiag_energy_strength)
905 : ! not yet fully implemented
906 6609 : CPASSERT(.NOT. settings%do_ener)
907 :
908 : ! write OT output
909 :
910 6609 : IF (output_unit > 0) THEN
911 3421 : WRITE (output_unit, '(/,A)') " ----------------------------------- OT ---------------------------------------"
912 3421 : IF (settings%do_rotation) THEN
913 100 : WRITE (output_unit, '(A)') " Allowing for rotations "
914 : END IF
915 3421 : IF (settings%do_ener) THEN
916 0 : WRITE (output_unit, '(A,L2)') " Optimizing orbital energies "
917 : END IF
918 5 : SELECT CASE (settings%OT_METHOD)
919 : CASE ("SD")
920 5 : WRITE (output_unit, '(A)') " Minimizer : SD : steepest descent"
921 : CASE ("CG")
922 1188 : WRITE (output_unit, '(A)') " Minimizer : CG : conjugate gradient"
923 : CASE ("DIIS")
924 2221 : WRITE (output_unit, '(A)') " Minimizer : DIIS : direct inversion"
925 2221 : WRITE (output_unit, '(A)') " in the iterative subspace"
926 2221 : WRITE (output_unit, '(A,I3,A)') " using ", settings%diis_m, " DIIS vectors"
927 2221 : IF (settings%safer_diis) THEN
928 2221 : WRITE (output_unit, '(A,I3,A)') " safer DIIS on"
929 : ELSE
930 0 : WRITE (output_unit, '(A,I3,A)') " safer DIIS off"
931 : END IF
932 : CASE ("BROY")
933 7 : WRITE (output_unit, '(A)') " Minimizer : BROYDEN : Broyden "
934 7 : WRITE (output_unit, '(A,F16.8)') " BETA : ", settings%broyden_beta
935 7 : WRITE (output_unit, '(A,F16.8)') " GAMMA : ", settings%broyden_gamma
936 7 : WRITE (output_unit, '(A,F16.8)') " SIGMA : ", settings%broyden_sigma
937 7 : WRITE (output_unit, '(A,I3,A)') " using : - ", &
938 14 : settings%diis_m, " BROYDEN vectors"
939 : CASE DEFAULT
940 3421 : WRITE (output_unit, '(3A)') " Minimizer : ", settings%OT_METHOD, " : UNKNOWN"
941 : END SELECT
942 15 : SELECT CASE (settings%preconditioner_name)
943 : CASE ("FULL_SINGLE")
944 15 : WRITE (output_unit, '(A)') " Preconditioner : FULL_SINGLE : diagonalization based"
945 : CASE ("FULL_SINGLE_INVERSE")
946 1445 : WRITE (output_unit, '(A,/,A)') " Preconditioner : FULL_SINGLE_INVERSE : inversion of ", &
947 2890 : " H + eS - 2*(Sc)(c^T*H*c+const)(Sc)^T"
948 : CASE ("FULL_ALL")
949 962 : WRITE (output_unit, '(A)') " Preconditioner : FULL_ALL : diagonalization, state selective"
950 : CASE ("FULL_KINETIC")
951 630 : WRITE (output_unit, '(A)') " Preconditioner : FULL_KINETIC : inversion of T + eS"
952 : CASE ("FULL_S_INVERSE")
953 18 : WRITE (output_unit, '(A)') " Preconditioner : FULL_S_INVERSE : cholesky inversion of S"
954 : CASE ("SPARSE_DIAG")
955 : WRITE (output_unit, '(A)') &
956 0 : " Preconditioner : SPARSE_DIAG : diagonal atomic block diagonalization"
957 : CASE ("SPARSE_KINETIC")
958 0 : WRITE (output_unit, '(A)') " Preconditioner : SPARSE_KINETIC : sparse linear solver for T + eS"
959 : CASE ("NONE")
960 351 : WRITE (output_unit, '(A)') " Preconditioner : NONE"
961 : CASE DEFAULT
962 3421 : WRITE (output_unit, '(3A)') " Preconditioner : ", settings%preconditioner_name, " : UNKNOWN"
963 : END SELECT
964 :
965 3421 : WRITE (output_unit, '(A)') " Precond_solver : "//TRIM(settings%precond_solver_name)
966 :
967 3421 : IF (settings%OT_METHOD .EQ. "SD" .OR. settings%OT_METHOD .EQ. "CG") THEN
968 1162 : SELECT CASE (settings%line_search_method)
969 : CASE ("2PNT")
970 1162 : WRITE (output_unit, '(A)') " Line search : 2PNT : 2 energies, one gradient"
971 : CASE ("3PNT")
972 19 : WRITE (output_unit, '(A)') " Line search : 3PNT : 3 energies"
973 : CASE ("GOLD")
974 6 : WRITE (output_unit, '(A)') " Line search : GOLD : bracketing and golden section search"
975 6 : WRITE (output_unit, '(A,F14.8)') " target rel accuracy : ", settings%gold_target
976 : CASE ("NONE")
977 1 : WRITE (output_unit, '(A)') " Line search : NONE"
978 : CASE DEFAULT
979 3421 : WRITE (output_unit, '(3A)') " Line search : ", settings%line_search_method, " : UNKNOWN"
980 : END SELECT
981 : END IF
982 3421 : WRITE (output_unit, '(A,F14.8,T49,A,F14.8)') " stepsize :", settings%ds_min, &
983 6842 : " energy_gap :", settings%energy_gap
984 3421 : IF (settings%ot_algorithm .EQ. 'TOD') THEN
985 3118 : WRITE (output_unit, '(A,E14.5,T49,A,I14)') " eps_taylor :", settings%eps_taylor, &
986 6236 : " max_taylor :", settings%max_taylor
987 : END IF
988 3421 : IF (settings%ot_algorithm .EQ. 'REF') THEN
989 303 : WRITE (output_unit, '(A,1X,A,T49,A,I14)') " ortho_irac :", settings%ortho_irac, &
990 606 : " irac_degree :", settings%irac_degree
991 303 : WRITE (output_unit, '(A,I14,T49,A,E14.5)') " max_irac :", settings%max_irac, &
992 606 : " eps_irac :", settings%eps_irac
993 303 : WRITE (output_unit, '(A,E14.5,T49,A,E10.3)') " eps_irac_switch:", settings%eps_irac_switch, &
994 606 : " eps_irac_quick_exit:", settings%eps_irac_quick_exit
995 303 : WRITE (output_unit, '(A,L2)') " on_the_fly_loc :", settings%on_the_fly_loc
996 : END IF
997 3421 : WRITE (output_unit, '(A)') " ----------------------------------- OT ---------------------------------------"
998 : WRITE (UNIT=output_unit, &
999 : FMT="(/,T3,A,T12,A,T31,A,T39,A,T59,A,T75,A,/,T3,A)") &
1000 3421 : "Step", "Update method", "Time", "Convergence", "Total energy", "Change", &
1001 6842 : REPEAT("-", 78)
1002 : END IF
1003 :
1004 6609 : CALL timestop(handle)
1005 :
1006 6609 : END SUBROUTINE ot_readwrite_input
1007 :
1008 0 : END MODULE qs_ot_types
|
