(git:08825b9)
Loading...
Searching...
No Matches
preconditioner.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief 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! **************************************************************************************************
17 USE cp_cfm_types, ONLY: cp_cfm_create,&
22 USE cp_dbcsr_api, ONLY: dbcsr_get_info,&
35 USE cp_fm_types, ONLY: cp_fm_create,&
44 USE input_constants, ONLY: &
54 USE kinds, ONLY: default_string_length,&
55 dp
93#include "./base/base_uses.f90"
94
95 IMPLICIT NONE
96
97 PRIVATE
98
99 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'preconditioner'
100 REAL(KIND=dp), PARAMETER, PRIVATE :: fermi_low_rank_gap_factor = 4.0_dp, &
101 fermi_low_rank_degeneracy_tolerance = 1.0e-8_dp, &
102 fermi_low_rank_min_window = 1.0_dp
103
113
114! The public interface for apply preconditioner, the routines can be found in preconditioner_apply.F
116 MODULE PROCEDURE apply_preconditioner_cfm_complex
117 MODULE PROCEDURE apply_preconditioner_dbcsr
119 MODULE PROCEDURE apply_preconditioner_fm
120 END INTERFACE
121
122! **************************************************************************************************
123
124CONTAINS
125
126! **************************************************************************************************
127
128! **************************************************************************************************
129!> \brief Build the block-circulant inverse of a replicated-cell SPD operator.
130!> \param preconditioner_env ...
131!> \param operator_matrix ...
132!> \param operator_name ...
133!> \param rhs_count number of occupied orbital columns to precondition
134! **************************************************************************************************
135 SUBROUTINE prepare_lattice_operator(preconditioner_env, operator_matrix, operator_name, rhs_count)
136
137 TYPE(preconditioner_type) :: preconditioner_env
138 TYPE(dbcsr_type), POINTER :: operator_matrix
139 CHARACTER(LEN=*), INTENT(IN) :: operator_name
140 INTEGER, INTENT(IN) :: rhs_count
141
142 INTEGER :: alpha, beta, block_size, col_cell, delta_index, i, info, j, local_info, n, ncell, &
143 ncol_local, nrow_local, output_unit, row_cell
144 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: index_to_cell
145 INTEGER, DIMENSION(3) :: delta
146 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
147 LOGICAL :: used_fft
148 REAL(KIND=dp) :: denominator, difference, numerator, value
149 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: dense_operator, xkp
150 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: blocks
151 REAL(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
152 POINTER :: local_data
153 TYPE(cp_fm_struct_type), POINTER :: fm_struct
154 TYPE(cp_fm_type) :: operator_fm, work
155 TYPE(cp_logger_type), POINTER :: logger
156
157 cpassert(ASSOCIATED(operator_matrix))
158 cpassert(preconditioner_env%lattice_fft_mode /= ot_lattice_fft_off)
159 preconditioner_env%lattice_fft_active = .false.
160 preconditioner_env%lattice_projection_error = 0.0_dp
161 preconditioner_env%lattice_corrected_projection_error = 0.0_dp
162 preconditioner_env%lattice_cost_ratio = 0.0_dp
163 preconditioner_env%lattice_storage_ratio = 0.0_dp
164
165 IF (any(preconditioner_env%lattice_dims <= 0)) cpabort("Lattice dimensions must be positive")
166
167 IF (preconditioner_env%lattice_fft_mode == ot_lattice_fft_auto .AND. &
168 product(preconditioner_env%lattice_dims) < lattice_fft_auto_min_cells) THEN
169 IF (preconditioner_env%para_env%is_source()) THEN
170 logger => cp_get_default_logger()
171 output_unit = cp_logger_get_default_io_unit(logger)
172 IF (output_unit > 0) THEN
173 WRITE (output_unit, '(T2,A,A,A,I0,A,I0)') &
174 "OT| Lattice FFT AUTO kept the conventional ", trim(operator_name), &
175 " solver: ", product(preconditioner_env%lattice_dims), &
176 " cells are fewer than ", lattice_fft_auto_min_cells
177 END IF
178 END IF
179 RETURN
180 END IF
181
182 NULLIFY (fm_struct)
183 CALL dbcsr_get_info(operator_matrix, nfullrows_total=n)
184 ncell = product(preconditioner_env%lattice_dims)
185 IF (modulo(n, ncell) /= 0) cpabort("AO dimension does not match replicated cells")
186 block_size = n/ncell
187 preconditioner_env%lattice_block_size = block_size
188 preconditioner_env%lattice_cost_ratio = lattice_fft_cost_ratio( &
189 preconditioner_env%lattice_dims, block_size, rhs_count, &
190 preconditioner_env%para_env%num_pe) + &
192 preconditioner_env%lattice_dims, &
193 preconditioner_env%lattice_local_cells)
194 preconditioner_env%lattice_storage_ratio = lattice_fft_storage_ratio( &
195 preconditioner_env%lattice_dims, 1) + &
197 preconditioner_env%lattice_dims, 1, &
198 preconditioner_env%lattice_local_cells)
199 IF (preconditioner_env%lattice_fft_mode == ot_lattice_fft_auto .AND. &
200 (preconditioner_env%lattice_cost_ratio > lattice_fft_auto_max_cost_ratio .OR. &
201 preconditioner_env%lattice_storage_ratio > lattice_fft_auto_max_storage_ratio)) THEN
202 IF (preconditioner_env%para_env%is_source()) THEN
203 logger => cp_get_default_logger()
204 output_unit = cp_logger_get_default_io_unit(logger)
205 IF (output_unit > 0) THEN
206 WRITE (output_unit, '(T2,A,A,A,F8.3,A,F8.3,A)') &
207 "OT| Lattice FFT AUTO kept the conventional ", trim(operator_name), &
208 " solver: estimated work ratio ", preconditioner_env%lattice_cost_ratio, &
209 " storage ratio ", preconditioner_env%lattice_storage_ratio, " selected: F"
210 END IF
211 END IF
212 RETURN
213 END IF
214 CALL cp_fm_struct_create(fm_struct, nrow_global=n, ncol_global=n, &
215 context=preconditioner_env%ctxt, &
216 para_env=preconditioner_env%para_env)
217 CALL cp_fm_create(operator_fm, fm_struct, name="lattice preconditioner operator")
218 CALL cp_fm_struct_release(fm_struct)
219 CALL cp_fm_set_all(operator_fm, alpha=0.0_dp, beta=0.0_dp)
220 CALL copy_dbcsr_to_fm(operator_matrix, operator_fm)
221 CALL cp_fm_create(work, operator_fm%matrix_struct, name="lattice operator work")
222 CALL cp_fm_uplo_to_full(operator_fm, work)
223 CALL cp_fm_release(work)
224 CALL cp_fm_get_info(operator_fm, nrow_local=nrow_local, ncol_local=ncol_local, &
225 row_indices=row_indices, col_indices=col_indices, local_data=local_data)
226 ALLOCATE (blocks(block_size, block_size, ncell), source=0.0_dp)
227 ALLOCATE (index_to_cell(3, ncell), xkp(3, ncell))
228 CALL lattice_grid(preconditioner_env%lattice_dims, index_to_cell, xkp)
229 IF (preconditioner_env%lattice_local_cells > 0) THEN
230 ALLOCATE (dense_operator(n, n), source=0.0_dp)
231 END IF
232 denominator = 0.0_dp
233 DO j = 1, ncol_local
234 DO i = 1, nrow_local
235 IF (row_indices(i) <= col_indices(j)) THEN
236 value = local_data(i, j)
237 row_cell = (row_indices(i) - 1)/block_size + 1
238 col_cell = (col_indices(j) - 1)/block_size + 1
239 alpha = modulo(row_indices(i) - 1, block_size) + 1
240 beta = modulo(col_indices(j) - 1, block_size) + 1
241 delta = modulo(index_to_cell(:, row_cell) - index_to_cell(:, col_cell), &
242 preconditioner_env%lattice_dims)
243 delta_index = delta(1) + preconditioner_env%lattice_dims(1)*(delta(2) + &
244 preconditioner_env%lattice_dims(2)*delta(3)) + 1
245 blocks(alpha, beta, delta_index) = blocks(alpha, beta, delta_index) + &
246 value/real(ncell, kind=dp)
247 denominator = denominator + value*value
248 IF (ALLOCATED(dense_operator)) THEN
249 dense_operator(row_indices(i), col_indices(j)) = value
250 END IF
251 IF (row_indices(i) < col_indices(j)) THEN
252 delta = modulo(index_to_cell(:, col_cell) - index_to_cell(:, row_cell), &
253 preconditioner_env%lattice_dims)
254 delta_index = delta(1) + preconditioner_env%lattice_dims(1)*(delta(2) + &
255 preconditioner_env%lattice_dims(2)*delta(3)) + 1
256 blocks(beta, alpha, delta_index) = blocks(beta, alpha, delta_index) + &
257 value/real(ncell, kind=dp)
258 denominator = denominator + value*value
259 END IF
260 END IF
261 END DO
262 END DO
263 CALL preconditioner_env%para_env%sum(blocks)
264 numerator = 0.0_dp
265 DO j = 1, ncol_local
266 DO i = 1, nrow_local
267 IF (row_indices(i) <= col_indices(j)) THEN
268 value = local_data(i, j)
269 row_cell = (row_indices(i) - 1)/block_size + 1
270 col_cell = (col_indices(j) - 1)/block_size + 1
271 alpha = modulo(row_indices(i) - 1, block_size) + 1
272 beta = modulo(col_indices(j) - 1, block_size) + 1
273 delta = modulo(index_to_cell(:, row_cell) - index_to_cell(:, col_cell), &
274 preconditioner_env%lattice_dims)
275 delta_index = delta(1) + preconditioner_env%lattice_dims(1)*(delta(2) + &
276 preconditioner_env%lattice_dims(2)*delta(3)) + 1
277 difference = value - blocks(alpha, beta, delta_index)
278 numerator = numerator + difference*difference
279 IF (row_indices(i) < col_indices(j)) THEN
280 delta = modulo(index_to_cell(:, col_cell) - index_to_cell(:, row_cell), &
281 preconditioner_env%lattice_dims)
282 delta_index = delta(1) + preconditioner_env%lattice_dims(1)*(delta(2) + &
283 preconditioner_env%lattice_dims(2)*delta(3)) + 1
284 difference = value - blocks(beta, alpha, delta_index)
285 numerator = numerator + difference*difference
286 END IF
287 END IF
288 END DO
289 END DO
290 CALL preconditioner_env%para_env%sum(numerator)
291 CALL preconditioner_env%para_env%sum(denominator)
292 preconditioner_env%lattice_projection_error = &
293 sqrt(numerator/max(denominator, tiny(denominator)))
295 blocks, preconditioner_env%lattice_dims, preconditioner_env%lattice_inverse_k, used_fft, info)
296 preconditioner_env%lattice_corrected_projection_error = &
297 preconditioner_env%lattice_projection_error
298 IF (ALLOCATED(dense_operator)) THEN
299 CALL preconditioner_env%para_env%sum(dense_operator)
300 DO j = 1, n
301 DO i = j + 1, n
302 dense_operator(i, j) = dense_operator(j, i)
303 END DO
304 END DO
306 dense_operator, preconditioner_env%lattice_dims, &
307 preconditioner_env%lattice_local_cells, &
308 preconditioner_env%lattice_local_indices, &
309 preconditioner_env%lattice_local_operator, &
310 preconditioner_env%lattice_local_coarse_inverse, &
311 preconditioner_env%lattice_corrected_projection_error, local_info)
312 IF (local_info /= 0) THEN
313 cpabort("LATTICE_FFT local correction found a non-positive coarse operator")
314 END IF
315 DEALLOCATE (dense_operator)
316 END IF
317 CALL cp_fm_release(operator_fm)
318 DEALLOCATE (blocks, index_to_cell, xkp)
319 SELECT CASE (preconditioner_env%lattice_fft_mode)
320 CASE (ot_lattice_fft_on)
321 IF (info /= 0) THEN
322 cpabort("LATTICE_FFT found a non-positive Fourier block for "//trim(operator_name))
323 END IF
324 IF (.NOT. used_fft) THEN
325 cpabort("LATTICE_FFT could not initialize its FFT for "//trim(operator_name))
326 END IF
327 preconditioner_env%lattice_fft_active = .true.
329 preconditioner_env%lattice_fft_active = lattice_fft_auto_select( &
330 preconditioner_env%lattice_dims, &
331 preconditioner_env%lattice_block_size, rhs_count, 1, &
332 preconditioner_env%para_env%num_pe, &
333 preconditioner_env%lattice_projection_error, used_fft, info)
334 CASE DEFAULT
335 cpabort("Unknown LATTICE_FFT mode")
336 END SELECT
337
338 IF (.NOT. preconditioner_env%lattice_fft_active .AND. &
339 ALLOCATED(preconditioner_env%lattice_inverse_k)) THEN
340 DEALLOCATE (preconditioner_env%lattice_inverse_k)
341 END IF
342 IF (preconditioner_env%para_env%is_source()) THEN
343 logger => cp_get_default_logger()
344 output_unit = cp_logger_get_default_io_unit(logger)
345 IF (output_unit > 0) THEN
346 WRITE (output_unit, '(T2,A,A,A,3(I0,1X),A,I0)') &
347 "OT| Lattice FFT ", trim(operator_name), " cells: ", &
348 preconditioner_env%lattice_dims, "AO block: ", preconditioner_env%lattice_block_size
349 WRITE (output_unit, '(T2,A,ES12.4,A,L1)') &
350 "OT| Block-circulant projection error: ", &
351 preconditioner_env%lattice_projection_error, &
352 " selected: ", preconditioner_env%lattice_fft_active
353 IF (preconditioner_env%lattice_local_cells > 0) THEN
354 WRITE (output_unit, '(T2,A,I0,A,ES12.4)') &
355 "OT| Balanced local correction cells: ", &
356 preconditioner_env%lattice_local_cells, " residual error: ", &
357 preconditioner_env%lattice_corrected_projection_error
358 END IF
359 WRITE (output_unit, '(T2,A,F8.3,A,F8.3)') &
360 "OT| Estimated lattice/dense work ratio: ", &
361 preconditioner_env%lattice_cost_ratio, " storage ratio: ", &
362 preconditioner_env%lattice_storage_ratio
363 END IF
364 END IF
365
366 END SUBROUTINE prepare_lattice_operator
367
368! **************************************************************************************************
369!> \brief Build state-dependent block-circulant inverses for rotation-covariant FULL_ALL.
370!> \param preconditioner_env ...
371!> \param matrix_s overlap matrix
372!> \param operator_name ...
373! **************************************************************************************************
374 SUBROUTINE prepare_lattice_state_operator(preconditioner_env, matrix_s, operator_name)
375
376 TYPE(preconditioner_type) :: preconditioner_env
377 TYPE(dbcsr_type), POINTER :: matrix_s
378 CHARACTER(LEN=*), INTENT(IN) :: operator_name
379
380 COMPLEX(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :) :: inverse_k
381 INTEGER :: i, info, istate, j, local_info, &
382 local_rank, n, ncell, ncol_local, &
383 nrow_local, nstate, output_unit
384 INTEGER, ALLOCATABLE, DIMENSION(:) :: local_indices
385 INTEGER, DIMENSION(:), POINTER :: col_indices, row_indices
386 LOGICAL :: all_selected, selected, used_fft
387 REAL(KIND=dp) :: corrected_projection_error, &
388 projection_error
389 REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: dense_operator, dense_sv, dense_vectors, &
390 local_coarse_inverse, local_operator, &
391 scaled_sv
392 REAL(KIND=dp), CONTIGUOUS, DIMENSION(:, :), &
393 POINTER :: local_data
394 TYPE(cp_fm_type) :: sv_fm
395 TYPE(cp_logger_type), POINTER :: logger
396
397 cpassert(ASSOCIATED(matrix_s))
398 cpassert(ASSOCIATED(preconditioner_env%fm))
399 cpassert(ASSOCIATED(preconditioner_env%full_evals))
400 cpassert(ASSOCIATED(preconditioner_env%occ_evals))
401 cpassert(preconditioner_env%lattice_fft_mode /= ot_lattice_fft_off)
402
403 preconditioner_env%lattice_fft_active = .false.
404 preconditioner_env%lattice_projection_error = 0.0_dp
405 preconditioner_env%lattice_corrected_projection_error = 0.0_dp
406 preconditioner_env%lattice_cost_ratio = 0.0_dp
407 preconditioner_env%lattice_storage_ratio = 0.0_dp
408 IF (ALLOCATED(preconditioner_env%lattice_state_inverse_k)) THEN
409 DEALLOCATE (preconditioner_env%lattice_state_inverse_k)
410 END IF
411 IF (ALLOCATED(preconditioner_env%lattice_occ_vectors)) THEN
412 DEALLOCATE (preconditioner_env%lattice_occ_vectors)
413 END IF
414 IF (ALLOCATED(preconditioner_env%lattice_occ_dual)) THEN
415 DEALLOCATE (preconditioner_env%lattice_occ_dual)
416 END IF
417 IF (ALLOCATED(preconditioner_env%lattice_local_indices)) THEN
418 DEALLOCATE (preconditioner_env%lattice_local_indices)
419 END IF
420 IF (ALLOCATED(preconditioner_env%lattice_local_operator)) THEN
421 DEALLOCATE (preconditioner_env%lattice_local_operator)
422 END IF
423 IF (ALLOCATED(preconditioner_env%lattice_local_coarse_inverse)) THEN
424 DEALLOCATE (preconditioner_env%lattice_local_coarse_inverse)
425 END IF
426 IF (ALLOCATED(preconditioner_env%lattice_state_local_indices)) THEN
427 DEALLOCATE (preconditioner_env%lattice_state_local_indices)
428 END IF
429 IF (ALLOCATED(preconditioner_env%lattice_state_local_operator)) THEN
430 DEALLOCATE (preconditioner_env%lattice_state_local_operator)
431 END IF
432 IF (ALLOCATED(preconditioner_env%lattice_state_local_coarse_inverse)) THEN
433 DEALLOCATE (preconditioner_env%lattice_state_local_coarse_inverse)
434 END IF
435 IF (any(preconditioner_env%lattice_dims <= 0)) cpabort("Lattice dimensions must be positive")
436 IF (preconditioner_env%lattice_fft_mode == ot_lattice_fft_auto .AND. &
437 product(preconditioner_env%lattice_dims) < lattice_fft_auto_min_cells) THEN
438 IF (preconditioner_env%para_env%is_source()) THEN
439 logger => cp_get_default_logger()
440 output_unit = cp_logger_get_default_io_unit(logger)
441 IF (output_unit > 0) THEN
442 WRITE (output_unit, '(T2,A,A,A,I0,A,I0)') &
443 "OT| Lattice FFT AUTO kept the conventional ", trim(operator_name), &
444 " solver: ", product(preconditioner_env%lattice_dims), &
445 " cells are fewer than ", lattice_fft_auto_min_cells
446 END IF
447 END IF
448 RETURN
449 END IF
450
451 CALL cp_fm_get_info(preconditioner_env%fm, nrow_global=n, ncol_global=j)
452 cpassert(j == n .AND. SIZE(preconditioner_env%full_evals) == n)
453 nstate = SIZE(preconditioner_env%occ_evals)
454 cpassert(nstate > 0)
455 ncell = product(preconditioner_env%lattice_dims)
456 IF (modulo(n, ncell) /= 0) cpabort("AO dimension does not match replicated cells")
457 preconditioner_env%lattice_block_size = n/ncell
458 preconditioner_env%lattice_cost_ratio = lattice_fft_cost_ratio( &
459 preconditioner_env%lattice_dims, &
460 preconditioner_env%lattice_block_size, nstate, &
461 preconditioner_env%para_env%num_pe)
462 preconditioner_env%lattice_storage_ratio = lattice_fft_storage_ratio( &
463 preconditioner_env%lattice_dims, nstate)
464 preconditioner_env%lattice_cost_ratio = preconditioner_env%lattice_cost_ratio + &
466 preconditioner_env%lattice_dims, &
467 preconditioner_env%lattice_local_cells)
468 preconditioner_env%lattice_storage_ratio = preconditioner_env%lattice_storage_ratio + &
470 preconditioner_env%lattice_dims, nstate, &
471 preconditioner_env%lattice_local_cells)
472 IF (preconditioner_env%lattice_fft_mode == ot_lattice_fft_auto .AND. &
473 (preconditioner_env%lattice_cost_ratio > lattice_fft_auto_max_cost_ratio .OR. &
474 preconditioner_env%lattice_storage_ratio > lattice_fft_auto_max_storage_ratio)) THEN
475 IF (preconditioner_env%para_env%is_source()) THEN
476 logger => cp_get_default_logger()
477 output_unit = cp_logger_get_default_io_unit(logger)
478 IF (output_unit > 0) THEN
479 WRITE (output_unit, '(T2,A,A,A,F8.3,A,F8.3,A)') &
480 "OT| Lattice FFT AUTO kept the conventional ", trim(operator_name), &
481 " solver: estimated work ratio ", preconditioner_env%lattice_cost_ratio, &
482 " storage ratio ", preconditioner_env%lattice_storage_ratio, " selected: F"
483 END IF
484 END IF
485 RETURN
486 END IF
487 IF (preconditioner_env%lattice_local_cells > 0) THEN
488 local_rank = preconditioner_env%lattice_local_cells* &
489 preconditioner_env%lattice_block_size
490 ALLOCATE (preconditioner_env%lattice_state_local_indices(local_rank, nstate), &
491 preconditioner_env%lattice_state_local_operator(n, local_rank, nstate), &
492 preconditioner_env%lattice_state_local_coarse_inverse(local_rank, local_rank, nstate))
493 END IF
494 CALL cp_fm_create(sv_fm, preconditioner_env%fm%matrix_struct, name="lattice FULL_ALL S V")
495 CALL cp_dbcsr_sm_fm_multiply(matrix_s, preconditioner_env%fm, sv_fm, n)
496 CALL cp_fm_get_info(sv_fm, nrow_local=nrow_local, ncol_local=ncol_local, &
497 row_indices=row_indices, col_indices=col_indices, local_data=local_data)
498 ALLOCATE (dense_sv(n, n), dense_vectors(n, nstate), source=0.0_dp)
499 DO j = 1, ncol_local
500 DO i = 1, nrow_local
501 dense_sv(row_indices(i), col_indices(j)) = local_data(i, j)
502 END DO
503 END DO
504 CALL preconditioner_env%para_env%sum(dense_sv)
505 CALL cp_fm_release(sv_fm)
506
507 CALL cp_fm_get_info(preconditioner_env%fm, nrow_local=nrow_local, ncol_local=ncol_local, &
508 row_indices=row_indices, col_indices=col_indices, local_data=local_data)
509 DO j = 1, ncol_local
510 IF (col_indices(j) > nstate) cycle
511 DO i = 1, nrow_local
512 dense_vectors(row_indices(i), col_indices(j)) = local_data(i, j)
513 END DO
514 END DO
515 CALL preconditioner_env%para_env%sum(dense_vectors)
516 ALLOCATE (preconditioner_env%lattice_occ_vectors(n, nstate), &
517 preconditioner_env%lattice_occ_dual(n, nstate))
518 preconditioner_env%lattice_occ_vectors(:, :) = dense_vectors
519 preconditioner_env%lattice_occ_dual(:, :) = dense_sv(:, 1:nstate)
520
521 ALLOCATE (dense_operator(n, n), scaled_sv(n, n))
522 all_selected = .true.
523 DO istate = 1, nstate
524 scaled_sv(:, :) = dense_sv
525 DO j = 1, n
526 scaled_sv(:, j) = scaled_sv(:, j)*max( &
527 preconditioner_env%energy_gap, &
528 preconditioner_env%full_evals(j) - preconditioner_env%occ_evals(istate))
529 END DO
530 dense_operator(:, :) = matmul(scaled_sv, transpose(dense_sv))
531 CALL build_lattice_inverse(dense_operator, preconditioner_env%lattice_dims, inverse_k, &
532 projection_error, used_fft, info)
533 preconditioner_env%lattice_projection_error = &
534 max(preconditioner_env%lattice_projection_error, projection_error)
535 corrected_projection_error = projection_error
536 IF (preconditioner_env%lattice_local_cells > 0) THEN
538 dense_operator, preconditioner_env%lattice_dims, &
539 preconditioner_env%lattice_local_cells, local_indices, &
540 local_operator, local_coarse_inverse, corrected_projection_error, local_info)
541 IF (local_info /= 0) THEN
542 cpabort("LATTICE_FFT local correction found a non-positive coarse operator")
543 END IF
544 preconditioner_env%lattice_state_local_indices(:, istate) = local_indices
545 preconditioner_env%lattice_state_local_operator(:, :, istate) = local_operator
546 preconditioner_env%lattice_state_local_coarse_inverse(:, :, istate) = local_coarse_inverse
547 DEALLOCATE (local_coarse_inverse, local_indices, local_operator)
548 END IF
549 preconditioner_env%lattice_corrected_projection_error = &
550 max(preconditioner_env%lattice_corrected_projection_error, corrected_projection_error)
551 selected = lattice_fft_auto_select(preconditioner_env%lattice_dims, &
552 preconditioner_env%lattice_block_size, nstate, nstate, &
553 preconditioner_env%para_env%num_pe, projection_error, used_fft, info)
554 IF (preconditioner_env%lattice_fft_mode == ot_lattice_fft_on) THEN
555 IF (info /= 0) THEN
556 cpabort("LATTICE_FFT found a non-positive Fourier block for "//trim(operator_name))
557 END IF
558 IF (.NOT. used_fft) THEN
559 cpabort("LATTICE_FFT could not initialize its FFT for "//trim(operator_name))
560 END IF
561 ELSE
562 all_selected = all_selected .AND. selected
563 END IF
564 IF (.NOT. ALLOCATED(preconditioner_env%lattice_state_inverse_k)) THEN
565 ALLOCATE (preconditioner_env%lattice_state_inverse_k( &
566 SIZE(inverse_k, 1), SIZE(inverse_k, 2), SIZE(inverse_k, 3), nstate))
567 END IF
568 IF (info == 0 .AND. used_fft) THEN
569 preconditioner_env%lattice_state_inverse_k(:, :, :, istate) = inverse_k
570 END IF
571 DEALLOCATE (inverse_k)
572 END DO
573
574 preconditioner_env%lattice_fft_active = &
575 preconditioner_env%lattice_fft_mode == ot_lattice_fft_on .OR. all_selected
576 IF (.NOT. preconditioner_env%lattice_fft_active .AND. &
577 ALLOCATED(preconditioner_env%lattice_state_inverse_k)) THEN
578 DEALLOCATE (preconditioner_env%lattice_state_inverse_k)
579 END IF
580 IF (.NOT. preconditioner_env%lattice_fft_active) THEN
581 DEALLOCATE (preconditioner_env%lattice_occ_vectors, &
582 preconditioner_env%lattice_occ_dual)
583 IF (ALLOCATED(preconditioner_env%lattice_state_local_indices)) THEN
584 DEALLOCATE (preconditioner_env%lattice_state_local_indices, &
585 preconditioner_env%lattice_state_local_operator, &
586 preconditioner_env%lattice_state_local_coarse_inverse)
587 END IF
588 END IF
589 IF (preconditioner_env%para_env%is_source()) THEN
590 logger => cp_get_default_logger()
591 output_unit = cp_logger_get_default_io_unit(logger)
592 IF (output_unit > 0) THEN
593 WRITE (output_unit, '(T2,A,A,A,3(I0,1X),A,I0,A,I0)') &
594 "OT| Lattice FFT ", trim(operator_name), " cells: ", &
595 preconditioner_env%lattice_dims, "AO block: ", &
596 preconditioner_env%lattice_block_size, " states: ", nstate
597 WRITE (output_unit, '(T2,A,ES12.4,A,L1)') &
598 "OT| Maximum block-circulant projection error: ", &
599 preconditioner_env%lattice_projection_error, &
600 " selected: ", preconditioner_env%lattice_fft_active
601 IF (preconditioner_env%lattice_local_cells > 0) THEN
602 WRITE (output_unit, '(T2,A,I0,A,ES12.4)') &
603 "OT| Balanced local correction cells/state: ", &
604 preconditioner_env%lattice_local_cells, " maximum residual error: ", &
605 preconditioner_env%lattice_corrected_projection_error
606 END IF
607 WRITE (output_unit, '(T2,A,F8.3,A,F8.3)') &
608 "OT| Estimated lattice/dense work ratio: ", &
609 preconditioner_env%lattice_cost_ratio, " storage ratio: ", &
610 preconditioner_env%lattice_storage_ratio
611 END IF
612 END IF
613
614 DEALLOCATE (dense_operator, dense_sv, dense_vectors, scaled_sv)
615
616 END SUBROUTINE prepare_lattice_state_operator
617
618! **************************************************************************************************
619!> \brief Build the overlap-inverse base of the bounded low-rank spectral preconditioner.
620!> \param preconditioner_env ...
621!> \param matrix_s ...
622! **************************************************************************************************
623 SUBROUTINE prepare_low_rank_base(preconditioner_env, matrix_s)
624
625 TYPE(preconditioner_type) :: preconditioner_env
626 TYPE(dbcsr_type), OPTIONAL, POINTER :: matrix_s
627
628 INTEGER :: info
629 TYPE(cp_fm_type) :: work
630
631 cpassert(ASSOCIATED(preconditioner_env%fm))
632 cpassert(.NOT. ASSOCIATED(preconditioner_env%base_fm))
633
634 IF (preconditioner_env%lattice_fft_mode /= ot_lattice_fft_off) THEN
635 IF (.NOT. PRESENT(matrix_s)) THEN
636 cpabort("LATTICE_FFT requires a non-unit overlap metric for FERMI_LOW_RANK")
637 END IF
638 cpassert(ASSOCIATED(matrix_s))
639 CALL prepare_lattice_operator(preconditioner_env, matrix_s, "FERMI_LOW_RANK overlap", &
640 SIZE(preconditioner_env%occ_evals))
641 IF (preconditioner_env%lattice_fft_active) THEN
642 preconditioner_env%low_rank_base = ot_low_rank_base_lattice_fft
643 RETURN
644 END IF
645 preconditioner_env%low_rank_base = ot_low_rank_base_overlap
646 END IF
647
648 IF (preconditioner_env%low_rank_base /= ot_low_rank_base_overlap) THEN
649 cpabort("Unknown FERMI_LOW_RANK base operator")
650 END IF
651 ALLOCATE (preconditioner_env%base_fm)
652 CALL cp_fm_create(preconditioner_env%base_fm, preconditioner_env%fm%matrix_struct, &
653 name="FERMI_LOW_RANK overlap inverse")
654
655 IF (PRESENT(matrix_s)) THEN
656 cpassert(ASSOCIATED(matrix_s))
657 CALL copy_dbcsr_to_fm(matrix_s, preconditioner_env%base_fm)
658 CALL cp_fm_cholesky_decompose(preconditioner_env%base_fm, info_out=info)
659 IF (info /= 0) cpabort("FERMI_LOW_RANK requires a positive-definite overlap matrix")
660 CALL cp_fm_cholesky_invert(preconditioner_env%base_fm)
661 CALL cp_fm_create(work, preconditioner_env%base_fm%matrix_struct, name="low-rank base work")
662 CALL cp_fm_uplo_to_full(preconditioner_env%base_fm, work)
663 CALL cp_fm_release(work)
664 ELSE
665 CALL cp_fm_set_all(preconditioner_env%base_fm, alpha=0.0_dp, beta=1.0_dp)
666 END IF
667
668 END SUBROUTINE prepare_low_rank_base
669
670! **************************************************************************************************
671
672! creates a preconditioner for the system (H-energy_homo S)
673! this preconditioner is (must be) symmetric positive definite.
674! currently uses a atom-block-diagonal form
675! each block will be ....
676! might overwrite matrix_h, matrix_t
677
678! **************************************************************************************************
679!> \brief ...
680!> \param preconditioner_env ...
681!> \param precon_type ...
682!> \param solver_type ...
683!> \param matrix_h ...
684!> \param matrix_s ...
685!> \param matrix_t ...
686!> \param mo_set ...
687!> \param energy_gap ...
688!> \param convert_precond_to_dbcsr ...
689!> \param chol_type ...
690!> \par History
691!> 09.2014 removed some unused or unfinished methods
692!> removed sparse preconditioners and the
693!> sparse approximate inverse at rev 14341 [Florian Schiffmann]
694! **************************************************************************************************
695 SUBROUTINE make_preconditioner(preconditioner_env, precon_type, solver_type, matrix_h, matrix_s, &
696 matrix_t, mo_set, energy_gap, convert_precond_to_dbcsr, chol_type)
697
698 TYPE(preconditioner_type) :: preconditioner_env
699 INTEGER, INTENT(IN) :: precon_type, solver_type
700 TYPE(dbcsr_type), POINTER :: matrix_h
701 TYPE(dbcsr_type), OPTIONAL, POINTER :: matrix_s, matrix_t
702 TYPE(mo_set_type), INTENT(IN) :: mo_set
703 REAL(kind=dp) :: energy_gap
704 LOGICAL, INTENT(IN), OPTIONAL :: convert_precond_to_dbcsr
705 INTEGER, INTENT(IN), OPTIONAL :: chol_type
706
707 CHARACTER(len=*), PARAMETER :: routinen = 'make_preconditioner'
708
709 INTEGER :: handle, k, max_rank, my_solver_type, &
710 nao, nhomo, output_unit
711 LOGICAL :: my_convert_precond_to_dbcsr, &
712 needs_full_spectrum, needs_homo, &
713 use_mo_coeff_b
714 REAL(kind=dp) :: energy_homo
715 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: eigenvalues_ot
716 TYPE(cp_fm_struct_type), POINTER :: fm_struct
717 TYPE(cp_fm_type) :: mo_occ
718 TYPE(cp_fm_type), POINTER :: mo_coeff
719 TYPE(cp_logger_type), POINTER :: logger
720 TYPE(dbcsr_type), POINTER :: mo_coeff_b
721
722 CALL timeset(routinen, handle)
723
724 CALL get_mo_set(mo_set=mo_set, mo_coeff=mo_coeff, mo_coeff_b=mo_coeff_b, homo=nhomo)
725 use_mo_coeff_b = mo_set%use_mo_coeff_b
726 CALL cp_fm_get_info(mo_coeff, ncol_global=k, nrow_global=nao)
727
728 ! Starting some matrix mess, check where to store the result in preconditioner_env, fm or dbcsr_matrix
729 my_convert_precond_to_dbcsr = .false.
730 IF (PRESENT(convert_precond_to_dbcsr)) my_convert_precond_to_dbcsr = convert_precond_to_dbcsr
731 IF (precon_type == ot_precond_fermi_low_rank .AND. &
732 preconditioner_env%low_rank_base == ot_low_rank_base_lattice_fft) THEN
733 ! The FFT base is an operator rather than a stored DBCSR matrix. Keep the
734 ! complementary vectors in their distributed full-matrix representation.
735 my_convert_precond_to_dbcsr = .false.
736 END IF
737
738 ! Thanks to the mess with the matrices we need to make sure in this case that the
739 ! Previous inverse is properly stored as a sparse matrix, fm gets deallocated here
740 ! if it wasn't anyway
741 IF (preconditioner_env%solver == ot_precond_solver_update) THEN
742 CALL transfer_fm_to_dbcsr(preconditioner_env%fm, preconditioner_env%dbcsr_matrix, matrix_h)
743 END IF
744
745 needs_full_spectrum = .false.
746 needs_homo = .false.
747
748 SELECT CASE (precon_type)
750 needs_full_spectrum = .true.
751 ! both of them need the coefficients as fm's, more matrix mess
752 IF (use_mo_coeff_b) THEN
753 CALL copy_dbcsr_to_fm(mo_coeff_b, mo_coeff)
754 END IF
756 needs_homo = .true.
757 ! XXXX to be removed if homo estimate only is implemented
758 needs_full_spectrum = .true.
760 ! these should be happy without an estimate for the homo energy
761 ! preconditioning can not depend on an absolute eigenvalue, only on eigenvalue differences
762 CASE DEFAULT
763 cpabort("The preconditioner is unknown ...")
764 END SELECT
765
766 ALLOCATE (eigenvalues_ot(k))
767 energy_homo = 0.0_dp
768 IF (needs_full_spectrum) THEN
769 ! XXXXXXXXXXXXXXXX do not touch the initial MOs, could be harmful for either
770 ! the case of non-equivalent MOs but also for the derivate
771 ! we could already have all eigenvalues e.g. full_all and we could skip this
772 ! to be optimised later.
773 ! one flaw is that not all SCF methods (i.e. that go over mo_derivs directly)
774 ! have a 'valid' matrix_h... (we even don't know what evals are in that case)
775 IF (use_mo_coeff_b) THEN
776 CALL calculate_subspace_eigenvalues(mo_coeff_b, matrix_h, &
777 eigenvalues_ot, do_rotation=.false., &
778 para_env=mo_coeff%matrix_struct%para_env, &
779 blacs_env=mo_coeff%matrix_struct%context)
780 ELSE
781 CALL calculate_subspace_eigenvalues(mo_coeff, matrix_h, &
782 eigenvalues_ot, do_rotation=.false.)
783 END IF
784 IF (k > 0) THEN
785 cpassert(nhomo > 0 .AND. nhomo <= k)
786 energy_homo = eigenvalues_ot(nhomo)
787 END IF
788 ELSE
789 IF (needs_homo) THEN
790 cpabort("Not yet implemented")
791 END IF
792 END IF
793
794 ! After all bits and pieces of checking and initialization, here comes the
795 ! part where the preconditioner matrix gets created and solved.
796 ! This will give the matrices for later use
797 my_solver_type = solver_type
798 preconditioner_env%in_use = precon_type
799 preconditioner_env%cholesky_use = cholesky_reduce
800 IF (PRESENT(chol_type)) preconditioner_env%cholesky_use = chol_type
801 preconditioner_env%in_use = precon_type
802 preconditioner_env%lattice_fft_active = .false.
803 preconditioner_env%lattice_corrected_projection_error = 0.0_dp
804 IF (ALLOCATED(preconditioner_env%lattice_inverse_k)) THEN
805 DEALLOCATE (preconditioner_env%lattice_inverse_k)
806 END IF
807 IF (ALLOCATED(preconditioner_env%lattice_state_inverse_k)) THEN
808 DEALLOCATE (preconditioner_env%lattice_state_inverse_k)
809 END IF
810 IF (ALLOCATED(preconditioner_env%lattice_occ_vectors)) THEN
811 DEALLOCATE (preconditioner_env%lattice_occ_vectors)
812 END IF
813 IF (ALLOCATED(preconditioner_env%lattice_occ_dual)) THEN
814 DEALLOCATE (preconditioner_env%lattice_occ_dual)
815 END IF
816 IF (ALLOCATED(preconditioner_env%lattice_local_indices)) THEN
817 DEALLOCATE (preconditioner_env%lattice_local_indices)
818 END IF
819 IF (ALLOCATED(preconditioner_env%lattice_local_operator)) THEN
820 DEALLOCATE (preconditioner_env%lattice_local_operator)
821 END IF
822 IF (ALLOCATED(preconditioner_env%lattice_local_coarse_inverse)) THEN
823 DEALLOCATE (preconditioner_env%lattice_local_coarse_inverse)
824 END IF
825 IF (ALLOCATED(preconditioner_env%lattice_state_local_indices)) THEN
826 DEALLOCATE (preconditioner_env%lattice_state_local_indices)
827 END IF
828 IF (ALLOCATED(preconditioner_env%lattice_state_local_operator)) THEN
829 DEALLOCATE (preconditioner_env%lattice_state_local_operator)
830 END IF
831 IF (ALLOCATED(preconditioner_env%lattice_state_local_coarse_inverse)) THEN
832 DEALLOCATE (preconditioner_env%lattice_state_local_coarse_inverse)
833 END IF
834 IF (nhomo == k) THEN
835 CALL make_preconditioner_matrix(preconditioner_env, matrix_h, matrix_s, matrix_t, mo_coeff, &
836 energy_homo, eigenvalues_ot, energy_gap, my_solver_type)
837 ELSE
838 CALL cp_fm_struct_create(fm_struct, nrow_global=nao, ncol_global=nhomo, &
839 context=preconditioner_env%ctxt, &
840 para_env=preconditioner_env%para_env)
841 CALL cp_fm_create(mo_occ, fm_struct)
842 CALL cp_fm_to_fm(mo_coeff, mo_occ, nhomo)
843 CALL cp_fm_struct_release(fm_struct)
844 !
845 CALL make_preconditioner_matrix(preconditioner_env, matrix_h, matrix_s, matrix_t, mo_occ, &
846 energy_homo, eigenvalues_ot(1:nhomo), energy_gap, my_solver_type)
847 !
848 CALL cp_fm_release(mo_occ)
849 END IF
850
851 IF (preconditioner_env%lattice_fft_mode /= ot_lattice_fft_off) THEN
852 SELECT CASE (precon_type)
854 cpassert(ASSOCIATED(preconditioner_env%sparse_matrix))
855 CALL prepare_lattice_operator(preconditioner_env, preconditioner_env%sparse_matrix, &
856 "FULL_KINETIC", nhomo)
858 cpassert(ASSOCIATED(preconditioner_env%sparse_matrix))
859 CALL prepare_lattice_operator(preconditioner_env, preconditioner_env%sparse_matrix, &
860 "FULL_SINGLE_INVERSE", nhomo)
862 cpassert(ASSOCIATED(preconditioner_env%sparse_matrix))
863 CALL prepare_lattice_operator(preconditioner_env, preconditioner_env%sparse_matrix, &
864 "FULL_S_INVERSE", nhomo)
866 cpassert(PRESENT(matrix_s))
867 cpassert(ASSOCIATED(matrix_s))
868 CALL prepare_lattice_state_operator(preconditioner_env, matrix_s, &
869 "FULL_ALL_COVARIANT")
870 END SELECT
871 END IF
872
873 IF (.NOT. preconditioner_env%lattice_fft_active) THEN
874 CALL solve_preconditioner(my_solver_type, preconditioner_env, matrix_s, matrix_h)
875 ELSE
876 preconditioner_env%solver = ot_precond_solver_inv_chol
877 END IF
878
879 IF (precon_type == ot_precond_fermi_low_rank) THEN
880 preconditioner_env%spectral_reference = energy_homo
881 preconditioner_env%spectral_window = &
882 max(fermi_low_rank_min_window, fermi_low_rank_gap_factor*preconditioner_env%energy_gap)
883 preconditioner_env%spectral_base_scale = 1.0_dp/preconditioner_env%spectral_window
884 max_rank = preconditioner_env%spectral_max_rank
885 IF (max_rank == -1) max_rank = SIZE(preconditioner_env%full_evals) - nhomo
886 preconditioner_env%spectral_rank = low_rank_select_rank( &
887 preconditioner_env%full_evals, nhomo, max_rank, &
888 fermi_low_rank_degeneracy_tolerance)
889 CALL prepare_low_rank_base(preconditioner_env, matrix_s)
890 IF (preconditioner_env%lattice_fft_active) THEN
891 ! With an approximate base, retain only non-negative spectral updates.
892 ! This preserves positive definiteness independently of the projection error.
893 DO WHILE (preconditioner_env%spectral_rank > 0)
895 preconditioner_env%full_evals(nhomo + preconditioner_env%spectral_rank), &
896 preconditioner_env%spectral_reference, preconditioner_env%energy_gap) >= &
897 preconditioner_env%spectral_base_scale) EXIT
898 preconditioner_env%spectral_rank = preconditioner_env%spectral_rank - 1
899 END DO
900 END IF
901 IF (preconditioner_env%para_env%is_source()) THEN
902 logger => cp_get_default_logger()
903 output_unit = cp_logger_get_default_io_unit(logger)
904 IF (output_unit > 0) THEN
905 WRITE (output_unit, '(T2,A,F14.8,A,I0)') &
906 "OT| FERMI_LOW_RANK common reference [hartree]: ", &
907 preconditioner_env%spectral_reference, " retained rank: ", &
908 preconditioner_env%spectral_rank
909 END IF
910 END IF
911 END IF
912
913 ! Here comes more matrix mess, make sure to output the correct matrix format,
914 ! A bit pointless to convert the cholesky factorized version as it doesn't work in
915 ! dbcsr form and will crash later,...
916 IF (.NOT. preconditioner_env%lattice_fft_active) THEN
917 IF (my_convert_precond_to_dbcsr) THEN
918 CALL transfer_fm_to_dbcsr(preconditioner_env%fm, preconditioner_env%dbcsr_matrix, matrix_h)
919 IF (precon_type == ot_precond_fermi_low_rank .AND. &
920 preconditioner_env%low_rank_base == ot_low_rank_base_overlap) THEN
921 CALL transfer_fm_to_dbcsr(preconditioner_env%base_fm, &
922 preconditioner_env%base_dbcsr_matrix, matrix_h)
923 END IF
924 ELSE
925 CALL transfer_dbcsr_to_fm(preconditioner_env%dbcsr_matrix, preconditioner_env%fm, &
926 preconditioner_env%para_env, preconditioner_env%ctxt)
927 END IF
928 END IF
929
930 DEALLOCATE (eigenvalues_ot)
931
932 CALL timestop(handle)
933
934 END SUBROUTINE make_preconditioner
935
936! **************************************************************************************************
937!> \brief Construct FULL_ALL directly from one complex H(k), S(k), and C(k) channel.
938!> \param preconditioner_env preconditioner storage
939!> \param matrix_c_re real part of the active OT reference orbitals
940!> \param matrix_c_im imaginary part of the active OT reference orbitals
941!> \param matrix_h_re real part of H(k)
942!> \param matrix_h_im imaginary part of H(k)
943!> \param matrix_s_re real part of S(k)
944!> \param matrix_s_im imaginary part of S(k)
945!> \param mo_set ...
946!> \param energy_gap denominator floor
947!> \param solver_type inverse application selected in the OT input
948! **************************************************************************************************
949 SUBROUTINE make_preconditioner_complex_full_all(preconditioner_env, &
950 matrix_c_re, matrix_c_im, &
951 matrix_h_re, matrix_h_im, &
952 matrix_s_re, matrix_s_im, &
953 mo_set, energy_gap, solver_type)
954
955 TYPE(preconditioner_type) :: preconditioner_env
956 TYPE(dbcsr_type), POINTER :: matrix_c_re, matrix_c_im, matrix_h_re, &
957 matrix_h_im, matrix_s_re, matrix_s_im
958 TYPE(mo_set_type), INTENT(IN) :: mo_set
959 REAL(kind=dp), INTENT(IN) :: energy_gap
960 INTEGER, INTENT(IN) :: solver_type
961
962 CHARACTER(len=*), PARAMETER :: routinen = 'make_preconditioner_complex_full_all'
963
964 INTEGER :: handle, k
965 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues
966 TYPE(cp_cfm_type) :: matrix_c, matrix_h, matrix_s
967
968 CALL timeset(routinen, handle)
969
970 NULLIFY (eigenvalues)
971 cpassert(ASSOCIATED(matrix_c_re))
972 cpassert(ASSOCIATED(matrix_c_im))
973 cpassert(ASSOCIATED(matrix_h_re))
974 cpassert(ASSOCIATED(matrix_h_im))
975 cpassert(ASSOCIATED(matrix_s_re))
976 cpassert(ASSOCIATED(matrix_s_im))
977 IF (solver_type /= ot_precond_solver_default) THEN
978 cpabort('Complex FULL_ALL supports only PRECOND_SOLVER DEFAULT')
979 END IF
980 CALL dbcsr_get_info(matrix_c_re, nfullcols_total=k)
981 cpassert(k > 0)
982 CALL get_mo_set(mo_set, eigenvalues=eigenvalues)
983 cpassert(ASSOCIATED(eigenvalues))
984 cpassert(SIZE(eigenvalues) >= k)
985
986 CALL dbcsr_pair_to_cfm(matrix_c_re, matrix_c_im, preconditioner_env, &
987 'complex FULL_ALL C', matrix_c)
988 CALL dbcsr_pair_to_cfm(matrix_h_re, matrix_h_im, preconditioner_env, &
989 'complex FULL_ALL H', matrix_h)
990 CALL dbcsr_pair_to_cfm(matrix_s_re, matrix_s_im, preconditioner_env, &
991 'complex FULL_ALL S', matrix_s)
992 CALL make_complex_full_all(preconditioner_env, matrix_c, matrix_h, matrix_s, &
993 eigenvalues(1:k), energy_gap)
994
995 CALL cp_cfm_release(matrix_c)
996 CALL cp_cfm_release(matrix_s)
997 CALL cp_cfm_release(matrix_h)
998
999 CALL timestop(handle)
1000
1002
1003! **************************************************************************************************
1004!> \brief Construct a complex FULL_SINGLE preconditioner from H(k) and S(k).
1005!> \param preconditioner_env preconditioner storage
1006!> \param matrix_h_re real part of H(k)
1007!> \param matrix_h_im imaginary part of H(k)
1008!> \param matrix_s_re real part of S(k)
1009!> \param matrix_s_im imaginary part of S(k)
1010!> \param mo_set orbital energy labels defining the occupied edge
1011!> \param energy_gap denominator floor
1012!> \param solver_type requested inverse solver
1013! **************************************************************************************************
1014 SUBROUTINE make_preconditioner_complex_full_single(preconditioner_env, &
1015 matrix_h_re, matrix_h_im, &
1016 matrix_s_re, matrix_s_im, &
1017 mo_set, energy_gap, solver_type)
1018
1019 TYPE(preconditioner_type) :: preconditioner_env
1020 TYPE(dbcsr_type), POINTER :: matrix_h_re, matrix_h_im, matrix_s_re, &
1021 matrix_s_im
1022 TYPE(mo_set_type), INTENT(IN) :: mo_set
1023 REAL(kind=dp), INTENT(IN) :: energy_gap
1024 INTEGER, INTENT(IN) :: solver_type
1025
1026 INTEGER :: homo
1027 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues
1028 TYPE(cp_cfm_type) :: matrix_h, matrix_s
1029
1030 NULLIFY (eigenvalues)
1031 cpassert(ASSOCIATED(matrix_h_re))
1032 cpassert(ASSOCIATED(matrix_h_im))
1033 cpassert(ASSOCIATED(matrix_s_re))
1034 cpassert(ASSOCIATED(matrix_s_im))
1035 IF (solver_type /= ot_precond_solver_default) THEN
1036 cpabort('Complex FULL_SINGLE supports only PRECOND_SOLVER DEFAULT')
1037 END IF
1038 CALL get_mo_set(mo_set, homo=homo, eigenvalues=eigenvalues)
1039 cpassert(ASSOCIATED(eigenvalues))
1040 cpassert(homo > 0 .AND. homo <= SIZE(eigenvalues))
1041
1042 CALL dbcsr_pair_to_cfm(matrix_h_re, matrix_h_im, preconditioner_env, &
1043 'complex FULL_SINGLE H', matrix_h)
1044 CALL dbcsr_pair_to_cfm(matrix_s_re, matrix_s_im, preconditioner_env, &
1045 'complex FULL_SINGLE S', matrix_s)
1046 CALL make_complex_full_single(preconditioner_env, matrix_h, matrix_s, &
1047 eigenvalues(homo), energy_gap)
1048 CALL cp_cfm_release(matrix_s)
1049 CALL cp_cfm_release(matrix_h)
1050
1052
1053! **************************************************************************************************
1054!> \brief Construct FERMI_LOW_RANK directly for one complex k-point channel.
1055!> \param preconditioner_env preconditioner storage
1056!> \param matrix_c_re real part of the active OT reference orbitals
1057!> \param matrix_c_im imaginary part of the active OT reference orbitals
1058!> \param matrix_h_re real part of H(k)
1059!> \param matrix_h_im imaginary part of H(k)
1060!> \param matrix_s_re real part of S(k)
1061!> \param matrix_s_im imaginary part of S(k)
1062!> \param energy_gap denominator floor
1063!> \param max_rank maximum complementary rank, or -1 for the full complement
1064!> \param solver_type inverse application selected in the OT input
1065! **************************************************************************************************
1067 preconditioner_env, matrix_c_re, matrix_c_im, matrix_h_re, matrix_h_im, &
1068 matrix_s_re, matrix_s_im, energy_gap, max_rank, solver_type)
1069
1070 TYPE(preconditioner_type) :: preconditioner_env
1071 TYPE(dbcsr_type), POINTER :: matrix_c_re, matrix_c_im, matrix_h_re, &
1072 matrix_h_im, matrix_s_re, matrix_s_im
1073 REAL(kind=dp), INTENT(IN) :: energy_gap
1074 INTEGER, INTENT(IN) :: max_rank, solver_type
1075
1076 TYPE(cp_cfm_type) :: matrix_c, matrix_h, matrix_s
1077
1078 IF (solver_type /= ot_precond_solver_default) THEN
1079 cpabort('Complex FERMI_LOW_RANK supports only PRECOND_SOLVER DEFAULT')
1080 END IF
1081
1082 CALL dbcsr_pair_to_cfm(matrix_c_re, matrix_c_im, preconditioner_env, &
1083 'complex FERMI_LOW_RANK C', matrix_c)
1084 CALL dbcsr_pair_to_cfm(matrix_h_re, matrix_h_im, preconditioner_env, &
1085 'complex FERMI_LOW_RANK H', matrix_h)
1086 CALL dbcsr_pair_to_cfm(matrix_s_re, matrix_s_im, preconditioner_env, &
1087 'complex FERMI_LOW_RANK S', matrix_s)
1089 preconditioner_env, matrix_c, matrix_h, matrix_s, energy_gap, &
1090 max_rank, max(fermi_low_rank_min_window, fermi_low_rank_gap_factor*energy_gap), &
1091 fermi_low_rank_degeneracy_tolerance)
1092
1093 CALL cp_cfm_release(matrix_s)
1094 CALL cp_cfm_release(matrix_h)
1095 CALL cp_cfm_release(matrix_c)
1096
1098
1099! **************************************************************************************************
1100!> \brief Construct rotation-covariant FULL_ALL for one complex k-point channel.
1101!> \param preconditioner_env preconditioner storage
1102!> \param matrix_c_re real part of the active OT reference orbitals
1103!> \param matrix_c_im imaginary part of the active OT reference orbitals
1104!> \param matrix_h_re real part of H(k)
1105!> \param matrix_h_im imaginary part of H(k)
1106!> \param matrix_s_re real part of S(k)
1107!> \param matrix_s_im imaginary part of S(k)
1108!> \param energy_gap denominator floor
1109!> \param solver_type inverse application selected in the OT input
1110!> \param occupation_signature occupations that define the admissible private gauge
1111! **************************************************************************************************
1113 preconditioner_env, matrix_c_re, matrix_c_im, matrix_h_re, matrix_h_im, &
1114 matrix_s_re, matrix_s_im, energy_gap, solver_type, occupation_signature)
1115
1116 TYPE(preconditioner_type) :: preconditioner_env
1117 TYPE(dbcsr_type), POINTER :: matrix_c_re, matrix_c_im, matrix_h_re, &
1118 matrix_h_im, matrix_s_re, matrix_s_im
1119 REAL(kind=dp), INTENT(IN) :: energy_gap
1120 INTEGER, INTENT(IN) :: solver_type
1121 REAL(kind=dp), DIMENSION(:, :), INTENT(IN), &
1122 OPTIONAL :: occupation_signature
1123
1124 TYPE(cp_cfm_type) :: matrix_c, matrix_h, matrix_s
1125
1126 IF (solver_type /= ot_precond_solver_default) THEN
1127 cpabort('Complex FULL_ALL_COVARIANT supports only PRECOND_SOLVER DEFAULT')
1128 END IF
1129 CALL dbcsr_pair_to_cfm(matrix_c_re, matrix_c_im, preconditioner_env, &
1130 'complex FULL_ALL_COVARIANT C', matrix_c)
1131 CALL dbcsr_pair_to_cfm(matrix_h_re, matrix_h_im, preconditioner_env, &
1132 'complex FULL_ALL_COVARIANT H', matrix_h)
1133 CALL dbcsr_pair_to_cfm(matrix_s_re, matrix_s_im, preconditioner_env, &
1134 'complex FULL_ALL_COVARIANT S', matrix_s)
1135 CALL make_complex_full_all_covariant(preconditioner_env, matrix_c, matrix_h, matrix_s, &
1136 energy_gap, occupation_signature)
1137
1138 CALL cp_cfm_release(matrix_s)
1139 CALL cp_cfm_release(matrix_h)
1140 CALL cp_cfm_release(matrix_c)
1141
1143
1144! **************************************************************************************************
1145!> \brief Construct a complex FULL_SINGLE_INVERSE preconditioner without discarding Im(H,S,C).
1146!> \param preconditioner_env preconditioner storage
1147!> \param matrix_c_re real part of the occupied reference orbitals
1148!> \param matrix_c_im imaginary part of the occupied reference orbitals
1149!> \param matrix_h_re real part of H(k)
1150!> \param matrix_h_im imaginary part of H(k)
1151!> \param matrix_s_re real part of S(k)
1152!> \param matrix_s_im imaginary part of S(k)
1153!> \param energy_gap lower spectral bound
1154!> \param solver_type requested inverse solver
1155! **************************************************************************************************
1157 matrix_c_re, matrix_c_im, &
1158 matrix_h_re, matrix_h_im, &
1159 matrix_s_re, matrix_s_im, &
1160 energy_gap, solver_type)
1161
1162 TYPE(preconditioner_type) :: preconditioner_env
1163 TYPE(dbcsr_type), POINTER :: matrix_c_re, matrix_c_im, matrix_h_re, &
1164 matrix_h_im, matrix_s_re, matrix_s_im
1165 REAL(kind=dp), INTENT(IN) :: energy_gap
1166 INTEGER, INTENT(IN) :: solver_type
1167
1168 INTEGER :: k, n
1169 TYPE(cp_cfm_type) :: matrix_c, matrix_h, matrix_s
1170 TYPE(cp_fm_struct_type), POINTER :: fm_struct_mo, fm_struct_square
1171 TYPE(cp_fm_type) :: matrix_c_im_fm, matrix_c_re_fm, &
1172 matrix_h_im_fm, matrix_h_re_fm, &
1173 matrix_s_im_fm, matrix_s_re_fm
1174
1175 NULLIFY (fm_struct_mo, fm_struct_square)
1176 cpassert(ASSOCIATED(matrix_c_re))
1177 cpassert(ASSOCIATED(matrix_c_im))
1178 cpassert(ASSOCIATED(matrix_h_re))
1179 cpassert(ASSOCIATED(matrix_h_im))
1180 cpassert(ASSOCIATED(matrix_s_re))
1181 cpassert(ASSOCIATED(matrix_s_im))
1182 IF (solver_type /= ot_precond_solver_default .AND. &
1183 solver_type /= ot_precond_solver_inv_chol) THEN
1184 cpabort("Complex FULL_SINGLE_INVERSE supports only PRECOND_SOLVER DEFAULT/INVERSE_CHOLESKY")
1185 END IF
1186 CALL dbcsr_get_info(matrix_c_re, nfullrows_total=n, nfullcols_total=k)
1187 cpassert(n > 0 .AND. k > 0 .AND. k <= n)
1188
1189 CALL cp_fm_struct_create(fm_struct_mo, nrow_global=n, ncol_global=k, &
1190 context=preconditioner_env%ctxt, &
1191 para_env=preconditioner_env%para_env)
1192 CALL cp_fm_struct_create(fm_struct_square, nrow_global=n, ncol_global=n, &
1193 context=preconditioner_env%ctxt, &
1194 para_env=preconditioner_env%para_env)
1195 CALL cp_fm_create(matrix_c_re_fm, fm_struct_mo, name='complex preconditioner C real')
1196 CALL cp_fm_create(matrix_c_im_fm, fm_struct_mo, name='complex preconditioner C imaginary')
1197 CALL cp_fm_create(matrix_h_re_fm, fm_struct_square, name='complex preconditioner H real')
1198 CALL cp_fm_create(matrix_h_im_fm, fm_struct_square, name='complex preconditioner H imaginary')
1199 CALL cp_fm_create(matrix_s_re_fm, fm_struct_square, name='complex preconditioner S real')
1200 CALL cp_fm_create(matrix_s_im_fm, fm_struct_square, name='complex preconditioner S imaginary')
1201 CALL cp_fm_struct_release(fm_struct_mo)
1202 CALL cp_fm_struct_release(fm_struct_square)
1203 CALL copy_dbcsr_to_fm(matrix_c_re, matrix_c_re_fm)
1204 CALL copy_dbcsr_to_fm(matrix_c_im, matrix_c_im_fm)
1205 CALL copy_dbcsr_to_fm(matrix_h_re, matrix_h_re_fm)
1206 CALL copy_dbcsr_to_fm(matrix_h_im, matrix_h_im_fm)
1207 CALL copy_dbcsr_to_fm(matrix_s_re, matrix_s_re_fm)
1208 CALL copy_dbcsr_to_fm(matrix_s_im, matrix_s_im_fm)
1209
1210 CALL cp_cfm_create(matrix_c, matrix_c_re_fm%matrix_struct, name='complex preconditioner C')
1211 CALL cp_cfm_create(matrix_h, matrix_h_re_fm%matrix_struct, name='complex preconditioner H')
1212 CALL cp_cfm_create(matrix_s, matrix_s_re_fm%matrix_struct, name='complex preconditioner S')
1213 CALL cp_fm_to_cfm(matrix_c_re_fm, matrix_c_im_fm, matrix_c)
1214 CALL cp_fm_to_cfm(matrix_h_re_fm, matrix_h_im_fm, matrix_h)
1215 CALL cp_fm_to_cfm(matrix_s_re_fm, matrix_s_im_fm, matrix_s)
1216 CALL make_complex_full_single_inverse(preconditioner_env, matrix_c, matrix_h, matrix_s, &
1217 energy_gap)
1218 preconditioner_env%solver = solver_type
1219
1220 CALL cp_cfm_release(matrix_s)
1221 CALL cp_cfm_release(matrix_h)
1222 CALL cp_cfm_release(matrix_c)
1223 CALL cp_fm_release(matrix_s_im_fm)
1224 CALL cp_fm_release(matrix_s_re_fm)
1225 CALL cp_fm_release(matrix_h_im_fm)
1226 CALL cp_fm_release(matrix_h_re_fm)
1227 CALL cp_fm_release(matrix_c_im_fm)
1228 CALL cp_fm_release(matrix_c_re_fm)
1229
1231
1232! **************************************************************************************************
1233!> \brief Construct a complex FULL_S_INVERSE preconditioner.
1234!> \param preconditioner_env preconditioner storage
1235!> \param matrix_s_re real part of S(k)
1236!> \param matrix_s_im imaginary part of S(k)
1237!> \param solver_type requested inverse solver
1238! **************************************************************************************************
1239 SUBROUTINE make_preconditioner_complex_full_s_inverse(preconditioner_env, &
1240 matrix_s_re, matrix_s_im, solver_type)
1241
1242 TYPE(preconditioner_type) :: preconditioner_env
1243 TYPE(dbcsr_type), POINTER :: matrix_s_re, matrix_s_im
1244 INTEGER, INTENT(IN) :: solver_type
1245
1246 TYPE(cp_cfm_type) :: matrix_s
1247
1248 IF (solver_type /= ot_precond_solver_default .AND. &
1249 solver_type /= ot_precond_solver_inv_chol) THEN
1250 cpabort("Complex FULL_S_INVERSE supports only PRECOND_SOLVER DEFAULT/INVERSE_CHOLESKY")
1251 END IF
1252 CALL dbcsr_pair_to_cfm(matrix_s_re, matrix_s_im, preconditioner_env, &
1253 'complex preconditioner S', matrix_s)
1254 CALL make_complex_full_s_inverse(preconditioner_env, matrix_s)
1255 preconditioner_env%solver = solver_type
1256 CALL cp_cfm_release(matrix_s)
1257
1259
1260! **************************************************************************************************
1261!> \brief Construct a complex FULL_KINETIC preconditioner.
1262!> \param preconditioner_env preconditioner storage
1263!> \param matrix_t_re real part of T(k)
1264!> \param matrix_t_im imaginary part of T(k)
1265!> \param matrix_s_re real part of S(k)
1266!> \param matrix_s_im imaginary part of S(k)
1267!> \param energy_gap non-negative overlap shift
1268!> \param solver_type requested inverse solver
1269! **************************************************************************************************
1270 SUBROUTINE make_preconditioner_complex_full_kinetic(preconditioner_env, &
1271 matrix_t_re, matrix_t_im, &
1272 matrix_s_re, matrix_s_im, &
1273 energy_gap, solver_type)
1274
1275 TYPE(preconditioner_type) :: preconditioner_env
1276 TYPE(dbcsr_type), POINTER :: matrix_t_re, matrix_t_im, matrix_s_re, &
1277 matrix_s_im
1278 REAL(kind=dp), INTENT(IN) :: energy_gap
1279 INTEGER, INTENT(IN) :: solver_type
1280
1281 TYPE(cp_cfm_type) :: matrix_s, matrix_t
1282
1283 IF (solver_type /= ot_precond_solver_default .AND. &
1284 solver_type /= ot_precond_solver_inv_chol) THEN
1285 cpabort("Complex FULL_KINETIC supports only PRECOND_SOLVER DEFAULT/INVERSE_CHOLESKY")
1286 END IF
1287 CALL dbcsr_pair_to_cfm(matrix_t_re, matrix_t_im, preconditioner_env, &
1288 'complex preconditioner T', matrix_t)
1289 CALL dbcsr_pair_to_cfm(matrix_s_re, matrix_s_im, preconditioner_env, &
1290 'complex preconditioner S', matrix_s)
1291 CALL make_complex_full_kinetic(preconditioner_env, matrix_t, matrix_s, energy_gap)
1292 preconditioner_env%solver = solver_type
1293 CALL cp_cfm_release(matrix_s)
1294 CALL cp_cfm_release(matrix_t)
1295
1297
1298! **************************************************************************************************
1299!> \brief Copy a real/imaginary DBCSR pair to one distributed complex full matrix.
1300!> \param matrix_re real matrix component
1301!> \param matrix_im imaginary matrix component
1302!> \param preconditioner_env source of the distribution context
1303!> \param matrix_name matrix label
1304!> \param matrix complex output matrix
1305! **************************************************************************************************
1306 SUBROUTINE dbcsr_pair_to_cfm(matrix_re, matrix_im, preconditioner_env, matrix_name, matrix)
1307
1308 TYPE(dbcsr_type), POINTER :: matrix_re, matrix_im
1309 TYPE(preconditioner_type) :: preconditioner_env
1310 CHARACTER(LEN=*), INTENT(IN) :: matrix_name
1311 TYPE(cp_cfm_type), INTENT(OUT) :: matrix
1312
1313 INTEGER :: ncol, nrow
1314 TYPE(cp_fm_struct_type), POINTER :: fm_struct
1315 TYPE(cp_fm_type) :: matrix_im_fm, matrix_re_fm
1316
1317 NULLIFY (fm_struct)
1318 cpassert(ASSOCIATED(matrix_re))
1319 cpassert(ASSOCIATED(matrix_im))
1320 CALL dbcsr_get_info(matrix_re, nfullrows_total=nrow, nfullcols_total=ncol)
1321 CALL cp_fm_struct_create(fm_struct, nrow_global=nrow, ncol_global=ncol, &
1322 context=preconditioner_env%ctxt, &
1323 para_env=preconditioner_env%para_env)
1324 CALL cp_fm_create(matrix_re_fm, fm_struct, name=trim(matrix_name)//' real')
1325 CALL cp_fm_create(matrix_im_fm, fm_struct, name=trim(matrix_name)//' imaginary')
1326 CALL cp_fm_struct_release(fm_struct)
1327 CALL copy_dbcsr_to_fm(matrix_re, matrix_re_fm)
1328 CALL copy_dbcsr_to_fm(matrix_im, matrix_im_fm)
1329 CALL cp_cfm_create(matrix, matrix_re_fm%matrix_struct, name=matrix_name)
1330 CALL cp_fm_to_cfm(matrix_re_fm, matrix_im_fm, matrix)
1331 CALL cp_fm_release(matrix_im_fm)
1332 CALL cp_fm_release(matrix_re_fm)
1333
1334 END SUBROUTINE dbcsr_pair_to_cfm
1335
1336! **************************************************************************************************
1337!> \brief Allows for a restart of the preconditioner
1338!> depending on the method it purges all arrays or keeps them
1339!> \param qs_env ...
1340!> \param preconditioner ...
1341!> \param prec_type ...
1342!> \param nspins ...
1343! **************************************************************************************************
1344 SUBROUTINE restart_preconditioner(qs_env, preconditioner, prec_type, nspins)
1345
1346 TYPE(qs_environment_type), POINTER :: qs_env
1347 TYPE(preconditioner_p_type), DIMENSION(:), POINTER :: preconditioner
1348 INTEGER, INTENT(IN) :: prec_type, nspins
1349
1350 INTEGER :: ispin
1351 TYPE(cp_blacs_env_type), POINTER :: blacs_env
1352 TYPE(mp_para_env_type), POINTER :: para_env
1353
1354 NULLIFY (para_env, blacs_env)
1355 CALL get_qs_env(qs_env, para_env=para_env, blacs_env=blacs_env)
1356
1357 IF (ASSOCIATED(preconditioner)) THEN
1358 SELECT CASE (prec_type)
1361 ot_precond_full_single) ! these depend on the ks matrix
1362 DO ispin = 1, SIZE(preconditioner)
1363 CALL destroy_preconditioner(preconditioner(ispin)%preconditioner)
1364 DEALLOCATE (preconditioner(ispin)%preconditioner)
1365 END DO
1366 DEALLOCATE (preconditioner)
1368 ot_precond_full_single_inverse) ! these are 'independent'
1369 ! do nothing
1370 CASE DEFAULT
1371 cpabort("Unknown preconditioner type")
1372 END SELECT
1373 END IF
1374
1375 ! add an OT preconditioner if none is present
1376 IF (.NOT. ASSOCIATED(preconditioner)) THEN
1377 SELECT CASE (prec_type)
1380 ALLOCATE (preconditioner(nspins))
1381 CASE DEFAULT
1382 ALLOCATE (preconditioner(1))
1383 END SELECT
1384 DO ispin = 1, SIZE(preconditioner)
1385 ALLOCATE (preconditioner(ispin)%preconditioner)
1386 CALL init_preconditioner(preconditioner(ispin)%preconditioner, &
1387 para_env=para_env, &
1388 blacs_env=blacs_env)
1389 END DO
1390 END IF
1391
1392 END SUBROUTINE restart_preconditioner
1393
1394! **************************************************************************************************
1395!> \brief ...
1396!> \param qs_env ...
1397!> \param mos ...
1398!> \param matrix_ks ...
1399!> \param matrix_s ...
1400!> \param ot_preconditioner ...
1401!> \param prec_type ...
1402!> \param solver_type ...
1403!> \param energy_gap ...
1404!> \param nspins ...
1405!> \param has_unit_metric ...
1406!> \param convert_to_dbcsr ...
1407!> \param chol_type ...
1408!> \param full_mo_set ...
1409!> \param chebyshev_degree ...
1410!> \param low_rank_base ...
1411!> \param fermi_low_rank_max_rank ...
1412!> \param lattice_fft ...
1413!> \param lattice_fft_local_cells ...
1414! **************************************************************************************************
1415 SUBROUTINE prepare_preconditioner(qs_env, mos, matrix_ks, matrix_s, &
1416 ot_preconditioner, prec_type, solver_type, &
1417 energy_gap, nspins, has_unit_metric, &
1418 convert_to_dbcsr, chol_type, full_mo_set, chebyshev_degree, &
1419 low_rank_base, fermi_low_rank_max_rank, lattice_fft, &
1420 lattice_fft_local_cells)
1421
1422 TYPE(qs_environment_type), POINTER :: qs_env
1423 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
1424 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
1425 TYPE(preconditioner_p_type), DIMENSION(:), POINTER :: ot_preconditioner
1426 INTEGER, INTENT(IN) :: prec_type, solver_type
1427 REAL(dp), INTENT(IN) :: energy_gap
1428 INTEGER, INTENT(IN) :: nspins
1429 LOGICAL, INTENT(IN), OPTIONAL :: has_unit_metric, convert_to_dbcsr
1430 INTEGER, INTENT(IN), OPTIONAL :: chol_type
1431 LOGICAL, INTENT(IN), OPTIONAL :: full_mo_set
1432 INTEGER, INTENT(IN), OPTIONAL :: chebyshev_degree, low_rank_base, &
1433 fermi_low_rank_max_rank, lattice_fft, &
1434 lattice_fft_local_cells
1435
1436 CHARACTER(LEN=*), PARAMETER :: routinen = 'prepare_preconditioner'
1437
1438 CHARACTER(LEN=default_string_length) :: msg
1439 INTEGER :: handle, icall, ispin, my_chebyshev_degree, my_fermi_low_rank_max_rank, &
1440 my_lattice_fft, my_lattice_fft_local_cells, my_low_rank_base, n_loops
1441 INTEGER, DIMENSION(5) :: nocc, norb
1442 INTEGER, DIMENSION(:), POINTER :: multiple_unit_cell
1443 LOGICAL :: do_co_rotate, do_kpoints, my_convert_to_dbcsr, my_full_mo_set, &
1444 my_has_unit_metric, use_covariant_full_all, use_mo_coeff_b
1445 TYPE(cp_blacs_env_type), POINTER :: blacs_env
1446 TYPE(cp_fm_type), POINTER :: mo_coeff
1447 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: kinetic
1448 TYPE(dbcsr_type), POINTER :: matrix_t, mo_coeff_b
1449 TYPE(dft_control_type), POINTER :: dft_control
1450 TYPE(mo_set_type), ALLOCATABLE, DIMENSION(:) :: canonical_mos
1451 TYPE(mp_para_env_type), POINTER :: para_env
1452 TYPE(section_vals_type), POINTER :: input, subsys_section
1453
1454 CALL timeset(routinen, handle)
1455 NULLIFY (matrix_t, mo_coeff_b, mo_coeff, kinetic, dft_control, para_env, blacs_env, &
1456 input, subsys_section, multiple_unit_cell)
1457 my_has_unit_metric = .false.
1458 IF (PRESENT(has_unit_metric)) my_has_unit_metric = has_unit_metric
1459 my_convert_to_dbcsr = .true.
1460 IF (PRESENT(convert_to_dbcsr)) my_convert_to_dbcsr = convert_to_dbcsr
1461 my_full_mo_set = .false.
1462 IF (PRESENT(full_mo_set)) my_full_mo_set = full_mo_set
1463 my_chebyshev_degree = 8
1464 IF (PRESENT(chebyshev_degree)) my_chebyshev_degree = chebyshev_degree
1465 my_low_rank_base = ot_low_rank_base_overlap
1466 IF (PRESENT(low_rank_base)) my_low_rank_base = low_rank_base
1467 my_fermi_low_rank_max_rank = 48
1468 IF (PRESENT(fermi_low_rank_max_rank)) my_fermi_low_rank_max_rank = fermi_low_rank_max_rank
1469 my_lattice_fft = ot_lattice_fft_off
1470 IF (PRESENT(lattice_fft)) my_lattice_fft = lattice_fft
1471 my_lattice_fft_local_cells = 0
1472 IF (PRESENT(lattice_fft_local_cells)) my_lattice_fft_local_cells = lattice_fft_local_cells
1473 IF (my_lattice_fft_local_cells < 0) cpabort("LATTICE_FFT_LOCAL_CELLS must not be negative")
1474 do_kpoints = .false.
1475 use_covariant_full_all = prec_type == ot_precond_full_all_covariant
1476
1477 CALL get_qs_env(qs_env, &
1478 dft_control=dft_control, &
1479 do_kpoints=do_kpoints, &
1480 input=input, &
1481 para_env=para_env, &
1482 blacs_env=blacs_env)
1483
1484 IF (my_low_rank_base == ot_low_rank_base_lattice_fft) THEN
1485 IF (prec_type /= ot_precond_fermi_low_rank) THEN
1486 cpabort("FERMI_LOW_RANK_BASE LATTICE_FFT requires PRECONDITIONER FERMI_LOW_RANK")
1487 END IF
1488 my_lattice_fft = ot_lattice_fft_on
1489 ELSE IF (my_low_rank_base /= ot_low_rank_base_overlap) THEN
1490 cpabort("Unknown FERMI_LOW_RANK_BASE setting")
1491 END IF
1492
1493 IF (my_lattice_fft /= ot_lattice_fft_off .AND. &
1494 my_lattice_fft /= ot_lattice_fft_on .AND. &
1495 my_lattice_fft /= ot_lattice_fft_auto) THEN
1496 cpabort("Unknown LATTICE_FFT setting")
1497 END IF
1498 IF (my_lattice_fft_local_cells > 0 .AND. my_lattice_fft /= ot_lattice_fft_on) THEN
1499 cpabort("LATTICE_FFT_LOCAL_CELLS requires LATTICE_FFT ON")
1500 END IF
1501 IF (my_lattice_fft /= ot_lattice_fft_off) THEN
1502 SELECT CASE (prec_type)
1505 CASE DEFAULT
1506 cpabort("LATTICE_FFT does not support the selected OT preconditioner")
1507 END SELECT
1508 IF (solver_type /= ot_precond_solver_default .AND. &
1509 solver_type /= ot_precond_solver_inv_chol) THEN
1510 cpabort("LATTICE_FFT requires PRECOND_SOLVER DEFAULT or INVERSE_CHOLESKY")
1511 END IF
1512 IF (my_has_unit_metric .OR. do_kpoints .OR. .NOT. ASSOCIATED(input)) THEN
1513 IF (my_lattice_fft == ot_lattice_fft_on) THEN
1514 cpabort("LATTICE_FFT ON requires a non-unit metric in a real Gamma-point calculation")
1515 END IF
1516 my_lattice_fft = ot_lattice_fft_off
1517 END IF
1518 END IF
1519 IF (my_lattice_fft /= ot_lattice_fft_off) THEN
1520 subsys_section => section_vals_get_subs_vals(input, "SUBSYS")
1521 IF (.NOT. ASSOCIATED(subsys_section)) THEN
1522 IF (my_lattice_fft == ot_lattice_fft_on) cpabort("LATTICE_FFT ON requires a SUBSYS section")
1523 my_lattice_fft = ot_lattice_fft_off
1524 END IF
1525 IF (my_lattice_fft /= ot_lattice_fft_off) THEN
1526 CALL section_vals_val_get(subsys_section, "CELL%MULTIPLE_UNIT_CELL", &
1527 i_vals=multiple_unit_cell)
1528 IF (.NOT. ASSOCIATED(multiple_unit_cell)) THEN
1529 IF (my_lattice_fft == ot_lattice_fft_on) THEN
1530 cpabort("LATTICE_FFT ON requires a valid CELL%MULTIPLE_UNIT_CELL")
1531 END IF
1532 my_lattice_fft = ot_lattice_fft_off
1533 ELSE
1534 IF (SIZE(multiple_unit_cell) /= 3 .OR. any(multiple_unit_cell <= 0)) THEN
1535 IF (my_lattice_fft == ot_lattice_fft_on) THEN
1536 cpabort("LATTICE_FFT ON requires a valid CELL%MULTIPLE_UNIT_CELL")
1537 END IF
1538 my_lattice_fft = ot_lattice_fft_off
1539 ELSE IF (product(multiple_unit_cell) <= 1 .AND. my_lattice_fft == ot_lattice_fft_on) THEN
1540 cpabort("LATTICE_FFT ON requires an explicitly replicated multiple unit cell")
1541 END IF
1542 IF (my_lattice_fft_local_cells >= product(multiple_unit_cell)) THEN
1543 cpabort("LATTICE_FFT_LOCAL_CELLS must be smaller than the replicated cell count")
1544 END IF
1545 END IF
1546 END IF
1547 END IF
1548
1549 IF (dft_control%qs_control%semi_empirical .OR. dft_control%qs_control%dftb .OR. &
1550 dft_control%qs_control%xtb) THEN
1551 IF (prec_type == ot_precond_full_kinetic) THEN
1552 msg = "Full_kinetic not available for semi-empirical methods"
1553 cpabort(trim(msg))
1554 END IF
1555 matrix_t => matrix_s(1)%matrix
1556 ELSE
1557 cpassert(.NOT. my_has_unit_metric)
1558 CALL get_qs_env(qs_env, kinetic=kinetic)
1559 matrix_t => kinetic(1)%matrix
1560 END IF
1561
1562 ! use full set of MOs or just occupied MOs
1563 nocc = 0
1564 norb = 0
1565 IF (my_full_mo_set) THEN
1566 DO ispin = 1, nspins
1567 CALL get_mo_set(mo_set=mos(ispin), homo=nocc(ispin), nmo=norb(ispin))
1568 CALL set_mo_set(mo_set=mos(ispin), homo=norb(ispin))
1569 END DO
1570 END IF
1571 !determines how often make preconditioner is called, spin dependent methods have to be called twice
1572 n_loops = 1
1573 IF (prec_type == ot_precond_full_single_inverse) n_loops = nspins
1574 ! check whether we need the ev and rotate the MOs
1575 SELECT CASE (prec_type)
1576 CASE (ot_precond_full_all)
1577 ! FULL_ALL needs canonical occupied eigenstates and one preconditioner per spin.
1578 n_loops = nspins
1579
1580 do_co_rotate = ASSOCIATED(qs_env%mo_derivs)
1581 DO ispin = 1, nspins
1582 CALL get_mo_set(mo_set=mos(ispin), mo_coeff_b=mo_coeff_b, mo_coeff=mo_coeff)
1583 use_mo_coeff_b = mos(ispin)%use_mo_coeff_b
1584 IF (use_mo_coeff_b .AND. do_co_rotate) THEN
1585 CALL calculate_subspace_eigenvalues(mo_coeff_b, matrix_ks(ispin)%matrix, &
1586 do_rotation=.true., &
1587 co_rotate=qs_env%mo_derivs(ispin)%matrix, &
1588 para_env=para_env, &
1589 blacs_env=blacs_env)
1590 ELSE IF (use_mo_coeff_b) THEN
1591 CALL calculate_subspace_eigenvalues(mo_coeff_b, matrix_ks(ispin)%matrix, &
1592 do_rotation=.true., &
1593 para_env=para_env, &
1594 blacs_env=blacs_env)
1595 ELSE
1596 CALL calculate_subspace_eigenvalues(mo_coeff, matrix_ks(ispin)%matrix, &
1597 do_rotation=.true.)
1598 END IF
1599 END DO
1601 ! The common-reference low-rank operator depends only on the occupied projector.
1602 ! Keep the input gauge unchanged so that OT ROTATION remains well defined.
1603 n_loops = nspins
1605 ! Diagonalize a private occupied copy. The physical MOs and their
1606 ! derivatives remain in the user-selected gauge.
1607 n_loops = nspins
1608 ALLOCATE (canonical_mos(nspins))
1609 DO ispin = 1, nspins
1610 CALL duplicate_mo_set(canonical_mos(ispin), mos(ispin))
1611 CALL get_mo_set(mo_set=canonical_mos(ispin), mo_coeff_b=mo_coeff_b, mo_coeff=mo_coeff)
1612 IF (canonical_mos(ispin)%use_mo_coeff_b) THEN
1613 CALL copy_dbcsr_to_fm(mo_coeff_b, mo_coeff)
1614 canonical_mos(ispin)%use_mo_coeff_b = .false.
1615 END IF
1616 CALL canonicalize_full_all_copy(ot_preconditioner(ispin)%preconditioner, &
1617 mo_coeff, matrix_ks(ispin)%matrix, &
1618 canonical_mos(ispin)%homo)
1619 END DO
1620 CASE DEFAULT
1621 ! No need to rotate the MOs
1622 END SELECT
1623
1624 ! check whether we have a preconditioner
1625 SELECT CASE (prec_type)
1626 CASE (ot_precond_none)
1627 DO ispin = 1, SIZE(ot_preconditioner)
1628 ot_preconditioner(ispin)%preconditioner%in_use = 0
1629 END DO
1630 CASE DEFAULT
1631 DO icall = 1, n_loops
1632 ot_preconditioner(icall)%preconditioner%polynomial_degree = my_chebyshev_degree
1633 ot_preconditioner(icall)%preconditioner%low_rank_base = my_low_rank_base
1634 ot_preconditioner(icall)%preconditioner%lattice_fft_mode = my_lattice_fft
1635 ot_preconditioner(icall)%preconditioner%lattice_local_cells = my_lattice_fft_local_cells
1636 ot_preconditioner(icall)%preconditioner%lattice_dims = [1, 1, 1]
1637 IF (my_lattice_fft /= ot_lattice_fft_off) THEN
1638 cpassert(ASSOCIATED(multiple_unit_cell))
1639 ot_preconditioner(icall)%preconditioner%lattice_dims = multiple_unit_cell
1640 END IF
1641 ot_preconditioner(icall)%preconditioner%spectral_max_rank = my_fermi_low_rank_max_rank
1642 IF (my_has_unit_metric) THEN
1643 IF (use_covariant_full_all) THEN
1644 CALL make_preconditioner(ot_preconditioner(icall)%preconditioner, &
1645 prec_type, solver_type, &
1646 matrix_h=matrix_ks(icall)%matrix, &
1647 mo_set=canonical_mos(icall), &
1648 energy_gap=energy_gap, &
1649 convert_precond_to_dbcsr=my_convert_to_dbcsr)
1650 ELSE
1651 CALL make_preconditioner(ot_preconditioner(icall)%preconditioner, &
1652 prec_type, solver_type, &
1653 matrix_h=matrix_ks(icall)%matrix, &
1654 mo_set=mos(icall), &
1655 energy_gap=energy_gap, &
1656 convert_precond_to_dbcsr=my_convert_to_dbcsr)
1657 END IF
1658 ELSE
1659 IF (use_covariant_full_all) THEN
1660 CALL make_preconditioner(ot_preconditioner(icall)%preconditioner, &
1661 prec_type, solver_type, &
1662 matrix_h=matrix_ks(icall)%matrix, &
1663 matrix_s=matrix_s(1)%matrix, matrix_t=matrix_t, &
1664 mo_set=canonical_mos(icall), &
1665 energy_gap=energy_gap, &
1666 convert_precond_to_dbcsr=my_convert_to_dbcsr, chol_type=chol_type)
1667 ELSE
1668 CALL make_preconditioner(ot_preconditioner(icall)%preconditioner, &
1669 prec_type, solver_type, &
1670 matrix_h=matrix_ks(icall)%matrix, &
1671 matrix_s=matrix_s(1)%matrix, matrix_t=matrix_t, &
1672 mo_set=mos(icall), energy_gap=energy_gap, &
1673 convert_precond_to_dbcsr=my_convert_to_dbcsr, chol_type=chol_type)
1674 END IF
1675 END IF
1676 END DO
1677 END SELECT
1678
1679 IF (ALLOCATED(canonical_mos)) THEN
1680 DO ispin = 1, SIZE(canonical_mos)
1681 CALL deallocate_mo_set(canonical_mos(ispin))
1682 END DO
1683 DEALLOCATE (canonical_mos)
1684 END IF
1685
1686 ! reset homo values
1687 IF (my_full_mo_set) THEN
1688 DO ispin = 1, nspins
1689 CALL set_mo_set(mo_set=mos(ispin), homo=nocc(ispin))
1690 END DO
1691 END IF
1692
1693 CALL timestop(handle)
1694
1695 END SUBROUTINE prepare_preconditioner
1696
1697! **************************************************************************************************
1698!> \brief Canonicalize an occupied copy and retain the occupied-space rotation.
1699!> \param preconditioner_env ...
1700!> \param mo_coeff ...
1701!> \param matrix_h ...
1702!> \param nocc ...
1703! **************************************************************************************************
1704 SUBROUTINE canonicalize_full_all_copy(preconditioner_env, mo_coeff, matrix_h, nocc)
1705
1706 TYPE(preconditioner_type) :: preconditioner_env
1707 TYPE(cp_fm_type), POINTER :: mo_coeff
1708 TYPE(dbcsr_type), POINTER :: matrix_h
1709 INTEGER, INTENT(IN) :: nocc
1710
1711 CHARACTER(LEN=*), PARAMETER :: routinen = 'canonicalize_full_all_copy'
1712
1713 INTEGER :: handle, nao
1714 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: occupied_evals
1715 TYPE(cp_fm_struct_type), POINTER :: fm_struct
1716 TYPE(cp_fm_type) :: canonical_coeff, h_block, weighted_coeff
1717
1718 CALL timeset(routinen, handle)
1719 CALL cp_fm_get_info(mo_coeff, nrow_global=nao)
1720 cpassert(nocc > 0)
1721
1722 CALL cp_fm_struct_create(fm_struct, nrow_global=nao, ncol_global=nocc, &
1723 context=preconditioner_env%ctxt, &
1724 para_env=preconditioner_env%para_env)
1725 CALL cp_fm_create(canonical_coeff, fm_struct, name=routinen)
1726 CALL cp_fm_create(weighted_coeff, fm_struct, name=routinen)
1727 CALL cp_fm_struct_release(fm_struct)
1728 CALL cp_fm_to_fm(mo_coeff, canonical_coeff, nocc)
1729
1730 CALL cp_fm_struct_create(fm_struct, nrow_global=nocc, ncol_global=nocc, &
1731 context=preconditioner_env%ctxt, &
1732 para_env=preconditioner_env%para_env)
1733 CALL cp_fm_create(h_block, fm_struct, name=routinen)
1734 IF (ASSOCIATED(preconditioner_env%occ_rotation_dbcsr)) THEN
1735 CALL dbcsr_release_p(preconditioner_env%occ_rotation_dbcsr)
1736 END IF
1737 IF (ASSOCIATED(preconditioner_env%occ_rotation)) THEN
1738 CALL cp_fm_release(preconditioner_env%occ_rotation)
1739 DEALLOCATE (preconditioner_env%occ_rotation)
1740 END IF
1741 ALLOCATE (preconditioner_env%occ_rotation)
1742 CALL cp_fm_create(preconditioner_env%occ_rotation, fm_struct, name="occupied rotation")
1743 CALL cp_fm_struct_release(fm_struct)
1744
1745 CALL cp_dbcsr_sm_fm_multiply(matrix_h, canonical_coeff, weighted_coeff, nocc)
1746 CALL parallel_gemm('T', 'N', nocc, nocc, nao, 1.0_dp, canonical_coeff, &
1747 weighted_coeff, 0.0_dp, h_block)
1748 ALLOCATE (occupied_evals(nocc))
1749 CALL choose_eigv_solver(h_block, preconditioner_env%occ_rotation, occupied_evals)
1750
1751 CALL parallel_gemm('N', 'N', nao, nocc, nocc, 1.0_dp, canonical_coeff, &
1752 preconditioner_env%occ_rotation, 0.0_dp, weighted_coeff)
1753 CALL cp_fm_to_fm(weighted_coeff, mo_coeff, nocc, 1, 1)
1754
1755 DEALLOCATE (occupied_evals)
1756 CALL cp_fm_release(canonical_coeff)
1757 CALL cp_fm_release(h_block)
1758 CALL cp_fm_release(weighted_coeff)
1759 CALL timestop(handle)
1760
1761 END SUBROUTINE canonicalize_full_all_copy
1762
1763END MODULE preconditioner
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
methods related to the blacs parallel environment
Represents a complex full matrix distributed on many processors.
subroutine, public cp_cfm_release(matrix)
Releases a full matrix.
subroutine, public cp_fm_to_cfm(msourcer, msourcei, mtarget)
Construct a complex full matrix by taking its real and imaginary parts from two separate real-value f...
subroutine, public cp_cfm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
Creates a new full matrix with the given structure.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_release_p(matrix)
...
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)
...
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
subroutine, public copy_dbcsr_to_fm(matrix, fm, plan, check_plan)
Copy a DBCSR matrix to a BLACS matrix.
Basic linear algebra operations for full matrices.
subroutine, public cp_fm_uplo_to_full(matrix, work, uplo)
given a triangular matrix according to uplo, computes the corresponding full matrix
various cholesky decomposition related routines
subroutine, public cp_fm_cholesky_invert(matrix, n, info_out, work)
Replace the upper Cholesky factor by the upper triangle of the inverse. When work is provided,...
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,...
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
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:262
represent the structure of a full matrix
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
subroutine, public cp_fm_struct_release(fmstruct)
releases a full matrix structure
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
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
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
subroutine, public cp_fm_create(matrix, matrix_struct, name, nrow, ncol, set_zero)
creates a new full matrix with the given structure
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public ot_low_rank_base_overlap
integer, parameter, public ot_precond_full_all_covariant
integer, parameter, public ot_low_rank_base_lattice_fft
integer, parameter, public ot_lattice_fft_on
integer, parameter, public ot_precond_full_kinetic
integer, parameter, public cholesky_reduce
integer, parameter, public ot_lattice_fft_auto
integer, parameter, public ot_precond_solver_default
integer, parameter, public ot_lattice_fft_off
integer, parameter, public ot_precond_full_single
integer, parameter, public ot_precond_solver_inv_chol
integer, parameter, public ot_precond_none
integer, parameter, public ot_precond_full_single_inverse
integer, parameter, public ot_precond_fermi_low_rank
integer, parameter, public ot_precond_s_inverse
integer, parameter, public ot_precond_solver_update
integer, parameter, public ot_precond_full_all
objects that represent the structure of input sections and the data contained in an input section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Block-circulant inverse operators for replicated-cell OT preconditioners.
real(kind=dp), parameter, public lattice_fft_auto_max_cost_ratio
subroutine, public build_lattice_local_correction(operator_matrix, lattice_dims, local_cell_count, correction_indices, operator_columns, coarse_inverse, corrected_projection_error, info)
Build a balanced local correction for selected non-circulant cell blocks.
subroutine, public build_lattice_inverse(operator_matrix, lattice_dims, inverse_k, projection_error, used_fft, info)
Project a cell-major SPD operator onto lattice translations and invert its Fourier blocks.
pure real(kind=dp) function, public lattice_local_storage_ratio(lattice_dims, operator_count, local_cell_count)
Estimate balanced-local-correction storage relative to one dense real operator.
pure real(kind=dp) function, public lattice_fft_storage_ratio(lattice_dims, operator_count)
Estimate reciprocal-block storage relative to one dense real operator.
pure logical function, public lattice_fft_auto_select(lattice_dims, block_size, rhs_count, operator_count, num_pe, projection_error, used_fft, info)
Conservative error-, work-, and storage-aware AUTO selector.
subroutine, public lattice_grid(lattice_dims, index_to_cell, xkp)
Enumerate a regular direct/reciprocal lattice with the first direction varying fastest.
integer, parameter, public lattice_fft_auto_min_cells
subroutine, public build_lattice_inverse_from_blocks(blocks, lattice_dims, inverse_k, used_fft, info)
Transform and invert already projected real-space lattice blocks.
pure real(kind=dp) function, public lattice_local_cost_ratio(lattice_dims, local_cell_count)
Estimate balanced-local-correction work relative to a dense inverse application.
pure real(kind=dp) function, public lattice_fft_cost_ratio(lattice_dims, block_size, rhs_count, num_pe)
Estimate lattice-FFT application work relative to a dense inverse application.
real(kind=dp), parameter, public lattice_fft_auto_max_storage_ratio
Small algebraic helpers for the rotationally covariant low-rank OT preconditioner.
pure elemental real(kind=dp) function, public low_rank_inverse_weight(eigenvalue, reference, energy_gap)
Spectral inverse weight relative to a common occupied reference level.
pure integer function, public low_rank_select_rank(eigenvalues, nocc, max_rank, degeneracy_tolerance)
Select a bounded spectral rank without cutting a degenerate boundary manifold.
Interface to the message passing library MPI.
basic linear algebra operations for full matrixes
computes preconditioners, and implements methods to apply them currently used in qs_ot
subroutine, public apply_preconditioner_dbcsr_complex(preconditioner_env, matrix_in_re, matrix_in_im, matrix_out_re, matrix_out_im)
Apply a complex orbital preconditioner to paired real/imaginary DBCSR matrices.
subroutine, public apply_preconditioner_fm(preconditioner_env, matrix_in, matrix_out)
applies a previously created preconditioner to a full matrix
subroutine, public apply_preconditioner_cfm_complex(preconditioner_env, matrix_in, matrix_out)
Apply a complex k-point orbital preconditioner.
subroutine, public apply_preconditioner_dbcsr(preconditioner_env, matrix_in, matrix_out)
...
computes preconditioners, and implements methods to apply them currently used in qs_ot
subroutine, public make_complex_full_single_inverse(preconditioner_env, matrix_c0, matrix_h, matrix_s, energy_gap)
Build a gauge-covariant FULL_SINGLE_INVERSE operator for a complex k-point channel.
subroutine, public make_complex_full_kinetic(preconditioner_env, matrix_t, matrix_s, energy_gap)
Build the inverse complex kinetic-plus-overlap preconditioner.
subroutine, public make_complex_full_s_inverse(preconditioner_env, matrix_s)
Build the inverse complex overlap preconditioner.
subroutine, public make_complex_full_all_covariant(preconditioner_env, matrix_c0, matrix_h, matrix_s, energy_gap, occupation_signature)
Build the gauge-covariant FULL_ALL inverse for a complex k-point channel.
subroutine, public make_complex_full_single(preconditioner_env, matrix_h, matrix_s, energy_homo, energy_gap)
Build the complex spectral FULL_SINGLE preconditioner.
subroutine, public make_complex_fermi_low_rank(preconditioner_env, matrix_c0, matrix_h, matrix_s, energy_gap, max_rank, spectral_window, degeneracy_tolerance)
Build the bounded low-rank spectral inverse for a complex k-point channel.
subroutine, public make_preconditioner_matrix(preconditioner_env, matrix_h, matrix_s, matrix_t, mo_coeff, energy_homo, eigenvalues_ot, energy_gap, my_solver_type)
...
subroutine, public make_complex_full_all(preconditioner_env, matrix_c0, matrix_h, matrix_s, c0_evals, energy_gap)
Build the state-selective FULL_ALL operator for a complex k-point channel. The occupied/reference sub...
solves the preconditioner, contains to utility function for fm<->dbcsr transfers, should be moved soo...
subroutine, public transfer_dbcsr_to_fm(dbcsr_matrix, fm_matrix, para_env, context)
transfers a dbcsr to a full matrix
subroutine, public solve_preconditioner(my_solver_type, preconditioner_env, matrix_s, matrix_h)
...
subroutine, public transfer_fm_to_dbcsr(fm_matrix, dbcsr_matrix, template_mat)
transfers a full matrix to a dbcsr
types of preconditioners
subroutine, public init_preconditioner(preconditioner_env, para_env, blacs_env)
...
subroutine, public destroy_preconditioner(preconditioner_env)
...
computes preconditioners, and implements methods to apply them currently used in qs_ot
subroutine, public make_preconditioner_complex_full_all(preconditioner_env, matrix_c_re, matrix_c_im, matrix_h_re, matrix_h_im, matrix_s_re, matrix_s_im, mo_set, energy_gap, solver_type)
Construct FULL_ALL directly from one complex H(k), S(k), and C(k) channel.
subroutine, public make_preconditioner_complex_full_s_inverse(preconditioner_env, matrix_s_re, matrix_s_im, solver_type)
Construct a complex FULL_S_INVERSE preconditioner.
subroutine, public make_preconditioner_complex_full_single(preconditioner_env, matrix_h_re, matrix_h_im, matrix_s_re, matrix_s_im, mo_set, energy_gap, solver_type)
Construct a complex FULL_SINGLE preconditioner from H(k) and S(k).
subroutine, public dbcsr_pair_to_cfm(matrix_re, matrix_im, preconditioner_env, matrix_name, matrix)
Copy a real/imaginary DBCSR pair to one distributed complex full matrix.
subroutine, public make_preconditioner_complex_full_single_inverse(preconditioner_env, matrix_c_re, matrix_c_im, matrix_h_re, matrix_h_im, matrix_s_re, matrix_s_im, energy_gap, solver_type)
Construct a complex FULL_SINGLE_INVERSE preconditioner without discarding Im(H,S,C).
subroutine, public make_preconditioner_complex_full_all_covariant(preconditioner_env, matrix_c_re, matrix_c_im, matrix_h_re, matrix_h_im, matrix_s_re, matrix_s_im, energy_gap, solver_type, occupation_signature)
Construct rotation-covariant FULL_ALL for one complex k-point channel.
subroutine, public restart_preconditioner(qs_env, preconditioner, prec_type, nspins)
Allows for a restart of the preconditioner depending on the method it purges all arrays or keeps them...
subroutine, public prepare_preconditioner(qs_env, mos, matrix_ks, matrix_s, ot_preconditioner, prec_type, solver_type, energy_gap, nspins, has_unit_metric, convert_to_dbcsr, chol_type, full_mo_set, chebyshev_degree, low_rank_base, fermi_low_rank_max_rank, lattice_fft, lattice_fft_local_cells)
...
subroutine, public make_preconditioner_complex_fermi_low_rank(preconditioner_env, matrix_c_re, matrix_c_im, matrix_h_re, matrix_h_im, matrix_s_re, matrix_s_im, energy_gap, max_rank, solver_type)
Construct FERMI_LOW_RANK directly for one complex k-point channel.
subroutine, public make_preconditioner_complex_full_kinetic(preconditioner_env, matrix_t_re, matrix_t_im, matrix_s_re, matrix_s_im, energy_gap, solver_type)
Construct a complex FULL_KINETIC preconditioner.
subroutine, public make_preconditioner(preconditioner_env, precon_type, solver_type, matrix_h, matrix_s, matrix_t, mo_set, energy_gap, convert_precond_to_dbcsr, chol_type)
...
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, mimic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_vhxc, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
collects routines that perform operations directly related to MOs
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
subroutine, public duplicate_mo_set(mo_set_new, mo_set_old)
allocate a new mo_set, and copy the old data
subroutine, public set_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, uniform_occupation, kts, mu, flexible_electron_count)
Set the components of a MO set data structure.
subroutine, public deallocate_mo_set(mo_set)
Deallocate a wavefunction data structure.
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count, cmo_coeff)
Get the components of a MO set data structure.
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
Represent a complex full matrix.
keeps the information about the structure of a full matrix
represent a full matrix
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment