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