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