(git:374b731)
Loading...
Searching...
No Matches
preconditioner_makes.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief computes preconditioners, and implements methods to apply them
10!> currently used in qs_ot
11!> \par History
12!> - [UB] 2009-05-13 Adding stable approximate inverse (full and sparse)
13!> \author Joost VandeVondele (09.2002)
14! **************************************************************************************************
15MODULE preconditioner_makes
16 USE arnoldi_api, ONLY: arnoldi_data_type,&
17 arnoldi_ev,&
18 deallocate_arnoldi_data,&
19 get_selected_ritz_val,&
20 get_selected_ritz_vector,&
21 set_arnoldi_initial_vector,&
22 setup_arnoldi_data
23 USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
24 cp_dbcsr_m_by_n_from_template,&
25 cp_dbcsr_sm_fm_multiply,&
26 cp_fm_to_dbcsr_row_template
27 USE cp_fm_basic_linalg, ONLY: cp_fm_column_scale,&
28 cp_fm_triangular_invert,&
29 cp_fm_triangular_multiply,&
30 cp_fm_upper_to_full
31 USE cp_fm_cholesky, ONLY: cp_fm_cholesky_decompose,&
32 cp_fm_cholesky_reduce,&
33 cp_fm_cholesky_restore
34 USE cp_fm_diag, ONLY: choose_eigv_solver
35 USE cp_fm_struct, ONLY: cp_fm_struct_create,&
36 cp_fm_struct_release,&
37 cp_fm_struct_type
38 USE cp_fm_types, ONLY: cp_fm_create,&
39 cp_fm_get_diag,&
40 cp_fm_get_info,&
41 cp_fm_release,&
42 cp_fm_to_fm,&
43 cp_fm_type
44 USE dbcsr_api, ONLY: &
45 dbcsr_add, dbcsr_add_on_diag, dbcsr_copy, dbcsr_create, dbcsr_get_info, dbcsr_multiply, &
46 dbcsr_p_type, dbcsr_release, dbcsr_type, dbcsr_type_symmetric
47 USE input_constants, ONLY: &
48 cholesky_inverse, cholesky_reduce, ot_precond_full_all, ot_precond_full_kinetic, &
49 ot_precond_full_single, ot_precond_full_single_inverse, ot_precond_s_inverse, &
50 ot_precond_solver_default, ot_precond_solver_inv_chol
51 USE kinds, ONLY: dp
52 USE parallel_gemm_api, ONLY: parallel_gemm
53 USE preconditioner_types, ONLY: preconditioner_type
54#include "./base/base_uses.f90"
55
56 IMPLICIT NONE
57
58 PRIVATE
59
60 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'preconditioner_makes'
61
62 PUBLIC :: make_preconditioner_matrix
63
64CONTAINS
65
66! **************************************************************************************************
67!> \brief ...
68!> \param preconditioner_env ...
69!> \param matrix_h ...
70!> \param matrix_s ...
71!> \param matrix_t ...
72!> \param mo_coeff ...
73!> \param energy_homo ...
74!> \param eigenvalues_ot ...
75!> \param energy_gap ...
76!> \param my_solver_type ...
77! **************************************************************************************************
78 SUBROUTINE make_preconditioner_matrix(preconditioner_env, matrix_h, matrix_s, matrix_t, mo_coeff, &
79 energy_homo, eigenvalues_ot, energy_gap, &
80 my_solver_type)
81 TYPE(preconditioner_type) :: preconditioner_env
82 TYPE(dbcsr_type), POINTER :: matrix_h
83 TYPE(dbcsr_type), OPTIONAL, POINTER :: matrix_s, matrix_t
84 TYPE(cp_fm_type), INTENT(IN) :: mo_coeff
85 REAL(kind=dp) :: energy_homo
86 REAL(kind=dp), DIMENSION(:) :: eigenvalues_ot
87 REAL(kind=dp) :: energy_gap
88 INTEGER :: my_solver_type
89
90 INTEGER :: precon_type
91
92 precon_type = preconditioner_env%in_use
93 SELECT CASE (precon_type)
94 CASE (ot_precond_full_single)
95 IF (my_solver_type .NE. ot_precond_solver_default) &
96 cpabort("Only PRECOND_SOLVER DEFAULT for the moment")
97 IF (PRESENT(matrix_s)) THEN
98 CALL make_full_single(preconditioner_env, preconditioner_env%fm, &
99 matrix_h, matrix_s, energy_homo, energy_gap)
100 ELSE
101 CALL make_full_single_ortho(preconditioner_env, preconditioner_env%fm, &
102 matrix_h, energy_homo, energy_gap)
103 END IF
104
105 CASE (ot_precond_s_inverse)
106 IF (my_solver_type .EQ. ot_precond_solver_default) my_solver_type = ot_precond_solver_inv_chol
107 IF (.NOT. PRESENT(matrix_s)) &
108 cpabort("Type for S=1 not implemented")
109 CALL make_full_s_inverse(preconditioner_env, matrix_s)
110
111 CASE (ot_precond_full_kinetic)
112 IF (my_solver_type .EQ. ot_precond_solver_default) my_solver_type = ot_precond_solver_inv_chol
113 IF (.NOT. (PRESENT(matrix_s) .AND. PRESENT(matrix_t))) &
114 cpabort("Type for S=1 not implemented")
115 CALL make_full_kinetic(preconditioner_env, matrix_t, matrix_s, energy_gap)
116 CASE (ot_precond_full_single_inverse)
117 IF (my_solver_type .EQ. ot_precond_solver_default) my_solver_type = ot_precond_solver_inv_chol
118 CALL make_full_single_inverse(preconditioner_env, mo_coeff, matrix_h, energy_gap, &
119 matrix_s=matrix_s)
120 CASE (ot_precond_full_all)
121 IF (my_solver_type .NE. ot_precond_solver_default) THEN
122 cpabort("Only PRECOND_SOLVER DEFAULT for the moment")
123 END IF
124 IF (PRESENT(matrix_s)) THEN
125 CALL make_full_all(preconditioner_env, mo_coeff, matrix_h, matrix_s, &
126 eigenvalues_ot, energy_gap)
127 ELSE
128 CALL make_full_all_ortho(preconditioner_env, mo_coeff, matrix_h, &
129 eigenvalues_ot, energy_gap)
130 END IF
131
132 CASE DEFAULT
133 cpabort("Type not implemented")
134 END SELECT
135
136 END SUBROUTINE make_preconditioner_matrix
137
138! **************************************************************************************************
139!> \brief Simply takes the overlap matrix as preconditioner
140!> \param preconditioner_env ...
141!> \param matrix_s ...
142! **************************************************************************************************
143 SUBROUTINE make_full_s_inverse(preconditioner_env, matrix_s)
144 TYPE(preconditioner_type) :: preconditioner_env
145 TYPE(dbcsr_type), POINTER :: matrix_s
146
147 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_s_inverse'
148
149 INTEGER :: handle
150
151 CALL timeset(routinen, handle)
152
153 cpassert(ASSOCIATED(matrix_s))
154
155 IF (.NOT. ASSOCIATED(preconditioner_env%sparse_matrix)) THEN
156 ALLOCATE (preconditioner_env%sparse_matrix)
157 END IF
158 CALL dbcsr_copy(preconditioner_env%sparse_matrix, matrix_s, name="full_kinetic")
159
160 CALL timestop(handle)
161
162 END SUBROUTINE make_full_s_inverse
163
164! **************************************************************************************************
165!> \brief kinetic matrix+shift*overlap as preconditioner. Cheap but could
166!> be better
167!> \param preconditioner_env ...
168!> \param matrix_t ...
169!> \param matrix_s ...
170!> \param energy_gap ...
171! **************************************************************************************************
172 SUBROUTINE make_full_kinetic(preconditioner_env, matrix_t, matrix_s, &
173 energy_gap)
174 TYPE(preconditioner_type) :: preconditioner_env
175 TYPE(dbcsr_type), POINTER :: matrix_t, matrix_s
176 REAL(kind=dp) :: energy_gap
177
178 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_kinetic'
179
180 INTEGER :: handle
181 REAL(kind=dp) :: shift
182
183 CALL timeset(routinen, handle)
184
185 cpassert(ASSOCIATED(matrix_t))
186 cpassert(ASSOCIATED(matrix_s))
187
188 IF (.NOT. ASSOCIATED(preconditioner_env%sparse_matrix)) THEN
189 ALLOCATE (preconditioner_env%sparse_matrix)
190 END IF
191 CALL dbcsr_copy(preconditioner_env%sparse_matrix, matrix_t, name="full_kinetic")
192
193 shift = max(0.0_dp, energy_gap)
194
195 CALL dbcsr_add(preconditioner_env%sparse_matrix, matrix_s, &
196 alpha_scalar=1.0_dp, beta_scalar=shift)
197
198 CALL timestop(handle)
199
200 END SUBROUTINE make_full_kinetic
201
202! **************************************************************************************************
203!> \brief full_single_preconditioner
204!> \param preconditioner_env ...
205!> \param fm ...
206!> \param matrix_h ...
207!> \param matrix_s ...
208!> \param energy_homo ...
209!> \param energy_gap ...
210! **************************************************************************************************
211 SUBROUTINE make_full_single(preconditioner_env, fm, matrix_h, matrix_s, &
212 energy_homo, energy_gap)
213 TYPE(preconditioner_type) :: preconditioner_env
214 TYPE(cp_fm_type), POINTER :: fm
215 TYPE(dbcsr_type), POINTER :: matrix_h, matrix_s
216 REAL(kind=dp) :: energy_homo, energy_gap
217
218 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_single'
219
220 INTEGER :: handle, i, n
221 REAL(kind=dp), DIMENSION(:), POINTER :: evals
222 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
223 TYPE(cp_fm_type) :: fm_h, fm_s
224
225 CALL timeset(routinen, handle)
226
227 NULLIFY (fm_struct_tmp, evals)
228
229 IF (ASSOCIATED(fm)) THEN
230 CALL cp_fm_release(fm)
231 DEALLOCATE (fm)
232 NULLIFY (fm)
233 END IF
234 CALL dbcsr_get_info(matrix_h, nfullrows_total=n)
235 ALLOCATE (evals(n))
236
237 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=n, ncol_global=n, &
238 context=preconditioner_env%ctxt, &
239 para_env=preconditioner_env%para_env)
240 ALLOCATE (fm)
241 CALL cp_fm_create(fm, fm_struct_tmp, name="preconditioner")
242 CALL cp_fm_create(fm_h, fm_struct_tmp, name="fm_h")
243 CALL cp_fm_create(fm_s, fm_struct_tmp, name="fm_s")
244 CALL cp_fm_struct_release(fm_struct_tmp)
245
246 CALL copy_dbcsr_to_fm(matrix_h, fm_h)
247 CALL copy_dbcsr_to_fm(matrix_s, fm_s)
248 CALL cp_fm_cholesky_decompose(fm_s)
249
250 SELECT CASE (preconditioner_env%cholesky_use)
251 CASE (cholesky_inverse)
252! if cho inverse
253 CALL cp_fm_triangular_invert(fm_s)
254 CALL cp_fm_upper_to_full(fm_h, fm)
255
256 CALL cp_fm_triangular_multiply(fm_s, fm_h, side="R", transpose_tr=.false., &
257 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=n, alpha=1.0_dp)
258 CALL cp_fm_triangular_multiply(fm_s, fm_h, side="L", transpose_tr=.true., &
259 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=n, alpha=1.0_dp)
260 CASE (cholesky_reduce)
261 CALL cp_fm_cholesky_reduce(fm_h, fm_s)
262 CASE DEFAULT
263 cpabort("cholesky type not implemented")
264 END SELECT
265
266 CALL choose_eigv_solver(fm_h, fm, evals)
267
268 SELECT CASE (preconditioner_env%cholesky_use)
269 CASE (cholesky_inverse)
270 CALL cp_fm_triangular_multiply(fm_s, fm, side="L", transpose_tr=.false., &
271 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=n, alpha=1.0_dp)
272 DO i = 1, n
273 evals(i) = 1.0_dp/max(evals(i) - energy_homo, energy_gap)
274 END DO
275 CALL cp_fm_to_fm(fm, fm_h)
276 CASE (cholesky_reduce)
277 CALL cp_fm_cholesky_restore(fm, n, fm_s, fm_h, "SOLVE")
278 DO i = 1, n
279 evals(i) = 1.0_dp/max(evals(i) - energy_homo, energy_gap)
280 END DO
281 CALL cp_fm_to_fm(fm_h, fm)
282 END SELECT
283
284 CALL cp_fm_column_scale(fm, evals)
285 CALL parallel_gemm('N', 'T', n, n, n, 1.0_dp, fm, fm_h, 0.0_dp, fm_s)
286 CALL cp_fm_to_fm(fm_s, fm)
287
288 DEALLOCATE (evals)
289 CALL cp_fm_release(fm_h)
290 CALL cp_fm_release(fm_s)
291
292 CALL timestop(handle)
293
294 END SUBROUTINE make_full_single
295
296! **************************************************************************************************
297!> \brief full single in the orthonormal basis
298!> \param preconditioner_env ...
299!> \param fm ...
300!> \param matrix_h ...
301!> \param energy_homo ...
302!> \param energy_gap ...
303! **************************************************************************************************
304 SUBROUTINE make_full_single_ortho(preconditioner_env, fm, matrix_h, &
305 energy_homo, energy_gap)
306 TYPE(preconditioner_type) :: preconditioner_env
307 TYPE(cp_fm_type), POINTER :: fm
308 TYPE(dbcsr_type), POINTER :: matrix_h
309 REAL(kind=dp) :: energy_homo, energy_gap
310
311 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_single_ortho'
312
313 INTEGER :: handle, i, n
314 REAL(kind=dp), DIMENSION(:), POINTER :: evals
315 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
316 TYPE(cp_fm_type) :: fm_h, fm_s
317
318 CALL timeset(routinen, handle)
319 NULLIFY (fm_struct_tmp, evals)
320
321 IF (ASSOCIATED(fm)) THEN
322 CALL cp_fm_release(fm)
323 DEALLOCATE (fm)
324 NULLIFY (fm)
325 END IF
326 CALL dbcsr_get_info(matrix_h, nfullrows_total=n)
327 ALLOCATE (evals(n))
328
329 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=n, ncol_global=n, &
330 context=preconditioner_env%ctxt, &
331 para_env=preconditioner_env%para_env)
332 ALLOCATE (fm)
333 CALL cp_fm_create(fm, fm_struct_tmp, name="preconditioner")
334 CALL cp_fm_create(fm_h, fm_struct_tmp, name="fm_h")
335 CALL cp_fm_create(fm_s, fm_struct_tmp, name="fm_s")
336 CALL cp_fm_struct_release(fm_struct_tmp)
337
338 CALL copy_dbcsr_to_fm(matrix_h, fm_h)
339
340 CALL choose_eigv_solver(fm_h, fm, evals)
341 DO i = 1, n
342 evals(i) = 1.0_dp/max(evals(i) - energy_homo, energy_gap)
343 END DO
344 CALL cp_fm_to_fm(fm, fm_h)
345 CALL cp_fm_column_scale(fm, evals)
346 CALL parallel_gemm('N', 'T', n, n, n, 1.0_dp, fm, fm_h, 0.0_dp, fm_s)
347 CALL cp_fm_to_fm(fm_s, fm)
348
349 DEALLOCATE (evals)
350 CALL cp_fm_release(fm_h)
351 CALL cp_fm_release(fm_s)
352
353 CALL timestop(handle)
354
355 END SUBROUTINE make_full_single_ortho
356
357! **************************************************************************************************
358!> \brief generates a state by state preconditioner based on the full hamiltonian matrix
359!> \param preconditioner_env ...
360!> \param matrix_c0 ...
361!> \param matrix_h ...
362!> \param matrix_s ...
363!> \param c0_evals ...
364!> \param energy_gap should be a slight underestimate of the physical energy gap for almost all systems
365!> the c0 are already ritz states of (h,s)
366!> \par History
367!> 10.2006 made more stable [Joost VandeVondele]
368!> \note
369!> includes error estimate on the hamiltonian matrix to result in a stable preconditioner
370!> a preconditioner for each eigenstate i is generated by keeping the factorized form
371!> U diag( something i ) U^T. It is important to only precondition in the subspace orthogonal to c0.
372!> not only is it the only part that matters, it also simplifies the computation of
373!> the lagrangian multipliers in the OT minimization (i.e. if the c0 here is different
374!> from the c0 used in the OT setup, there will be a bug).
375! **************************************************************************************************
376 SUBROUTINE make_full_all(preconditioner_env, matrix_c0, matrix_h, matrix_s, c0_evals, energy_gap)
377 TYPE(preconditioner_type) :: preconditioner_env
378 TYPE(cp_fm_type), INTENT(IN) :: matrix_c0
379 TYPE(dbcsr_type), POINTER :: matrix_h, matrix_s
380 REAL(kind=dp), DIMENSION(:) :: c0_evals
381 REAL(kind=dp) :: energy_gap
382
383 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_all'
384 REAL(kind=dp), PARAMETER :: fudge_factor = 0.25_dp, &
385 lambda_base = 10.0_dp
386
387 INTEGER :: handle, k, n
388 REAL(kind=dp) :: error_estimate, lambda
389 REAL(kind=dp), DIMENSION(:), POINTER :: diag, norms, shifted_evals
390 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
391 TYPE(cp_fm_type) :: matrix_hc0, matrix_left, matrix_s1, &
392 matrix_s2, matrix_sc0, matrix_shc0, &
393 matrix_tmp, ortho
394 TYPE(cp_fm_type), POINTER :: matrix_pre
395
396 CALL timeset(routinen, handle)
397
398 IF (ASSOCIATED(preconditioner_env%fm)) THEN
399 CALL cp_fm_release(preconditioner_env%fm)
400 DEALLOCATE (preconditioner_env%fm)
401 NULLIFY (preconditioner_env%fm)
402 END IF
403 CALL cp_fm_get_info(matrix_c0, nrow_global=n, ncol_global=k)
404 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=n, ncol_global=n, &
405 context=preconditioner_env%ctxt, &
406 para_env=preconditioner_env%para_env)
407 ALLOCATE (preconditioner_env%fm)
408 CALL cp_fm_create(preconditioner_env%fm, fm_struct_tmp, name="preconditioner_env%fm")
409 CALL cp_fm_create(ortho, fm_struct_tmp, name="ortho")
410 CALL cp_fm_create(matrix_tmp, fm_struct_tmp, name="matrix_tmp")
411 CALL cp_fm_struct_release(fm_struct_tmp)
412 ALLOCATE (preconditioner_env%full_evals(n))
413 ALLOCATE (preconditioner_env%occ_evals(k))
414
415 ! 0) cholesky decompose the overlap matrix, if this fails the basis is singular,
416 ! more than EPS_DEFAULT
417 CALL copy_dbcsr_to_fm(matrix_s, ortho)
418 CALL cp_fm_cholesky_decompose(ortho)
419! if cho inverse
420 IF (preconditioner_env%cholesky_use == cholesky_inverse) THEN
421 CALL cp_fm_triangular_invert(ortho)
422 END IF
423 ! 1) Construct a new H matrix, which has the current C0 as eigenvectors,
424 ! possibly shifted by an amount lambda,
425 ! and the same spectrum as the original H matrix in the space orthogonal to the C0
426 ! with P=C0 C0 ^ T
427 ! (1 - PS)^T H (1-PS) + (PS)^T (H - lambda S ) (PS)
428 ! we exploit that the C0 are already the ritz states of H
429 CALL cp_fm_create(matrix_sc0, matrix_c0%matrix_struct, name="sc0")
430 CALL cp_dbcsr_sm_fm_multiply(matrix_s, matrix_c0, matrix_sc0, k)
431 CALL cp_fm_create(matrix_hc0, matrix_c0%matrix_struct, name="hc0")
432 CALL cp_dbcsr_sm_fm_multiply(matrix_h, matrix_c0, matrix_hc0, k)
433
434 ! An aside, try to estimate the error on the ritz values, we'll need it later on
435 CALL cp_fm_create(matrix_shc0, matrix_c0%matrix_struct, name="shc0")
436
437 SELECT CASE (preconditioner_env%cholesky_use)
438 CASE (cholesky_inverse)
439! if cho inverse
440 CALL cp_fm_to_fm(matrix_hc0, matrix_shc0)
441 CALL cp_fm_triangular_multiply(ortho, matrix_shc0, side="L", transpose_tr=.true., &
442 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=k, alpha=1.0_dp)
443 CASE (cholesky_reduce)
444 CALL cp_fm_cholesky_restore(matrix_hc0, k, ortho, matrix_shc0, "SOLVE", transa="T")
445 CASE DEFAULT
446 cpabort("cholesky type not implemented")
447 END SELECT
448 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=k, ncol_global=k, &
449 context=preconditioner_env%ctxt, &
450 para_env=preconditioner_env%para_env)
451 CALL cp_fm_create(matrix_s1, fm_struct_tmp, name="matrix_s1")
452 CALL cp_fm_struct_release(fm_struct_tmp)
453 ! since we only use diagonal elements this is a bit of a waste
454 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, matrix_shc0, matrix_shc0, 0.0_dp, matrix_s1)
455 ALLOCATE (diag(k))
456 CALL cp_fm_get_diag(matrix_s1, diag)
457 error_estimate = maxval(sqrt(abs(diag - c0_evals**2)))
458 DEALLOCATE (diag)
459 CALL cp_fm_release(matrix_s1)
460 CALL cp_fm_release(matrix_shc0)
461 ! we'll only use the energy gap, if our estimate of the error on the eigenvalues
462 ! is small enough. A large error combined with a small energy gap would otherwise lead to
463 ! an aggressive but bad preconditioner. Only when the error is small (MD), we can precondition
464 ! aggressively
465 preconditioner_env%energy_gap = max(energy_gap, error_estimate*fudge_factor)
466 CALL copy_dbcsr_to_fm(matrix_h, matrix_tmp)
467 matrix_pre => preconditioner_env%fm
468 CALL cp_fm_upper_to_full(matrix_tmp, matrix_pre)
469 ! tmp = H ( 1 - PS )
470 CALL parallel_gemm('N', 'T', n, n, k, -1.0_dp, matrix_hc0, matrix_sc0, 1.0_dp, matrix_tmp)
471
472 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=k, ncol_global=n, &
473 context=preconditioner_env%ctxt, &
474 para_env=preconditioner_env%para_env)
475 CALL cp_fm_create(matrix_left, fm_struct_tmp, name="matrix_left")
476 CALL cp_fm_struct_release(fm_struct_tmp)
477 CALL parallel_gemm('T', 'N', k, n, n, 1.0_dp, matrix_c0, matrix_tmp, 0.0_dp, matrix_left)
478 ! tmp = (1 - PS)^T H (1-PS)
479 CALL parallel_gemm('N', 'N', n, n, k, -1.0_dp, matrix_sc0, matrix_left, 1.0_dp, matrix_tmp)
480 CALL cp_fm_release(matrix_left)
481
482 ALLOCATE (shifted_evals(k))
483 lambda = lambda_base + error_estimate
484 shifted_evals = c0_evals - lambda
485 CALL cp_fm_to_fm(matrix_sc0, matrix_hc0)
486 CALL cp_fm_column_scale(matrix_hc0, shifted_evals)
487 CALL parallel_gemm('N', 'T', n, n, k, 1.0_dp, matrix_hc0, matrix_sc0, 1.0_dp, matrix_tmp)
488
489 ! 2) diagonalize this operator
490 SELECT CASE (preconditioner_env%cholesky_use)
491 CASE (cholesky_inverse)
492 CALL cp_fm_triangular_multiply(ortho, matrix_tmp, side="R", transpose_tr=.false., &
493 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=n, alpha=1.0_dp)
494 CALL cp_fm_triangular_multiply(ortho, matrix_tmp, side="L", transpose_tr=.true., &
495 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=n, alpha=1.0_dp)
496 CASE (cholesky_reduce)
497 CALL cp_fm_cholesky_reduce(matrix_tmp, ortho)
498 END SELECT
499 CALL choose_eigv_solver(matrix_tmp, matrix_pre, preconditioner_env%full_evals)
500 SELECT CASE (preconditioner_env%cholesky_use)
501 CASE (cholesky_inverse)
502 CALL cp_fm_triangular_multiply(ortho, matrix_pre, side="L", transpose_tr=.false., &
503 invert_tr=.false., uplo_tr="U", n_rows=n, n_cols=n, alpha=1.0_dp)
504 CALL cp_fm_to_fm(matrix_pre, matrix_tmp)
505 CASE (cholesky_reduce)
506 CALL cp_fm_cholesky_restore(matrix_pre, n, ortho, matrix_tmp, "SOLVE")
507 CALL cp_fm_to_fm(matrix_tmp, matrix_pre)
508 END SELECT
509
510 ! test that the subspace remained conserved
511 IF (.false.) THEN
512 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=k, ncol_global=k, &
513 context=preconditioner_env%ctxt, &
514 para_env=preconditioner_env%para_env)
515 CALL cp_fm_create(matrix_s1, fm_struct_tmp, name="matrix_s1")
516 CALL cp_fm_create(matrix_s2, fm_struct_tmp, name="matrix_s2")
517 CALL cp_fm_struct_release(fm_struct_tmp)
518 ALLOCATE (norms(k))
519 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, matrix_sc0, matrix_tmp, 0.0_dp, matrix_s1)
520 CALL choose_eigv_solver(matrix_s1, matrix_s2, norms)
521 WRITE (*, *) "matrix norm deviation (should be close to zero): ", maxval(abs(abs(norms) - 1.0_dp))
522 DEALLOCATE (norms)
523 CALL cp_fm_release(matrix_s1)
524 CALL cp_fm_release(matrix_s2)
525 END IF
526
527 ! 3) replace the lowest k evals and evecs with what they should be
528 preconditioner_env%occ_evals = c0_evals
529 ! notice, this choice causes the preconditioner to be constant when applied to sc0 (see apply_full_all)
530 preconditioner_env%full_evals(1:k) = c0_evals
531 CALL cp_fm_to_fm(matrix_c0, matrix_pre, k, 1, 1)
532
533 CALL cp_fm_release(matrix_sc0)
534 CALL cp_fm_release(matrix_hc0)
535 CALL cp_fm_release(ortho)
536 CALL cp_fm_release(matrix_tmp)
537 DEALLOCATE (shifted_evals)
538 CALL timestop(handle)
539
540 END SUBROUTINE make_full_all
541
542! **************************************************************************************************
543!> \brief full all in the orthonormal basis
544!> \param preconditioner_env ...
545!> \param matrix_c0 ...
546!> \param matrix_h ...
547!> \param c0_evals ...
548!> \param energy_gap ...
549! **************************************************************************************************
550 SUBROUTINE make_full_all_ortho(preconditioner_env, matrix_c0, matrix_h, c0_evals, energy_gap)
551
552 TYPE(preconditioner_type) :: preconditioner_env
553 TYPE(cp_fm_type), INTENT(IN) :: matrix_c0
554 TYPE(dbcsr_type), POINTER :: matrix_h
555 REAL(kind=dp), DIMENSION(:) :: c0_evals
556 REAL(kind=dp) :: energy_gap
557
558 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_all_ortho'
559 REAL(kind=dp), PARAMETER :: fudge_factor = 0.25_dp, &
560 lambda_base = 10.0_dp
561
562 INTEGER :: handle, k, n
563 REAL(kind=dp) :: error_estimate, lambda
564 REAL(kind=dp), DIMENSION(:), POINTER :: diag, norms, shifted_evals
565 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
566 TYPE(cp_fm_type) :: matrix_hc0, matrix_left, matrix_s1, &
567 matrix_s2, matrix_sc0, matrix_tmp
568 TYPE(cp_fm_type), POINTER :: matrix_pre
569
570 CALL timeset(routinen, handle)
571
572 IF (ASSOCIATED(preconditioner_env%fm)) THEN
573 CALL cp_fm_release(preconditioner_env%fm)
574 DEALLOCATE (preconditioner_env%fm)
575 NULLIFY (preconditioner_env%fm)
576 END IF
577 CALL cp_fm_get_info(matrix_c0, nrow_global=n, ncol_global=k)
578 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=n, ncol_global=n, &
579 context=preconditioner_env%ctxt, &
580 para_env=preconditioner_env%para_env)
581 ALLOCATE (preconditioner_env%fm)
582 CALL cp_fm_create(preconditioner_env%fm, fm_struct_tmp, name="preconditioner_env%fm")
583 CALL cp_fm_create(matrix_tmp, fm_struct_tmp, name="matrix_tmp")
584 CALL cp_fm_struct_release(fm_struct_tmp)
585 ALLOCATE (preconditioner_env%full_evals(n))
586 ALLOCATE (preconditioner_env%occ_evals(k))
587
588 ! 1) Construct a new H matrix, which has the current C0 as eigenvectors,
589 ! possibly shifted by an amount lambda,
590 ! and the same spectrum as the original H matrix in the space orthogonal to the C0
591 ! with P=C0 C0 ^ T
592 ! (1 - PS)^T H (1-PS) + (PS)^T (H - lambda S ) (PS)
593 ! we exploit that the C0 are already the ritz states of H
594 CALL cp_fm_create(matrix_sc0, matrix_c0%matrix_struct, name="sc0")
595 CALL cp_fm_to_fm(matrix_c0, matrix_sc0)
596 CALL cp_fm_create(matrix_hc0, matrix_c0%matrix_struct, name="hc0")
597 CALL cp_dbcsr_sm_fm_multiply(matrix_h, matrix_c0, matrix_hc0, k)
598
599 ! An aside, try to estimate the error on the ritz values, we'll need it later on
600 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=k, ncol_global=k, &
601 context=preconditioner_env%ctxt, &
602 para_env=preconditioner_env%para_env)
603 CALL cp_fm_create(matrix_s1, fm_struct_tmp, name="matrix_s1")
604 CALL cp_fm_struct_release(fm_struct_tmp)
605 ! since we only use diagonal elements this is a bit of a waste
606 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, matrix_hc0, matrix_hc0, 0.0_dp, matrix_s1)
607 ALLOCATE (diag(k))
608 CALL cp_fm_get_diag(matrix_s1, diag)
609 error_estimate = maxval(sqrt(abs(diag - c0_evals**2)))
610 DEALLOCATE (diag)
611 CALL cp_fm_release(matrix_s1)
612 ! we'll only use the energy gap, if our estimate of the error on the eigenvalues
613 ! is small enough. A large error combined with a small energy gap would otherwise lead to
614 ! an aggressive but bad preconditioner. Only when the error is small (MD), we can precondition
615 ! aggressively
616 preconditioner_env%energy_gap = max(energy_gap, error_estimate*fudge_factor)
617
618 matrix_pre => preconditioner_env%fm
619 CALL copy_dbcsr_to_fm(matrix_h, matrix_tmp)
620 CALL cp_fm_upper_to_full(matrix_tmp, matrix_pre)
621 ! tmp = H ( 1 - PS )
622 CALL parallel_gemm('N', 'T', n, n, k, -1.0_dp, matrix_hc0, matrix_sc0, 1.0_dp, matrix_tmp)
623
624 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=k, ncol_global=n, &
625 context=preconditioner_env%ctxt, &
626 para_env=preconditioner_env%para_env)
627 CALL cp_fm_create(matrix_left, fm_struct_tmp, name="matrix_left")
628 CALL cp_fm_struct_release(fm_struct_tmp)
629 CALL parallel_gemm('T', 'N', k, n, n, 1.0_dp, matrix_c0, matrix_tmp, 0.0_dp, matrix_left)
630 ! tmp = (1 - PS)^T H (1-PS)
631 CALL parallel_gemm('N', 'N', n, n, k, -1.0_dp, matrix_sc0, matrix_left, 1.0_dp, matrix_tmp)
632 CALL cp_fm_release(matrix_left)
633
634 ALLOCATE (shifted_evals(k))
635 lambda = lambda_base + error_estimate
636 shifted_evals = c0_evals - lambda
637 CALL cp_fm_to_fm(matrix_sc0, matrix_hc0)
638 CALL cp_fm_column_scale(matrix_hc0, shifted_evals)
639 CALL parallel_gemm('N', 'T', n, n, k, 1.0_dp, matrix_hc0, matrix_sc0, 1.0_dp, matrix_tmp)
640
641 ! 2) diagonalize this operator
642 CALL choose_eigv_solver(matrix_tmp, matrix_pre, preconditioner_env%full_evals)
643
644 ! test that the subspace remained conserved
645 IF (.false.) THEN
646 CALL cp_fm_to_fm(matrix_pre, matrix_tmp)
647 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=k, ncol_global=k, &
648 context=preconditioner_env%ctxt, &
649 para_env=preconditioner_env%para_env)
650 CALL cp_fm_create(matrix_s1, fm_struct_tmp, name="matrix_s1")
651 CALL cp_fm_create(matrix_s2, fm_struct_tmp, name="matrix_s2")
652 CALL cp_fm_struct_release(fm_struct_tmp)
653 ALLOCATE (norms(k))
654 CALL parallel_gemm('T', 'N', k, k, n, 1.0_dp, matrix_sc0, matrix_tmp, 0.0_dp, matrix_s1)
655 CALL choose_eigv_solver(matrix_s1, matrix_s2, norms)
656
657 WRITE (*, *) "matrix norm deviation (should be close to zero): ", maxval(abs(abs(norms) - 1.0_dp))
658 DEALLOCATE (norms)
659 CALL cp_fm_release(matrix_s1)
660 CALL cp_fm_release(matrix_s2)
661 END IF
662
663 ! 3) replace the lowest k evals and evecs with what they should be
664 preconditioner_env%occ_evals = c0_evals
665 ! notice, this choice causes the preconditioner to be constant when applied to sc0 (see apply_full_all)
666 preconditioner_env%full_evals(1:k) = c0_evals
667 CALL cp_fm_to_fm(matrix_c0, matrix_pre, k, 1, 1)
668
669 CALL cp_fm_release(matrix_sc0)
670 CALL cp_fm_release(matrix_hc0)
671 CALL cp_fm_release(matrix_tmp)
672 DEALLOCATE (shifted_evals)
673
674 CALL timestop(handle)
675
676 END SUBROUTINE make_full_all_ortho
677
678! **************************************************************************************************
679!> \brief generates a preconditioner matrix H-lambda S+(SC)(2.0*CT*H*C+delta)(SC)^T
680!> for later inversion.
681!> H is the Kohn Sham matrix
682!> lambda*S shifts the spectrum of the generalized form up by lambda
683!> the last term only shifts the occupied space (reversing them in energy order)
684!> This form is implicitly multiplied from both sides by S^0.5
685!> This ensures we precondition the correct quantity
686!> Before this reads S^-0.5 H S^-0.5 + lambda + (S^0.5 C)shifts(S^0.5 C)T
687!> which might be a bit more obvious
688!> Replaced the old full_single_inverse at revision 14616
689!> \param preconditioner_env the preconditioner env
690!> \param matrix_c0 the MO coefficient matrix (fm)
691!> \param matrix_h Kohn-Sham matrix (dbcsr)
692!> \param energy_gap an additional shift in lambda=-E_homo+energy_gap
693!> \param matrix_s the overlap matrix if not orthonormal (dbcsr, optional)
694! **************************************************************************************************
695 SUBROUTINE make_full_single_inverse(preconditioner_env, matrix_c0, matrix_h, energy_gap, matrix_s)
696 TYPE(preconditioner_type) :: preconditioner_env
697 TYPE(cp_fm_type), INTENT(IN) :: matrix_c0
698 TYPE(dbcsr_type), POINTER :: matrix_h
699 REAL(kind=dp) :: energy_gap
700 TYPE(dbcsr_type), OPTIONAL, POINTER :: matrix_s
701
702 CHARACTER(len=*), PARAMETER :: routinen = 'make_full_single_inverse'
703
704 INTEGER :: handle, k, n
705 REAL(kind=dp) :: max_ev, min_ev, pre_shift
706 TYPE(arnoldi_data_type) :: my_arnoldi
707 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrices
708 TYPE(dbcsr_type), TARGET :: dbcsr_cthc, dbcsr_hc, dbcsr_sc, mo_dbcsr
709
710 CALL timeset(routinen, handle)
711
712 ! Allocate all working matrices needed
713 CALL cp_fm_get_info(matrix_c0, nrow_global=n, ncol_global=k)
714 ! copy the fm MO's to a sparse matrix, can be solved better if the sparse version is already present
715 ! but for the time beeing this will do
716 CALL cp_fm_to_dbcsr_row_template(mo_dbcsr, matrix_c0, matrix_h)
717 CALL dbcsr_create(dbcsr_sc, template=mo_dbcsr)
718 CALL dbcsr_create(dbcsr_hc, template=mo_dbcsr)
719 CALL cp_dbcsr_m_by_n_from_template(dbcsr_cthc, matrix_h, k, k, sym=dbcsr_type_symmetric)
720
721 ! Check whether the output matrix was already created, if not do it now
722 IF (.NOT. ASSOCIATED(preconditioner_env%sparse_matrix)) THEN
723 ALLOCATE (preconditioner_env%sparse_matrix)
724 END IF
725
726 ! Put the first term of the preconditioner (H) into the output matrix
727 CALL dbcsr_copy(preconditioner_env%sparse_matrix, matrix_h)
728
729 ! Precompute some matrices
730 ! S*C, if orthonormal this will be simply C so a copy will do
731 IF (PRESENT(matrix_s)) THEN
732 CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_s, mo_dbcsr, 0.0_dp, dbcsr_sc)
733 ELSE
734 CALL dbcsr_copy(dbcsr_sc, mo_dbcsr)
735 END IF
736
737!----------------------------compute the occupied subspace and shift it ------------------------------------
738 ! cT*H*C which will be used to shift the occupied states to 0
739 CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_h, mo_dbcsr, 0.0_dp, dbcsr_hc)
740 CALL dbcsr_multiply("T", "N", 1.0_dp, mo_dbcsr, dbcsr_hc, 0.0_dp, dbcsr_cthc)
741
742 ! Compute the Energy of the HOMO. We will use this as a reference energy
743 ALLOCATE (matrices(1))
744 matrices(1)%matrix => dbcsr_cthc
745 CALL setup_arnoldi_data(my_arnoldi, matrices, max_iter=20, threshold=1.0e-3_dp, selection_crit=2, &
746 nval_request=1, nrestarts=8, generalized_ev=.false., iram=.false.)
747 IF (ASSOCIATED(preconditioner_env%max_ev_vector)) &
748 CALL set_arnoldi_initial_vector(my_arnoldi, preconditioner_env%max_ev_vector)
749 CALL arnoldi_ev(matrices, my_arnoldi)
750 max_ev = real(get_selected_ritz_val(my_arnoldi, 1), dp)
751
752 ! save the ev as guess for the next time
753 IF (.NOT. ASSOCIATED(preconditioner_env%max_ev_vector)) ALLOCATE (preconditioner_env%max_ev_vector)
754 CALL get_selected_ritz_vector(my_arnoldi, 1, matrices(1)%matrix, preconditioner_env%max_ev_vector)
755 CALL deallocate_arnoldi_data(my_arnoldi)
756 DEALLOCATE (matrices)
757
758 ! Lets shift the occupied states a bit further up, -1.0 because we gonna subtract it from H
759 CALL dbcsr_add_on_diag(dbcsr_cthc, -0.5_dp)
760 ! Get the AO representation of the shift (see above why S is needed), W-matrix like object
761 CALL dbcsr_multiply("N", "N", 2.0_dp, dbcsr_sc, dbcsr_cthc, 0.0_dp, dbcsr_hc)
762 CALL dbcsr_multiply("N", "T", -1.0_dp, dbcsr_hc, dbcsr_sc, 1.0_dp, preconditioner_env%sparse_matrix)
763
764!-------------------------------------compute eigenvalues of H ----------------------------------------------
765 ! Setup the arnoldi procedure to compute the lowest ev. if S is present this has to be the generalized ev
766 IF (PRESENT(matrix_s)) THEN
767 ALLOCATE (matrices(2))
768 matrices(1)%matrix => preconditioner_env%sparse_matrix
769 matrices(2)%matrix => matrix_s
770 CALL setup_arnoldi_data(my_arnoldi, matrices, max_iter=20, threshold=2.0e-2_dp, selection_crit=3, &
771 nval_request=1, nrestarts=21, generalized_ev=.true., iram=.false.)
772 ELSE
773 ALLOCATE (matrices(1))
774 matrices(1)%matrix => preconditioner_env%sparse_matrix
775 CALL setup_arnoldi_data(my_arnoldi, matrices, max_iter=20, threshold=2.0e-2_dp, selection_crit=3, &
776 nval_request=1, nrestarts=8, generalized_ev=.false., iram=.false.)
777 END IF
778 IF (ASSOCIATED(preconditioner_env%min_ev_vector)) &
779 CALL set_arnoldi_initial_vector(my_arnoldi, preconditioner_env%min_ev_vector)
780
781 ! compute the LUMO energy
782 CALL arnoldi_ev(matrices, my_arnoldi)
783 min_ev = real(get_selected_ritz_val(my_arnoldi, 1), dp)
784
785 ! save the lumo vector for restarting in the next step
786 IF (.NOT. ASSOCIATED(preconditioner_env%min_ev_vector)) ALLOCATE (preconditioner_env%min_ev_vector)
787 CALL get_selected_ritz_vector(my_arnoldi, 1, matrices(1)%matrix, preconditioner_env%min_ev_vector)
788 CALL deallocate_arnoldi_data(my_arnoldi)
789 DEALLOCATE (matrices)
790
791!-------------------------------------compute eigenvalues of H ----------------------------------------------
792 ! Shift the Lumo to the 1.5*the computed energy_gap or the external energy gap value
793 ! The factor 1.5 is determined by trying. If the LUMO is positive, enough, just leave it alone
794 pre_shift = max(1.5_dp*(min_ev - max_ev), energy_gap)
795 IF (min_ev .LT. pre_shift) THEN
796 pre_shift = pre_shift - min_ev
797 ELSE
798 pre_shift = 0.0_dp
799 END IF
800 IF (PRESENT(matrix_s)) THEN
801 CALL dbcsr_add(preconditioner_env%sparse_matrix, matrix_s, 1.0_dp, pre_shift)
802 ELSE
803 CALL dbcsr_add_on_diag(preconditioner_env%sparse_matrix, pre_shift)
804 END IF
805
806 CALL dbcsr_release(mo_dbcsr)
807 CALL dbcsr_release(dbcsr_hc)
808 CALL dbcsr_release(dbcsr_sc)
809 CALL dbcsr_release(dbcsr_cthc)
810
811 CALL timestop(handle)
812
813 END SUBROUTINE make_full_single_inverse
814
815END MODULE preconditioner_makes
816
cp_fm_types::cp_fm_release
Definition cp_fm_types.F:90
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:85
parallel_gemm_api::parallel_gemm
Definition parallel_gemm_api.F:38
arnoldi_api
arnoldi iteration using dbcsr
Definition arnoldi_api.F:15
arnoldi_api::arnoldi_ev
subroutine, public arnoldi_ev(matrix, arnoldi_data)
Driver routine for different arnoldi eigenvalue methods the selection which one is to be taken is mad...
Definition arnoldi_api.F:69
cp_dbcsr_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_dbcsr_operations::cp_dbcsr_sm_fm_multiply
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
Definition cp_dbcsr_operations.F:744
cp_dbcsr_operations::cp_dbcsr_m_by_n_from_template
subroutine, public cp_dbcsr_m_by_n_from_template(matrix, template, m, n, sym, data_type)
Utility function to create an arbitrary shaped dbcsr matrix with the same processor grid as the templ...
Definition cp_dbcsr_operations.F:1105
cp_dbcsr_operations::copy_dbcsr_to_fm
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:208
cp_dbcsr_operations::cp_fm_to_dbcsr_row_template
subroutine, public cp_fm_to_dbcsr_row_template(matrix, fm_in, template)
Utility function to copy a specially shaped fm to dbcsr_matrix The result matrix will be the matrix i...
Definition cp_dbcsr_operations.F:1070
cp_fm_basic_linalg
basic linear algebra operations for full matrices
Definition cp_fm_basic_linalg.F:14
cp_fm_basic_linalg::cp_fm_column_scale
subroutine, public cp_fm_column_scale(matrixa, scaling)
scales column i of matrix a with scaling(i)
Definition cp_fm_basic_linalg.F:1882
cp_fm_basic_linalg::cp_fm_upper_to_full
subroutine, public cp_fm_upper_to_full(matrix, work)
given an upper triangular matrix computes the corresponding full matrix
Definition cp_fm_basic_linalg.F:1759
cp_fm_basic_linalg::cp_fm_cholesky_restore
subroutine, public cp_fm_cholesky_restore(fm_matrix, neig, fm_matrixb, fm_matrixout, op, pos, transa)
...
Definition cp_fm_basic_linalg.F:3162
cp_fm_basic_linalg::cp_fm_triangular_invert
subroutine, public cp_fm_triangular_invert(matrix_a, uplo_tr)
inverts a triangular matrix
Definition cp_fm_basic_linalg.F:2193
cp_fm_basic_linalg::cp_fm_triangular_multiply
subroutine, public cp_fm_triangular_multiply(triangular_matrix, matrix_b, side, transpose_tr, invert_tr, uplo_tr, unit_diag_tr, n_rows, n_cols, alpha)
multiplies in place by a triangular matrix: matrix_b = alpha op(triangular_matrix) matrix_b or (if si...
Definition cp_fm_basic_linalg.F:1592
cp_fm_cholesky
various cholesky decomposition related routines
Definition cp_fm_cholesky.F:14
cp_fm_cholesky::cp_fm_cholesky_decompose
subroutine, public cp_fm_cholesky_decompose(matrix, n, info_out)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U,...
Definition cp_fm_cholesky.F:50
cp_fm_cholesky::cp_fm_cholesky_reduce
subroutine, public cp_fm_cholesky_reduce(matrix, matrixb, itype)
reduce a matrix pencil A,B to normal form B has to be cholesky decomposed with cp_fm_cholesky_decompo...
Definition cp_fm_cholesky.F:192
cp_fm_diag
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
cp_fm_diag::choose_eigv_solver
subroutine, public choose_eigv_solver(matrix, eigenvectors, eigenvalues, info)
Choose the Eigensolver depending on which library is available ELPA seems to be unstable for small sy...
Definition cp_fm_diag.F:208
cp_fm_struct
represent the structure of a full matrix
Definition cp_fm_struct.F:14
cp_fm_struct::cp_fm_struct_create
subroutine, public cp_fm_struct_create(fmstruct, para_env, context, nrow_global, ncol_global, nrow_block, ncol_block, descriptor, first_p_pos, local_leading_dimension, template_fmstruct, square_blocks, force_block)
allocates and initializes a full matrix structure
Definition cp_fm_struct.F:125
cp_fm_struct::cp_fm_struct_release
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
Definition cp_fm_struct.F:320
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_diag
subroutine, public cp_fm_get_diag(matrix, diag)
returns the diagonal elements of a fm
Definition cp_fm_types.F:570
cp_fm_types::cp_fm_get_info
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
Definition cp_fm_types.F:1016
cp_fm_types::cp_fm_create
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
creates a new full matrix with the given structure
Definition cp_fm_types.F:167
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::ot_precond_full_kinetic
integer, parameter, public ot_precond_full_kinetic
Definition input_constants.F:513
input_constants::cholesky_reduce
integer, parameter, public cholesky_reduce
Definition input_constants.F:195
input_constants::cholesky_inverse
integer, parameter, public cholesky_inverse
Definition input_constants.F:195
input_constants::ot_precond_solver_default
integer, parameter, public ot_precond_solver_default
Definition input_constants.F:520
input_constants::ot_precond_full_single
integer, parameter, public ot_precond_full_single
Definition input_constants.F:513
input_constants::ot_precond_solver_inv_chol
integer, parameter, public ot_precond_solver_inv_chol
Definition input_constants.F:520
input_constants::ot_precond_full_single_inverse
integer, parameter, public ot_precond_full_single_inverse
Definition input_constants.F:513
input_constants::ot_precond_s_inverse
integer, parameter, public ot_precond_s_inverse
Definition input_constants.F:513
input_constants::ot_precond_full_all
integer, parameter, public ot_precond_full_all
Definition input_constants.F:513
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
parallel_gemm_api
basic linear algebra operations for full matrixes
Definition parallel_gemm_api.F:14
preconditioner_makes
computes preconditioners, and implements methods to apply them currently used in qs_ot
Definition preconditioner_makes.F:15
preconditioner_makes::make_preconditioner_matrix
subroutine, public make_preconditioner_matrix(preconditioner_env, matrix_h, matrix_s, matrix_t, mo_coeff, energy_homo, eigenvalues_ot, energy_gap, my_solver_type)
...
Definition preconditioner_makes.F:81
preconditioner_types
types of preconditioners
Definition preconditioner_types.F:14
cp_fm_struct::cp_fm_struct_type
keeps the information about the structure of a full matrix
Definition cp_fm_struct.F:82
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:116
preconditioner_types::preconditioner_type
Definition preconditioner_types.F:42