(git:1f9fd2c)
Loading...
Searching...
No Matches
almo_scf_qs.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 Interface between ALMO SCF and QS
10!> \par History
11!> 2011.05 created [Rustam Z Khaliullin]
12!> \author Rustam Z Khaliullin
13! **************************************************************************************************
22 USE cell_types, ONLY: cell_type,&
23 pbc
25 USE cp_dbcsr_api, ONLY: &
30 dbcsr_put_block, dbcsr_release, dbcsr_set, dbcsr_type, dbcsr_type_no_symmetry, &
37 USE cp_fm_types, ONLY: cp_fm_create,&
43 USE cp_units, ONLY: cp_unit_to_cp2k
52 USE kinds, ONLY: dp
62 USE qs_ks_types, ONLY: qs_ks_did_change,&
65 USE qs_mo_types, ONLY: allocate_mo_set,&
76 USE qs_rho_types, ONLY: qs_rho_get,&
78 USE qs_scf_types, ONLY: qs_scf_env_type,&
80#include "./base/base_uses.f90"
81
82 IMPLICIT NONE
83
84 PRIVATE
85
86 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'almo_scf_qs'
87
88 PUBLIC :: matrix_almo_create, &
97
98CONTAINS
99
100! **************************************************************************************************
101!> \brief create the ALMO matrix templates
102!> \param matrix_new ...
103!> \param matrix_qs ...
104!> \param almo_scf_env ...
105!> \param name_new ...
106!> \param size_keys ...
107!> \param symmetry_new ...
108!> \param spin_key ...
109!> \param init_domains ...
110!> \par History
111!> 2011.05 created [Rustam Z Khaliullin]
112!> \author Rustam Z Khaliullin
113! **************************************************************************************************
114 SUBROUTINE matrix_almo_create(matrix_new, matrix_qs, almo_scf_env, &
115 name_new, size_keys, symmetry_new, &
116 spin_key, init_domains)
117
118 TYPE(dbcsr_type) :: matrix_new, matrix_qs
119 TYPE(almo_scf_env_type), INTENT(IN) :: almo_scf_env
120 CHARACTER(len=*), INTENT(IN) :: name_new
121 INTEGER, DIMENSION(2), INTENT(IN) :: size_keys
122 CHARACTER, INTENT(IN) :: symmetry_new
123 INTEGER, INTENT(IN) :: spin_key
124 LOGICAL, INTENT(IN) :: init_domains
125
126 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_almo_create'
127
128 INTEGER :: dimen, handle, hold, iatom, iblock_col, &
129 iblock_row, imol, mynode, natoms, &
130 nblkrows_tot, nlength, nmols, row
131 INTEGER, DIMENSION(:), POINTER :: blk_distr, blk_sizes, block_sizes_new, col_blk_size, &
132 col_distr_new, col_sizes_new, distr_new_array, row_blk_size, row_distr_new, row_sizes_new
133 LOGICAL :: active, one_dim_is_mo, tr
134 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: new_block
135 TYPE(dbcsr_distribution_type) :: dist_new, dist_qs
136
137! dimension size: AO, MO, etc
138! almo_mat_dim_aobasis - no. of AOs,
139! almo_mat_dim_occ - no. of occupied MOs
140! almo_mat_dim_domains - no. of domains
141! symmetry type: dbcsr_type_no_symmetry, dbcsr_type_symmetric,
142! dbcsr_type_antisymmetric, dbcsr_type_hermitian, dbcsr_type_antihermitian
143! (see dbcsr_lib/dbcsr_types.F for other values)
144! spin_key: either 1 or 2 (0 is allowed for matrics in the AO basis)
145! TYPE(dbcsr_iterator_type) :: iter
146! REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: allones
147!-----------------------------------------------------------------------
148
149 CALL timeset(routinen, handle)
150
151 ! RZK-warning The structure of the matrices can be optimized:
152 ! 1. Diagonal matrices must be distributed evenly over the processes.
153 ! This can be achieved by distributing cpus: 012012-rows and 001122-cols
154 ! block_diagonal_flag is introduced but not used
155 ! 2. Multiplication of diagonally dominant matrices will be faster
156 ! if the diagonal blocks are local to the same processes.
157 ! 3. Systems of molecules of drastically different sizes might need
158 ! better distribution.
159
160 ! obtain distribution from the qs matrix - it might be useful
161 ! to get the structure of the AO dimensions
162 CALL dbcsr_get_info(matrix_qs, distribution=dist_qs)
163
164 natoms = almo_scf_env%natoms
165 nmols = almo_scf_env%nmolecules
166
167 DO dimen = 1, 2 ! 1 - row, 2 - column dimension
168
169 ! distribution pattern is the same for all matrix types (ao, occ, virt)
170 IF (dimen == 1) THEN !rows
171 CALL dbcsr_distribution_get(dist_qs, row_dist=blk_distr)
172 ELSE !columns
173 CALL dbcsr_distribution_get(dist_qs, col_dist=blk_distr)
174 END IF
175
176 IF (size_keys(dimen) == almo_mat_dim_aobasis) THEN ! this dimension is AO
177
178 ! structure of an AO dimension can be copied from matrix_qs
179 CALL dbcsr_get_info(matrix_qs, row_blk_size=blk_sizes)
180
181 ! atomic clustering of AOs
182 IF (almo_scf_env%mat_distr_aos == almo_mat_distr_atomic) THEN
183 ALLOCATE (block_sizes_new(natoms), distr_new_array(natoms))
184 block_sizes_new(:) = blk_sizes(:)
185 distr_new_array(:) = blk_distr(:)
186 ! molecular clustering of AOs
187 ELSE IF (almo_scf_env%mat_distr_aos == almo_mat_distr_molecular) THEN
188 ALLOCATE (block_sizes_new(nmols), distr_new_array(nmols))
189 block_sizes_new(:) = 0
190 DO iatom = 1, natoms
191 block_sizes_new(almo_scf_env%domain_index_of_atom(iatom)) = &
192 block_sizes_new(almo_scf_env%domain_index_of_atom(iatom)) + &
193 blk_sizes(iatom)
194 END DO
195 DO imol = 1, nmols
196 distr_new_array(imol) = &
197 blk_distr(almo_scf_env%first_atom_of_domain(imol))
198 END DO
199 ELSE
200 cpabort("Illegal distribution")
201 END IF
202
203 ELSE ! this dimension is not AO
204
205 IF (size_keys(dimen) == almo_mat_dim_occ .OR. &
206 size_keys(dimen) == almo_mat_dim_virt .OR. &
207 size_keys(dimen) == almo_mat_dim_virt_disc .OR. &
208 size_keys(dimen) == almo_mat_dim_virt_full) THEN ! this dim is MO
209
210 ! atomic clustering of MOs
211 IF (almo_scf_env%mat_distr_mos == almo_mat_distr_atomic) THEN
212 nlength = natoms
213 ALLOCATE (block_sizes_new(nlength))
214 block_sizes_new(:) = 0
215 IF (size_keys(dimen) == almo_mat_dim_occ) THEN
216 ! currently distributing atomic distr of mos is not allowed
217 ! RZK-warning define nocc_of_atom and nvirt_atom to implement it
218 !block_sizes_new(:)=almo_scf_env%nocc_of_atom(:,spin_key)
219 ELSE IF (size_keys(dimen) == almo_mat_dim_virt) THEN
220 !block_sizes_new(:)=almo_scf_env%nvirt_of_atom(:,spin_key)
221 END IF
222 ! molecular clustering of MOs
223 ELSE IF (almo_scf_env%mat_distr_mos == almo_mat_distr_molecular) THEN
224 nlength = nmols
225 ALLOCATE (block_sizes_new(nlength))
226 IF (size_keys(dimen) == almo_mat_dim_occ) THEN
227 block_sizes_new(:) = almo_scf_env%nocc_of_domain(:, spin_key)
228 ! Handle zero-electron fragments by adding one-orbital that
229 ! must remain zero at all times
230 WHERE (block_sizes_new == 0) block_sizes_new = 1
231 ELSE IF (size_keys(dimen) == almo_mat_dim_virt_disc) THEN
232 block_sizes_new(:) = almo_scf_env%nvirt_disc_of_domain(:, spin_key)
233 ELSE IF (size_keys(dimen) == almo_mat_dim_virt_full) THEN
234 block_sizes_new(:) = almo_scf_env%nvirt_full_of_domain(:, spin_key)
235 ELSE IF (size_keys(dimen) == almo_mat_dim_virt) THEN
236 block_sizes_new(:) = almo_scf_env%nvirt_of_domain(:, spin_key)
237 END IF
238 ELSE
239 cpabort("Illegal distribution")
240 END IF
241
242 ELSE
243
244 cpabort("Illegal dimension")
245
246 END IF ! end choosing dim size (occ, virt)
247
248 ! distribution for MOs is copied from AOs
249 ALLOCATE (distr_new_array(nlength))
250 ! atomic clustering
251 IF (almo_scf_env%mat_distr_mos == almo_mat_distr_atomic) THEN
252 distr_new_array(:) = blk_distr(:)
253 ! molecular clustering
254 ELSE IF (almo_scf_env%mat_distr_mos == almo_mat_distr_molecular) THEN
255 DO imol = 1, nmols
256 distr_new_array(imol) = &
257 blk_distr(almo_scf_env%first_atom_of_domain(imol))
258 END DO
259 END IF
260 END IF ! end choosing dimension size (AOs vs .NOT.AOs)
261
262 ! create final arrays
263 IF (dimen == 1) THEN !rows
264 row_sizes_new => block_sizes_new
265 row_distr_new => distr_new_array
266 ELSE !columns
267 col_sizes_new => block_sizes_new
268 col_distr_new => distr_new_array
269 END IF
270 END DO ! both rows and columns are done
271
272 ! Create the distribution
273 CALL dbcsr_distribution_new(dist_new, template=dist_qs, &
274 row_dist=row_distr_new, col_dist=col_distr_new, &
275 reuse_arrays=.true.)
276
277 ! Create the matrix
278 CALL dbcsr_create(matrix_new, name_new, &
279 dist_new, symmetry_new, &
280 row_sizes_new, col_sizes_new, reuse_arrays=.true.)
281 CALL dbcsr_distribution_release(dist_new)
282
283 ! fill out reqired blocks with 1.0_dp to tell the dbcsr library
284 ! which blocks to keep
285 IF (init_domains) THEN
286
287 CALL dbcsr_distribution_get(dist_new, mynode=mynode)
288 CALL dbcsr_work_create(matrix_new, work_mutable=.true.)
289 CALL dbcsr_get_info(matrix_new, nblkrows_total=nblkrows_tot, &
290 row_blk_size=row_blk_size, col_blk_size=col_blk_size)
291 ! startQQQ - this part of the code scales quadratically
292 ! therefore it is replaced with a less general but linear scaling algorithm below
293 ! the quadratic algorithm is kept to be re-written later
294 !QQQCALL dbcsr_get_info(matrix_new, nblkrows_total=nblkrows_tot, nblkcols_total=nblkcols_tot)
295 !QQQDO row = 1, nblkrows_tot
296 !QQQ DO col = 1, nblkcols_tot
297 !QQQ tr = .FALSE.
298 !QQQ iblock_row = row
299 !QQQ iblock_col = col
300 !QQQ CALL dbcsr_get_stored_coordinates(matrix_new, iblock_row, iblock_col, tr, hold)
301
302 !QQQ IF(hold==mynode) THEN
303 !QQQ
304 !QQQ ! RZK-warning replace with a function which says if this
305 !QQQ ! distribution block is active or not
306 !QQQ ! Translate indeces of distribution blocks to domain blocks
307 !QQQ if (size_keys(1)==almo_mat_dim_aobasis) then
308 !QQQ domain_row=almo_scf_env%domain_index_of_ao_block(iblock_row)
309 !QQQ else if (size_keys(2)==almo_mat_dim_occ .OR. &
310 !QQQ size_keys(2)==almo_mat_dim_virt .OR. &
311 !QQQ size_keys(2)==almo_mat_dim_virt_disc .OR. &
312 !QQQ size_keys(2)==almo_mat_dim_virt_full) then
313 !QQQ domain_row=almo_scf_env%domain_index_of_mo_block(iblock_row)
314 !QQQ else
315 !QQQ CPErrorMessage(cp_failure_level,routineP,"Illegal dimension")
316 !QQQ CPPrecondition(.FALSE.,cp_failure_level,routineP,failure)
317 !QQQ endif
318
319 !QQQ if (size_keys(2)==almo_mat_dim_aobasis) then
320 !QQQ domain_col=almo_scf_env%domain_index_of_ao_block(iblock_col)
321 !QQQ else if (size_keys(2)==almo_mat_dim_occ .OR. &
322 !QQQ size_keys(2)==almo_mat_dim_virt .OR. &
323 !QQQ size_keys(2)==almo_mat_dim_virt_disc .OR. &
324 !QQQ size_keys(2)==almo_mat_dim_virt_full) then
325 !QQQ domain_col=almo_scf_env%domain_index_of_mo_block(iblock_col)
326 !QQQ else
327 !QQQ CPErrorMessage(cp_failure_level,routineP,"Illegal dimension")
328 !QQQ CPPrecondition(.FALSE.,cp_failure_level,routineP,failure)
329 !QQQ endif
330
331 !QQQ ! Finds if we need this block
332 !QQQ ! only the block-diagonal constraint is implemented here
333 !QQQ active=.false.
334 !QQQ if (domain_row==domain_col) active=.true.
335
336 !QQQ IF (active) THEN
337 !QQQ ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
338 !QQQ new_block(:, :) = 1.0_dp
339 !QQQ CALL dbcsr_put_block(matrix_new, iblock_row, iblock_col, new_block)
340 !QQQ DEALLOCATE (new_block)
341 !QQQ ENDIF
342
343 !QQQ ENDIF ! mynode
344 !QQQ ENDDO
345 !QQQENDDO
346 !QQQtake care of zero-electron fragments
347 ! endQQQ - end of the quadratic part
348 ! start linear-scaling replacement:
349 ! works only for molecular blocks AND molecular distributions
350 DO row = 1, nblkrows_tot
351 tr = .false.
352 iblock_row = row
353 iblock_col = row
354 CALL dbcsr_get_stored_coordinates(matrix_new, iblock_row, iblock_col, hold)
355
356 IF (hold == mynode) THEN
357
358 active = .true.
359
360 one_dim_is_mo = .false.
361 DO dimen = 1, 2 ! 1 - row, 2 - column dimension
362 IF (size_keys(dimen) == almo_mat_dim_occ) one_dim_is_mo = .true.
363 END DO
364 IF (one_dim_is_mo) THEN
365 IF (almo_scf_env%nocc_of_domain(row, spin_key) == 0) active = .false.
366 END IF
367
368 one_dim_is_mo = .false.
369 DO dimen = 1, 2
370 IF (size_keys(dimen) == almo_mat_dim_virt) one_dim_is_mo = .true.
371 END DO
372 IF (one_dim_is_mo) THEN
373 IF (almo_scf_env%nvirt_of_domain(row, spin_key) == 0) active = .false.
374 END IF
375
376 one_dim_is_mo = .false.
377 DO dimen = 1, 2
378 IF (size_keys(dimen) == almo_mat_dim_virt_disc) one_dim_is_mo = .true.
379 END DO
380 IF (one_dim_is_mo) THEN
381 IF (almo_scf_env%nvirt_disc_of_domain(row, spin_key) == 0) active = .false.
382 END IF
383
384 one_dim_is_mo = .false.
385 DO dimen = 1, 2
386 IF (size_keys(dimen) == almo_mat_dim_virt_full) one_dim_is_mo = .true.
387 END DO
388 IF (one_dim_is_mo) THEN
389 IF (almo_scf_env%nvirt_full_of_domain(row, spin_key) == 0) active = .false.
390 END IF
391
392 IF (active) THEN
393 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
394 new_block(:, :) = 1.0_dp
395 CALL dbcsr_put_block(matrix_new, iblock_row, iblock_col, new_block)
396 DEALLOCATE (new_block)
397 END IF
398
399 END IF ! mynode
400 END DO
401 ! end lnear-scaling replacement
402
403 END IF ! init_domains
404
405 CALL dbcsr_finalize(matrix_new)
406
407 CALL timestop(handle)
408
409 END SUBROUTINE matrix_almo_create
410
411! **************************************************************************************************
412!> \brief convert between two types of matrices: QS style to ALMO style
413!> \param matrix_qs ...
414!> \param matrix_almo ...
415!> \param mat_distr_aos ...
416!> \par History
417!> 2011.06 created [Rustam Z Khaliullin]
418!> \author Rustam Z Khaliullin
419! **************************************************************************************************
420 SUBROUTINE matrix_qs_to_almo(matrix_qs, matrix_almo, mat_distr_aos)
421
422 TYPE(dbcsr_type) :: matrix_qs, matrix_almo
423 INTEGER :: mat_distr_aos
424
425 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_qs_to_almo'
426
427 INTEGER :: handle
428 TYPE(dbcsr_type) :: matrix_qs_nosym
429
430 CALL timeset(routinen, handle)
431 !RZK-warning if it's not a N(AO)xN(AO) matrix then stop
432
433 SELECT CASE (mat_distr_aos)
435 ! automatic data_type conversion
436 CALL dbcsr_copy(matrix_almo, matrix_qs)
438 ! desymmetrize the qs matrix
439 CALL dbcsr_create(matrix_qs_nosym, template=matrix_qs, matrix_type=dbcsr_type_no_symmetry)
440 CALL dbcsr_desymmetrize(matrix_qs, matrix_qs_nosym)
441
442 ! perform the magic complete_redistribute
443 ! before calling complete_redistribute set all blocks to zero
444 ! otherwise the non-zero elements of the redistributed matrix,
445 ! which are in zero-blocks of the original matrix, will remain
446 ! in the final redistributed matrix. this is a bug in
447 ! complete_redistribute. RZK-warning it should be later corrected by calling
448 ! dbcsr_set to 0.0 from within complete_redistribute
449 CALL dbcsr_set(matrix_almo, 0.0_dp)
450 CALL dbcsr_complete_redistribute(matrix_qs_nosym, matrix_almo)
451 CALL dbcsr_release(matrix_qs_nosym)
452
453 CASE DEFAULT
454 cpabort("")
455 END SELECT
456
457 CALL timestop(handle)
458
459 END SUBROUTINE matrix_qs_to_almo
460
461! **************************************************************************************************
462!> \brief convert between two types of matrices: ALMO style to QS style
463!> \param matrix_almo ...
464!> \param matrix_qs ...
465!> \param mat_distr_aos ...
466!> \par History
467!> 2011.06 created [Rustam Z Khaliullin]
468!> \author Rustam Z Khaliullin
469! **************************************************************************************************
470 SUBROUTINE matrix_almo_to_qs(matrix_almo, matrix_qs, mat_distr_aos)
471 TYPE(dbcsr_type) :: matrix_almo, matrix_qs
472 INTEGER, INTENT(IN) :: mat_distr_aos
473
474 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_almo_to_qs'
475
476 INTEGER :: handle
477 TYPE(dbcsr_type) :: matrix_almo_redist
478
479 CALL timeset(routinen, handle)
480 ! RZK-warning if it's not a N(AO)xN(AO) matrix then stop
481
482 SELECT CASE (mat_distr_aos)
484 CALL dbcsr_copy(matrix_qs, matrix_almo, keep_sparsity=.true.)
486 CALL dbcsr_create(matrix_almo_redist, template=matrix_qs)
487 CALL dbcsr_complete_redistribute(matrix_almo, matrix_almo_redist)
488 CALL dbcsr_set(matrix_qs, 0.0_dp)
489 CALL dbcsr_copy(matrix_qs, matrix_almo_redist, keep_sparsity=.true.)
490 CALL dbcsr_release(matrix_almo_redist)
491 CASE DEFAULT
492 cpabort("")
493 END SELECT
494
495 CALL timestop(handle)
496
497 END SUBROUTINE matrix_almo_to_qs
498
499! **************************************************************************************************
500!> \brief Initialization of the QS and ALMO KS matrix
501!> \param qs_env ...
502!> \param matrix_ks ...
503!> \param mat_distr_aos ...
504!> \param eps_filter ...
505!> \par History
506!> 2011.05 created [Rustam Z Khaliullin]
507!> \author Rustam Z Khaliullin
508! **************************************************************************************************
509 SUBROUTINE init_almo_ks_matrix_via_qs(qs_env, matrix_ks, mat_distr_aos, eps_filter)
510
511 TYPE(qs_environment_type), POINTER :: qs_env
512 TYPE(dbcsr_type), DIMENSION(:) :: matrix_ks
513 INTEGER :: mat_distr_aos
514 REAL(kind=dp) :: eps_filter
515
516 CHARACTER(len=*), PARAMETER :: routinen = 'init_almo_ks_matrix_via_qs'
517
518 INTEGER :: handle, ispin, nspin
519 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_qs_ks, matrix_qs_s
520 TYPE(dft_control_type), POINTER :: dft_control
521 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
522 POINTER :: sab_orb
523 TYPE(qs_ks_env_type), POINTER :: ks_env
524
525 CALL timeset(routinen, handle)
526
527 NULLIFY (sab_orb)
528
529 ! get basic quantities from the qs_env
530 CALL get_qs_env(qs_env, &
531 dft_control=dft_control, &
532 matrix_s=matrix_qs_s, &
533 matrix_ks=matrix_qs_ks, &
534 ks_env=ks_env, &
535 sab_orb=sab_orb)
536
537 nspin = dft_control%nspins
538
539 ! create matrix_ks in the QS env if necessary
540 IF (.NOT. ASSOCIATED(matrix_qs_ks)) THEN
541 CALL dbcsr_allocate_matrix_set(matrix_qs_ks, nspin)
542 DO ispin = 1, nspin
543 ALLOCATE (matrix_qs_ks(ispin)%matrix)
544 CALL dbcsr_create(matrix_qs_ks(ispin)%matrix, &
545 template=matrix_qs_s(1)%matrix)
546 CALL cp_dbcsr_alloc_block_from_nbl(matrix_qs_ks(ispin)%matrix, sab_orb)
547 CALL dbcsr_set(matrix_qs_ks(ispin)%matrix, 0.0_dp)
548 END DO
549 CALL set_ks_env(ks_env, matrix_ks=matrix_qs_ks)
550 END IF
551
552 ! copy to ALMO
553 DO ispin = 1, nspin
554 CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, matrix_ks(ispin), mat_distr_aos)
555 CALL dbcsr_filter(matrix_ks(ispin), eps_filter)
556 END DO
557
558 CALL timestop(handle)
559
560 END SUBROUTINE init_almo_ks_matrix_via_qs
561
562! **************************************************************************************************
563!> \brief Create MOs in the QS env to be able to return ALMOs to QS
564!> \param qs_env ...
565!> \param almo_scf_env ...
566!> \par History
567!> 2016.12 created [Yifei Shi]
568!> \author Yifei Shi
569! **************************************************************************************************
570 SUBROUTINE construct_qs_mos(qs_env, almo_scf_env)
571
572 TYPE(qs_environment_type), POINTER :: qs_env
573 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env
574
575 CHARACTER(len=*), PARAMETER :: routinen = 'construct_qs_mos'
576
577 INTEGER :: handle, ispin, ncol_fm, nrow_fm
578 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
579 TYPE(cp_fm_type) :: mo_fm_copy
580 TYPE(dft_control_type), POINTER :: dft_control
581 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
582 TYPE(qs_scf_env_type), POINTER :: scf_env
583
584 CALL timeset(routinen, handle)
585
586 ! create and init scf_env (this is necessary to return MOs to qs)
587 NULLIFY (mos, fm_struct_tmp, scf_env)
588 ALLOCATE (scf_env)
589 CALL scf_env_create(scf_env)
590
591 !CALL qs_scf_env_initialize(qs_env, scf_env)
592 CALL set_qs_env(qs_env, scf_env=scf_env)
593 CALL get_qs_env(qs_env, dft_control=dft_control, mos=mos)
594
595 CALL dbcsr_get_info(almo_scf_env%matrix_t(1), nfullrows_total=nrow_fm, nfullcols_total=ncol_fm)
596
597 ! allocate and init mo_set
598 DO ispin = 1, almo_scf_env%nspins
599 CALL dbcsr_get_info(almo_scf_env%matrix_t(ispin), nfullrows_total=nrow_fm, nfullcols_total=ncol_fm)
600
601 ! Currently only fm version of mo_set is usable.
602 ! First transform the matrix_t to fm version
603 ! Empty the containers to prevent memory leaks
604 CALL deallocate_mo_set(mos(ispin))
605
606 IF (almo_scf_env%nspins == 1) THEN
607 CALL allocate_mo_set(mo_set=mos(ispin), &
608 nao=nrow_fm, &
609 nmo=ncol_fm, &
610 nelectron=almo_scf_env%nelectrons_total, &
611 n_el_f=real(almo_scf_env%nelectrons_total, dp), &
612 maxocc=2.0_dp, &
613 flexible_electron_count=dft_control%relax_multiplicity)
614 ELSEIF (almo_scf_env%nspins == 2) THEN
615 CALL allocate_mo_set(mo_set=mos(ispin), &
616 nao=nrow_fm, &
617 nmo=ncol_fm, &
618 nelectron=sum(almo_scf_env%nocc_of_domain(:, ispin)), &
619 n_el_f=real(sum(almo_scf_env%nocc_of_domain(:, ispin)), dp), &
620 maxocc=1.0_dp, &
621 flexible_electron_count=dft_control%relax_multiplicity)
622 END IF
623
624 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=nrow_fm, ncol_global=ncol_fm, &
625 context=almo_scf_env%blacs_env, &
626 para_env=almo_scf_env%para_env)
627
628 CALL cp_fm_create(mo_fm_copy, fm_struct_tmp, name="t_orthogonal_converted_to_fm")
629 CALL cp_fm_struct_release(fm_struct_tmp)
630 !CALL copy_dbcsr_to_fm(almo_scf_env%matrix_t(ispin), mo_fm_copy)
631
632 CALL init_mo_set(mos(ispin), fm_ref=mo_fm_copy, name='fm_mo')
633
634 CALL cp_fm_release(mo_fm_copy)
635
636 END DO
637
638 CALL timestop(handle)
639
640 END SUBROUTINE construct_qs_mos
641
642! **************************************************************************************************
643!> \brief return density matrix to the qs_env
644!> \param qs_env ...
645!> \param matrix_p ...
646!> \param mat_distr_aos ...
647!> \par History
648!> 2011.05 created [Rustam Z Khaliullin]
649!> \author Rustam Z Khaliullin
650! **************************************************************************************************
651 SUBROUTINE almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
652 TYPE(qs_environment_type), POINTER :: qs_env
653 TYPE(dbcsr_type), DIMENSION(:) :: matrix_p
654 INTEGER, INTENT(IN) :: mat_distr_aos
655
656 CHARACTER(len=*), PARAMETER :: routinen = 'almo_dm_to_qs_env'
657
658 INTEGER :: handle, ispin, nspins
659 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
660 TYPE(qs_rho_type), POINTER :: rho
661
662 CALL timeset(routinen, handle)
663
664 NULLIFY (rho, rho_ao)
665 nspins = SIZE(matrix_p)
666 CALL get_qs_env(qs_env, rho=rho)
667 CALL qs_rho_get(rho, rho_ao=rho_ao)
668
669 ! set the new density matrix
670 DO ispin = 1, nspins
671 CALL matrix_almo_to_qs(matrix_p(ispin), &
672 rho_ao(ispin)%matrix, &
673 mat_distr_aos)
674 END DO
675 CALL qs_rho_update_rho(rho, qs_env=qs_env)
676 CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.true.)
677
678 CALL timestop(handle)
679
680 END SUBROUTINE almo_dm_to_qs_env
681
682! **************************************************************************************************
683!> \brief uses the ALMO density matrix
684!> to compute KS matrix (inside QS environment) and the new energy
685!> \param qs_env ...
686!> \param matrix_p ...
687!> \param energy_total ...
688!> \param mat_distr_aos ...
689!> \param smear ...
690!> \param kTS_sum ...
691!> \par History
692!> 2011.05 created [Rustam Z Khaliullin]
693!> 2018.09 smearing support [Ruben Staub]
694!> \author Rustam Z Khaliullin
695! **************************************************************************************************
696 SUBROUTINE almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos, smear, kTS_sum)
697 TYPE(qs_environment_type), POINTER :: qs_env
698 TYPE(dbcsr_type), DIMENSION(:) :: matrix_p
699 REAL(kind=dp) :: energy_total
700 INTEGER, INTENT(IN) :: mat_distr_aos
701 LOGICAL, INTENT(IN), OPTIONAL :: smear
702 REAL(kind=dp), INTENT(IN), OPTIONAL :: kts_sum
703
704 CHARACTER(len=*), PARAMETER :: routinen = 'almo_dm_to_qs_ks'
705
706 INTEGER :: handle
707 LOGICAL :: smearing
708 REAL(kind=dp) :: entropic_term
709 TYPE(qs_energy_type), POINTER :: energy
710
711 CALL timeset(routinen, handle)
712
713 IF (PRESENT(smear)) THEN
714 smearing = smear
715 ELSE
716 smearing = .false.
717 END IF
718
719 IF (PRESENT(kts_sum)) THEN
720 entropic_term = kts_sum
721 ELSE
722 entropic_term = 0.0_dp
723 END IF
724
725 NULLIFY (energy)
726 CALL get_qs_env(qs_env, energy=energy)
727 CALL almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
728 CALL qs_ks_update_qs_env(qs_env, calculate_forces=.false., just_energy=.false., &
729 print_active=.true.)
730
731 !! Add electronic entropy contribution if smearing is requested
732 !! Previous QS entropy is replaced by the sum of the entropy for each spin
733 IF (smearing) THEN
734 energy%total = energy%total - energy%kTS + entropic_term
735 END IF
736
737 energy_total = energy%total
738
739 CALL timestop(handle)
740
741 END SUBROUTINE almo_dm_to_qs_ks
742
743! **************************************************************************************************
744!> \brief uses the ALMO density matrix
745!> to compute ALMO KS matrix and the new energy
746!> \param qs_env ...
747!> \param matrix_p ...
748!> \param matrix_ks ...
749!> \param energy_total ...
750!> \param eps_filter ...
751!> \param mat_distr_aos ...
752!> \param smear ...
753!> \param kTS_sum ...
754!> \par History
755!> 2011.05 created [Rustam Z Khaliullin]
756!> 2018.09 smearing support [Ruben Staub]
757!> \author Rustam Z Khaliullin
758! **************************************************************************************************
759 SUBROUTINE almo_dm_to_almo_ks(qs_env, matrix_p, matrix_ks, energy_total, eps_filter, &
760 mat_distr_aos, smear, kTS_sum)
761
762 TYPE(qs_environment_type), POINTER :: qs_env
763 TYPE(dbcsr_type), DIMENSION(:) :: matrix_p, matrix_ks
764 REAL(kind=dp) :: energy_total, eps_filter
765 INTEGER, INTENT(IN) :: mat_distr_aos
766 LOGICAL, INTENT(IN), OPTIONAL :: smear
767 REAL(kind=dp), INTENT(IN), OPTIONAL :: kts_sum
768
769 CHARACTER(len=*), PARAMETER :: routinen = 'almo_dm_to_almo_ks'
770
771 INTEGER :: handle, ispin, nspins
772 LOGICAL :: smearing
773 REAL(kind=dp) :: entropic_term
774 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_qs_ks
775
776 CALL timeset(routinen, handle)
777
778 IF (PRESENT(smear)) THEN
779 smearing = smear
780 ELSE
781 smearing = .false.
782 END IF
783
784 IF (PRESENT(kts_sum)) THEN
785 entropic_term = kts_sum
786 ELSE
787 entropic_term = 0.0_dp
788 END IF
789
790 ! update KS matrix in the QS env
791 CALL almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos, &
792 smear=smearing, &
793 kts_sum=entropic_term)
794
795 nspins = SIZE(matrix_ks)
796
797 ! get KS matrix from the QS env and convert to the ALMO format
798 CALL get_qs_env(qs_env, matrix_ks=matrix_qs_ks)
799 DO ispin = 1, nspins
800 CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, matrix_ks(ispin), mat_distr_aos)
801 CALL dbcsr_filter(matrix_ks(ispin), eps_filter)
802 END DO
803
804 CALL timestop(handle)
805
806 END SUBROUTINE almo_dm_to_almo_ks
807
808! **************************************************************************************************
809!> \brief update qs_env total energy
810!> \param qs_env ...
811!> \param energy ...
812!> \param energy_singles_corr ...
813!> \par History
814!> 2013.03 created [Rustam Z Khaliullin]
815!> \author Rustam Z Khaliullin
816! **************************************************************************************************
817 SUBROUTINE almo_scf_update_ks_energy(qs_env, energy, energy_singles_corr)
818
819 TYPE(qs_environment_type), POINTER :: qs_env
820 REAL(kind=dp), INTENT(IN), OPTIONAL :: energy, energy_singles_corr
821
822 TYPE(qs_energy_type), POINTER :: qs_energy
823
824 CALL get_qs_env(qs_env, energy=qs_energy)
825
826 IF (PRESENT(energy_singles_corr)) THEN
827 qs_energy%singles_corr = energy_singles_corr
828 ELSE
829 qs_energy%singles_corr = 0.0_dp
830 END IF
831
832 IF (PRESENT(energy)) THEN
833 qs_energy%total = energy
834 END IF
835
836 qs_energy%total = qs_energy%total + qs_energy%singles_corr
837
838 END SUBROUTINE almo_scf_update_ks_energy
839
840! **************************************************************************************************
841!> \brief Creates the matrix that imposes absolute locality on MOs
842!> \param qs_env ...
843!> \param almo_scf_env ...
844!> \par History
845!> 2011.11 created [Rustam Z. Khaliullin]
846!> \author Rustam Z. Khaliullin
847! **************************************************************************************************
848 SUBROUTINE almo_scf_construct_quencher(qs_env, almo_scf_env)
849
850 TYPE(qs_environment_type), POINTER :: qs_env
851 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env
852
853 CHARACTER(len=*), PARAMETER :: routinen = 'almo_scf_construct_quencher'
854
855 CHARACTER :: sym
856 INTEGER :: col, contact_atom_1, contact_atom_2, domain_col, domain_map_local_entries, &
857 domain_row, global_entries, global_list_length, grid1, groupid, handle, hold, iatom, &
858 iatom2, iblock_col, iblock_row, idomain, idomain2, ientry, igrid, ineig, ineighbor, &
859 inode, inode2, ipair, ispin, jatom, jatom2, jdomain2, local_list_length, &
860 max_domain_neighbors, max_neig, mynode, nblkcols_tot, nblkrows_tot, nblks, ndomains, &
861 neig_temp, nnode2, nnodes, row, unit_nr
862 INTEGER, ALLOCATABLE, DIMENSION(:) :: current_number_neighbors, domain_entries_cpu, &
863 domain_map_global, domain_map_local, first_atom_of_molecule, global_list, &
864 last_atom_of_molecule, list_length_cpu, list_offset_cpu, local_list, offset_for_cpu
865 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: domain_grid, domain_neighbor_list, &
866 domain_neighbor_list_excessive
867 INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
868 LOGICAL :: already_listed, block_active, &
869 delayed_increment, found, &
870 max_neig_fails, tr
871 REAL(kind=dp) :: contact1_radius, contact2_radius, &
872 distance, distance_squared, overlap, &
873 r0, r1, s0, s1, trial_distance_squared
874 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: new_block
875 REAL(kind=dp), DIMENSION(3) :: rab
876 REAL(kind=dp), DIMENSION(:, :), POINTER :: p_old_block
877 TYPE(cell_type), POINTER :: cell
878 TYPE(cp_logger_type), POINTER :: logger
879 TYPE(dbcsr_distribution_type) :: dist
880 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
881 TYPE(dbcsr_type) :: matrix_s_sym
882 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
883 TYPE(mp_comm_type) :: group
885 DIMENSION(:), POINTER :: nl_iterator, nl_iterator2
886 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
887 POINTER :: sab_almo
888 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
889
890 CALL timeset(routinen, handle)
891
892 ! get a useful output_unit
893 logger => cp_get_default_logger()
894 IF (logger%para_env%is_source()) THEN
895 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
896 ELSE
897 unit_nr = -1
898 END IF
899
900 ndomains = almo_scf_env%ndomains
901
902 CALL get_qs_env(qs_env=qs_env, &
903 particle_set=particle_set, &
904 molecule_set=molecule_set, &
905 cell=cell, &
906 matrix_s=matrix_s, &
907 sab_almo=sab_almo)
908
909 ! if we are dealing with molecules get info about them
910 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular .OR. &
911 almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
912 ALLOCATE (first_atom_of_molecule(almo_scf_env%nmolecules))
913 ALLOCATE (last_atom_of_molecule(almo_scf_env%nmolecules))
914 CALL get_molecule_set_info(molecule_set, &
915 mol_to_first_atom=first_atom_of_molecule, &
916 mol_to_last_atom=last_atom_of_molecule)
917 END IF
918
919 ! create a symmetrized copy of the ao overlap
920 CALL dbcsr_create(matrix_s_sym, &
921 template=almo_scf_env%matrix_s(1), &
922 matrix_type=dbcsr_type_no_symmetry)
923 CALL dbcsr_get_info(almo_scf_env%matrix_s(1), &
924 matrix_type=sym)
925 IF (sym == dbcsr_type_no_symmetry) THEN
926 CALL dbcsr_copy(matrix_s_sym, almo_scf_env%matrix_s(1))
927 ELSE
928 CALL dbcsr_desymmetrize(almo_scf_env%matrix_s(1), &
929 matrix_s_sym)
930 END IF
931
932 ALLOCATE (almo_scf_env%quench_t(almo_scf_env%nspins))
933 ALLOCATE (almo_scf_env%domain_map(almo_scf_env%nspins))
934
935 DO ispin = 1, almo_scf_env%nspins
936
937 ! create the sparsity template for the occupied orbitals
938 CALL matrix_almo_create(matrix_new=almo_scf_env%quench_t(ispin), &
939 matrix_qs=matrix_s(1)%matrix, &
940 almo_scf_env=almo_scf_env, &
941 name_new="T_QUENCHER", &
943 symmetry_new=dbcsr_type_no_symmetry, &
944 spin_key=ispin, &
945 init_domains=.false.)
946
947 ! initialize distance quencher
948 CALL dbcsr_work_create(almo_scf_env%quench_t(ispin), work_mutable=.true.)
949 CALL dbcsr_get_info(almo_scf_env%quench_t(ispin), distribution=dist, &
950 nblkrows_total=nblkrows_tot, nblkcols_total=nblkcols_tot)
951 CALL dbcsr_distribution_get(dist, numnodes=nnodes, group=groupid, mynode=mynode)
952 CALL group%set_handle(groupid)
953
954 ! create global atom neighbor list from the local lists
955 ! first, calculate number of local pairs
956 local_list_length = 0
957 CALL neighbor_list_iterator_create(nl_iterator, sab_almo)
958 DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
959 ! nnode - total number of neighbors for iatom
960 ! inode - current neighbor count
961 CALL get_iterator_info(nl_iterator, &
962 iatom=iatom2, jatom=jatom2, inode=inode2, nnode=nnode2)
963 IF (inode2 == 1) THEN
964 local_list_length = local_list_length + nnode2
965 END IF
966 END DO
967 CALL neighbor_list_iterator_release(nl_iterator)
968
969 ! second, extract the local list to an array
970 ALLOCATE (local_list(2*local_list_length))
971 local_list(:) = 0
972 local_list_length = 0
973 CALL neighbor_list_iterator_create(nl_iterator2, sab_almo)
974 DO WHILE (neighbor_list_iterate(nl_iterator2) == 0)
975 CALL get_iterator_info(nl_iterator2, &
976 iatom=iatom2, jatom=jatom2)
977 local_list(2*local_list_length + 1) = iatom2
978 local_list(2*local_list_length + 2) = jatom2
979 local_list_length = local_list_length + 1
980 END DO ! end loop over pairs of atoms
981 CALL neighbor_list_iterator_release(nl_iterator2)
982
983 ! third, communicate local length to the other nodes
984 ALLOCATE (list_length_cpu(nnodes), list_offset_cpu(nnodes))
985 CALL group%allgather(2*local_list_length, list_length_cpu)
986
987 ! fourth, create a global list
988 list_offset_cpu(1) = 0
989 DO inode = 2, nnodes
990 list_offset_cpu(inode) = list_offset_cpu(inode - 1) + &
991 list_length_cpu(inode - 1)
992 END DO
993 global_list_length = list_offset_cpu(nnodes) + list_length_cpu(nnodes)
994
995 ! fifth, communicate all list data
996 ALLOCATE (global_list(global_list_length))
997 CALL group%allgatherv(local_list, global_list, &
998 list_length_cpu, list_offset_cpu)
999 DEALLOCATE (list_length_cpu, list_offset_cpu)
1000 DEALLOCATE (local_list)
1001
1002 ! calculate maximum number of atoms surrounding the domain
1003 ALLOCATE (current_number_neighbors(almo_scf_env%ndomains))
1004 current_number_neighbors(:) = 0
1005 global_list_length = global_list_length/2
1006 DO ipair = 1, global_list_length
1007 iatom2 = global_list(2*(ipair - 1) + 1)
1008 jatom2 = global_list(2*(ipair - 1) + 2)
1009 idomain2 = almo_scf_env%domain_index_of_atom(iatom2)
1010 jdomain2 = almo_scf_env%domain_index_of_atom(jatom2)
1011 ! add to the list
1012 current_number_neighbors(idomain2) = current_number_neighbors(idomain2) + 1
1013 ! add j,i with i,j
1014 IF (idomain2 /= jdomain2) THEN
1015 current_number_neighbors(jdomain2) = current_number_neighbors(jdomain2) + 1
1016 END IF
1017 END DO
1018 max_domain_neighbors = maxval(current_number_neighbors)
1019
1020 ! use the global atom neighbor list to create a global domain neighbor list
1021 ALLOCATE (domain_neighbor_list_excessive(ndomains, max_domain_neighbors))
1022 current_number_neighbors(:) = 1
1023 DO ipair = 1, ndomains
1024 domain_neighbor_list_excessive(ipair, 1) = ipair
1025 END DO
1026 DO ipair = 1, global_list_length
1027 iatom2 = global_list(2*(ipair - 1) + 1)
1028 jatom2 = global_list(2*(ipair - 1) + 2)
1029 idomain2 = almo_scf_env%domain_index_of_atom(iatom2)
1030 jdomain2 = almo_scf_env%domain_index_of_atom(jatom2)
1031 already_listed = .false.
1032 DO ineighbor = 1, current_number_neighbors(idomain2)
1033 IF (domain_neighbor_list_excessive(idomain2, ineighbor) == jdomain2) THEN
1034 already_listed = .true.
1035 EXIT
1036 END IF
1037 END DO
1038 IF (.NOT. already_listed) THEN
1039 ! add to the list
1040 current_number_neighbors(idomain2) = current_number_neighbors(idomain2) + 1
1041 domain_neighbor_list_excessive(idomain2, current_number_neighbors(idomain2)) = jdomain2
1042 ! add j,i with i,j
1043 IF (idomain2 /= jdomain2) THEN
1044 current_number_neighbors(jdomain2) = current_number_neighbors(jdomain2) + 1
1045 domain_neighbor_list_excessive(jdomain2, current_number_neighbors(jdomain2)) = idomain2
1046 END IF
1047 END IF
1048 END DO ! end loop over pairs of atoms
1049 DEALLOCATE (global_list)
1050
1051 max_domain_neighbors = maxval(current_number_neighbors)
1052 ALLOCATE (domain_neighbor_list(ndomains, max_domain_neighbors))
1053 domain_neighbor_list(:, :) = 0
1054 domain_neighbor_list(:, :) = domain_neighbor_list_excessive(:, 1:max_domain_neighbors)
1055 DEALLOCATE (domain_neighbor_list_excessive)
1056
1057 ALLOCATE (almo_scf_env%domain_map(ispin)%index1(ndomains))
1058 ALLOCATE (almo_scf_env%domain_map(ispin)%pairs(max_domain_neighbors*ndomains, 2))
1059 almo_scf_env%domain_map(ispin)%pairs(:, :) = 0
1060 almo_scf_env%domain_map(ispin)%index1(:) = 0
1061 domain_map_local_entries = 0
1062
1063 ! RZK-warning intermediate [0,1] quencher values are ill-defined
1064 ! for molecules (not continuous and conceptually inadequate)
1065
1066 CALL dbcsr_get_info(almo_scf_env%quench_t(ispin), &
1067 row_blk_size=row_blk_size, col_blk_size=col_blk_size)
1068 ! O(N) loop over domain pairs
1069 DO row = 1, nblkrows_tot
1070 DO col = 1, current_number_neighbors(row)
1071 tr = .false.
1072 iblock_row = row
1073 iblock_col = domain_neighbor_list(row, col)
1074 CALL dbcsr_get_stored_coordinates(almo_scf_env%quench_t(ispin), &
1075 iblock_row, iblock_col, hold)
1076
1077 IF (hold == mynode) THEN
1078
1079 ! Translate indices of distribution blocks to indices of domain blocks
1080 ! Rows are AOs
1081 domain_row = almo_scf_env%domain_index_of_ao_block(iblock_row)
1082 ! Columns are electrons (i.e. MOs)
1083 domain_col = almo_scf_env%domain_index_of_mo_block(iblock_col)
1084
1085 SELECT CASE (almo_scf_env%constraint_type)
1087
1088 block_active = .false.
1089 ! type of electron groups
1090 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
1091
1092 ! type of ao domains
1093 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1094
1095 ! ao domains are molecular / electron groups are molecular
1096 IF (domain_row == domain_col) THEN
1097 block_active = .true.
1098 END IF
1099
1100 ELSE ! ao domains are atomic
1101
1102 ! ao domains are atomic / electron groups are molecular
1103 cpabort("Illegal: atomic domains and molecular groups")
1104
1105 END IF
1106
1107 ELSE ! electron groups are atomic
1108
1109 ! type of ao domains
1110 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1111
1112 ! ao domains are molecular / electron groups are atomic
1113 cpabort("Illegal: molecular domains and atomic groups")
1114
1115 ELSE
1116
1117 ! ao domains are atomic / electron groups are atomic
1118 IF (domain_row == domain_col) THEN
1119 block_active = .true.
1120 END IF
1121
1122 END IF
1123
1124 END IF ! end type of electron groups
1125
1126 IF (block_active) THEN
1127
1128 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1129 new_block(:, :) = 1.0_dp
1130 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1131 DEALLOCATE (new_block)
1132
1133 IF (domain_map_local_entries >= max_domain_neighbors*almo_scf_env%ndomains) THEN
1134 cpabort("weird... max_domain_neighbors is exceeded")
1135 END IF
1136 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1137 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1138 domain_map_local_entries = domain_map_local_entries + 1
1139
1140 END IF
1141
1143
1144 ! type of electron groups
1145 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
1146
1147 ! type of ao domains
1148 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1149
1150 ! ao domains are molecular / electron groups are molecular
1151
1152 ! compute the maximum overlap between the atoms of the two molecules
1153 CALL dbcsr_get_block_p(matrix_s_sym, iblock_row, iblock_col, p_old_block, found)
1154 IF (found) THEN
1155 overlap = maxval(abs(p_old_block))
1156 ELSE
1157 overlap = 0.0_dp
1158 END IF
1159
1160 ELSE ! ao domains are atomic
1161
1162 ! ao domains are atomic / electron groups are molecular
1163 ! overlap_between_atom_and_molecule(atom=domain_row,molecule=domain_col)
1164 cpabort("atomic domains and molecular groups - NYI")
1165
1166 END IF
1167
1168 ELSE ! electron groups are atomic
1169
1170 ! type of ao domains
1171 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1172
1173 ! ao domains are molecular / electron groups are atomic
1174 ! overlap_between_atom_and_molecule(atom=domain_col,molecule=domain_row)
1175 cpabort("molecular domains and atomic groups - NYI")
1176
1177 ELSE
1178
1179 ! ao domains are atomic / electron groups are atomic
1180 ! compute max overlap between atoms: domain_row and domain_col
1181 CALL dbcsr_get_block_p(matrix_s_sym, iblock_row, iblock_col, p_old_block, found)
1182 IF (found) THEN
1183 overlap = maxval(abs(p_old_block))
1184 ELSE
1185 overlap = 0.0_dp
1186 END IF
1187
1188 END IF
1189
1190 END IF ! end type of electron groups
1191
1192 s0 = -log10(abs(almo_scf_env%quencher_s0))
1193 s1 = -log10(abs(almo_scf_env%quencher_s1))
1194 IF (overlap == 0.0_dp) THEN
1195 overlap = -log10(abs(almo_scf_env%eps_filter)) + 100.0_dp
1196 ELSE
1197 overlap = -log10(overlap)
1198 END IF
1199 IF (s0 < 0.0_dp) THEN
1200 cpabort("S0 is less than zero")
1201 END IF
1202 IF (s1 <= 0.0_dp) THEN
1203 cpabort("S1 is less than or equal to zero")
1204 END IF
1205 IF (s0 >= s1) THEN
1206 cpabort("S0 is greater than or equal to S1")
1207 END IF
1208
1209 ! Fill in non-zero blocks if AOs are close to the electron center
1210 IF (overlap < s1) THEN
1211 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1212 IF (overlap <= s0) THEN
1213 new_block(:, :) = 1.0_dp
1214 ELSE
1215 new_block(:, :) = 1.0_dp/(1.0_dp + exp(-(s0 - s1)/(s0 - overlap) - (s0 - s1)/(overlap - s1)))
1216 END IF
1217
1218 IF (abs(new_block(1, 1)) > abs(almo_scf_env%eps_filter)) THEN
1219 IF (domain_map_local_entries >= max_domain_neighbors*almo_scf_env%ndomains) THEN
1220 cpabort("weird... max_domain_neighbors is exceeded")
1221 END IF
1222 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1223 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1224 domain_map_local_entries = domain_map_local_entries + 1
1225 END IF
1226
1227 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1228 DEALLOCATE (new_block)
1229
1230 END IF
1231
1233
1234 ! type of electron groups
1235 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
1236
1237 ! type of ao domains
1238 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1239
1240 ! ao domains are molecular / electron groups are molecular
1241
1242 ! compute distance between molecules: domain_row and domain_col
1243 ! distance between molecules is defined as the smallest
1244 ! distance among all atom pairs
1245 IF (domain_row == domain_col) THEN
1246 distance = 0.0_dp
1247 contact_atom_1 = first_atom_of_molecule(domain_row)
1248 contact_atom_2 = first_atom_of_molecule(domain_col)
1249 ELSE
1250 distance_squared = 1.0e+100_dp
1251 contact_atom_1 = -1
1252 contact_atom_2 = -1
1253 DO iatom = first_atom_of_molecule(domain_row), last_atom_of_molecule(domain_row)
1254 DO jatom = first_atom_of_molecule(domain_col), last_atom_of_molecule(domain_col)
1255 rab(:) = pbc(particle_set(iatom)%r(:), particle_set(jatom)%r(:), cell)
1256 trial_distance_squared = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
1257 IF (trial_distance_squared < distance_squared) THEN
1258 distance_squared = trial_distance_squared
1259 contact_atom_1 = iatom
1260 contact_atom_2 = jatom
1261 END IF
1262 END DO ! jatom
1263 END DO ! iatom
1264 cpassert(contact_atom_1 > 0)
1265 distance = sqrt(distance_squared)
1266 END IF
1267
1268 ELSE ! ao domains are atomic
1269
1270 ! ao domains are atomic / electron groups are molecular
1271 !distance_between_atom_and_molecule(atom=domain_row,molecule=domain_col)
1272 cpabort("atomic domains and molecular groups - NYI")
1273
1274 END IF
1275
1276 ELSE ! electron groups are atomic
1277
1278 ! type of ao domains
1279 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1280
1281 ! ao domains are molecular / electron groups are atomic
1282 !distance_between_atom_and_molecule(atom=domain_col,molecule=domain_row)
1283 cpabort("molecular domains and atomic groups - NYI")
1284
1285 ELSE
1286
1287 ! ao domains are atomic / electron groups are atomic
1288 ! compute distance between atoms: domain_row and domain_col
1289 rab(:) = pbc(particle_set(domain_row)%r(:), particle_set(domain_col)%r(:), cell)
1290 distance = sqrt(rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3))
1291 contact_atom_1 = domain_row
1292 contact_atom_2 = domain_col
1293
1294 END IF
1295
1296 END IF ! end type of electron groups
1297
1298 ! get atomic radii to compute distance cutoff threshold
1299 IF (almo_scf_env%quencher_radius_type == do_bondparm_covalent) THEN
1300 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_1)%atomic_kind, &
1301 rcov=contact1_radius)
1302 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_2)%atomic_kind, &
1303 rcov=contact2_radius)
1304 ELSE IF (almo_scf_env%quencher_radius_type == do_bondparm_vdw) THEN
1305 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_1)%atomic_kind, &
1306 rvdw=contact1_radius)
1307 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_2)%atomic_kind, &
1308 rvdw=contact2_radius)
1309 ELSE
1310 cpabort("Illegal quencher_radius_type")
1311 END IF
1312 contact1_radius = cp_unit_to_cp2k(contact1_radius, "angstrom")
1313 contact2_radius = cp_unit_to_cp2k(contact2_radius, "angstrom")
1314
1315 !RZK-warning the procedure is faulty for molecules:
1316 ! the closest contacts should be found using
1317 ! the element specific radii
1318
1319 ! compute inner and outer cutoff radii
1320 r0 = almo_scf_env%quencher_r0_factor*(contact1_radius + contact2_radius)
1321 !+almo_scf_env%quencher_r0_shift
1322 r1 = almo_scf_env%quencher_r1_factor*(contact1_radius + contact2_radius)
1323 !+almo_scf_env%quencher_r1_shift
1324
1325 IF (r0 < 0.0_dp) THEN
1326 cpabort("R0 is less than zero")
1327 END IF
1328 IF (r1 <= 0.0_dp) THEN
1329 cpabort("R1 is less than or equal to zero")
1330 END IF
1331 IF (r0 > r1) THEN
1332 cpabort("R0 is greater than or equal to R1")
1333 END IF
1334
1335 ! Fill in non-zero blocks if AOs are close to the electron center
1336 IF (distance < r1) THEN
1337 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1338 IF (distance <= r0) THEN
1339 new_block(:, :) = 1.0_dp
1340 ELSE
1341 ! remove the intermediate values from the quencher temporarily
1342 cpabort("")
1343 new_block(:, :) = 1.0_dp/(1.0_dp + exp((r1 - r0)/(r0 - distance) + (r1 - r0)/(r1 - distance)))
1344 END IF
1345
1346 IF (abs(new_block(1, 1)) > abs(almo_scf_env%eps_filter)) THEN
1347 IF (domain_map_local_entries >= max_domain_neighbors*almo_scf_env%ndomains) THEN
1348 cpabort("weird... max_domain_neighbors is exceeded")
1349 END IF
1350 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1351 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1352 domain_map_local_entries = domain_map_local_entries + 1
1353 END IF
1354
1355 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1356 DEALLOCATE (new_block)
1357 END IF
1358
1359 CASE DEFAULT
1360 cpabort("Illegal constraint type")
1361 END SELECT
1362
1363 END IF ! mynode
1364
1365 END DO
1366 END DO ! end O(N) loop over pairs
1367
1368 DEALLOCATE (domain_neighbor_list)
1369 DEALLOCATE (current_number_neighbors)
1370
1371 CALL dbcsr_finalize(almo_scf_env%quench_t(ispin))
1372
1373 CALL dbcsr_filter(almo_scf_env%quench_t(ispin), &
1374 almo_scf_env%eps_filter)
1375
1376 ! check that both domain_map and quench_t have the same number of entries
1377 nblks = dbcsr_get_num_blocks(almo_scf_env%quench_t(ispin))
1378 IF (nblks /= domain_map_local_entries) THEN
1379 cpabort("number of blocks is wrong")
1380 END IF
1381
1382 ! first, communicate map sizes on the other nodes
1383 ALLOCATE (domain_entries_cpu(nnodes), offset_for_cpu(nnodes))
1384 CALL group%allgather(2*domain_map_local_entries, domain_entries_cpu)
1385
1386 ! second, create
1387 offset_for_cpu(1) = 0
1388 DO inode = 2, nnodes
1389 offset_for_cpu(inode) = offset_for_cpu(inode - 1) + &
1390 domain_entries_cpu(inode - 1)
1391 END DO
1392 global_entries = offset_for_cpu(nnodes) + domain_entries_cpu(nnodes)
1393
1394 ! communicate all entries
1395 ALLOCATE (domain_map_global(global_entries))
1396 ALLOCATE (domain_map_local(2*domain_map_local_entries))
1397 DO ientry = 1, domain_map_local_entries
1398 domain_map_local(2*(ientry - 1) + 1) = almo_scf_env%domain_map(ispin)%pairs(ientry, 1)
1399 domain_map_local(2*ientry) = almo_scf_env%domain_map(ispin)%pairs(ientry, 2)
1400 END DO
1401 CALL group%allgatherv(domain_map_local, domain_map_global, &
1402 domain_entries_cpu, offset_for_cpu)
1403 DEALLOCATE (domain_entries_cpu, offset_for_cpu)
1404 DEALLOCATE (domain_map_local)
1405
1406 DEALLOCATE (almo_scf_env%domain_map(ispin)%index1)
1407 DEALLOCATE (almo_scf_env%domain_map(ispin)%pairs)
1408 ALLOCATE (almo_scf_env%domain_map(ispin)%index1(ndomains))
1409 ALLOCATE (almo_scf_env%domain_map(ispin)%pairs(global_entries/2, 2))
1410 almo_scf_env%domain_map(ispin)%pairs(:, :) = 0
1411 almo_scf_env%domain_map(ispin)%index1(:) = 0
1412
1413 ! unpack the received data into a local variable
1414 ! since we do not know the maximum global number of neighbors
1415 ! try one. if fails increase the maximum number and try again
1416 ! until it succeeds
1417 max_neig = max_domain_neighbors
1418 max_neig_fails = .true.
1419 max_neig_loop: DO WHILE (max_neig_fails)
1420 ALLOCATE (domain_grid(almo_scf_env%ndomains, 0:max_neig))
1421 domain_grid(:, :) = 0
1422 ! init the number of collected neighbors
1423 domain_grid(:, 0) = 1
1424 ! loop over the records
1425 global_entries = global_entries/2
1426 DO ientry = 1, global_entries
1427 ! get the center
1428 grid1 = domain_map_global(2*ientry)
1429 ! get the neighbor
1430 ineig = domain_map_global(2*(ientry - 1) + 1)
1431 ! check boundaries
1432 IF (domain_grid(grid1, 0) > max_neig) THEN
1433 ! this neighbor will overstep the boundaries
1434 ! stop the trial and increase the max number of neighbors
1435 DEALLOCATE (domain_grid)
1436 max_neig = max_neig*2
1437 cycle max_neig_loop
1438 END IF
1439 ! for the current center loop over the collected neighbors
1440 ! to insert the current record in a numerical order
1441 delayed_increment = .false.
1442 DO igrid = 1, domain_grid(grid1, 0)
1443 ! compare the current neighbor with that already in the 'book'
1444 IF (ineig < domain_grid(grid1, igrid)) THEN
1445 ! if this one is smaller then insert it here and pick up the one
1446 ! from the book to continue inserting
1447 neig_temp = ineig
1448 ineig = domain_grid(grid1, igrid)
1449 domain_grid(grid1, igrid) = neig_temp
1450 ELSE
1451 IF (domain_grid(grid1, igrid) == 0) THEN
1452 ! got the empty slot now - insert the record
1453 domain_grid(grid1, igrid) = ineig
1454 ! increase the record counter but do it outside the loop
1455 delayed_increment = .true.
1456 END IF
1457 END IF
1458 END DO
1459 IF (delayed_increment) THEN
1460 domain_grid(grid1, 0) = domain_grid(grid1, 0) + 1
1461 ELSE
1462 ! should not be here - all records must be inserted
1463 cpabort("all records must be inserted")
1464 END IF
1465 END DO
1466 max_neig_fails = .false.
1467 END DO max_neig_loop
1468 DEALLOCATE (domain_map_global)
1469
1470 ientry = 1
1471 DO idomain = 1, almo_scf_env%ndomains
1472 DO ineig = 1, domain_grid(idomain, 0) - 1
1473 almo_scf_env%domain_map(ispin)%pairs(ientry, 1) = domain_grid(idomain, ineig)
1474 almo_scf_env%domain_map(ispin)%pairs(ientry, 2) = idomain
1475 ientry = ientry + 1
1476 END DO
1477 almo_scf_env%domain_map(ispin)%index1(idomain) = ientry
1478 END DO
1479 DEALLOCATE (domain_grid)
1480
1481 END DO ! ispin
1482 IF (almo_scf_env%nspins == 2) THEN
1483 CALL dbcsr_copy(almo_scf_env%quench_t(2), &
1484 almo_scf_env%quench_t(1))
1485 almo_scf_env%domain_map(2)%pairs(:, :) = &
1486 almo_scf_env%domain_map(1)%pairs(:, :)
1487 almo_scf_env%domain_map(2)%index1(:) = &
1488 almo_scf_env%domain_map(1)%index1(:)
1489 END IF
1490
1491 CALL dbcsr_release(matrix_s_sym)
1492
1493 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular .OR. &
1494 almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1495 DEALLOCATE (first_atom_of_molecule)
1496 DEALLOCATE (last_atom_of_molecule)
1497 END IF
1498
1499 !mynode = dbcsr_mp_mynode(dbcsr_distribution_mp(&
1500 ! dbcsr_distribution(almo_scf_env%quench_t(ispin))))
1501 !CALL dbcsr_get_info(almo_scf_env%quench_t(ispin), distribution=dist, &
1502 ! nblkrows_total=nblkrows_tot, nblkcols_total=nblkcols_tot)
1503 !DO row = 1, nblkrows_tot
1504 ! DO col = 1, nblkcols_tot
1505 ! tr = .FALSE.
1506 ! iblock_row = row
1507 ! iblock_col = col
1508 ! CALL dbcsr_get_stored_coordinates(almo_scf_env%quench_t(ispin),&
1509 ! iblock_row, iblock_col, tr, hold)
1510 ! CALL dbcsr_get_block_p(almo_scf_env%quench_t(ispin),&
1511 ! row, col, p_old_block, found)
1512 ! write(*,*) "RST_NOTE:", mynode, row, col, hold, found
1513 ! ENDDO
1514 !ENDDO
1515
1516 CALL timestop(handle)
1517
1518 END SUBROUTINE almo_scf_construct_quencher
1519
1520! *****************************************************************************
1521!> \brief Compute matrix W (energy-weighted density matrix) that is needed
1522!> for the evaluation of forces
1523!> \param matrix_w ...
1524!> \param almo_scf_env ...
1525!> \par History
1526!> 2015.03 created [Rustam Z. Khaliullin]
1527!> \author Rustam Z. Khaliullin
1528! **************************************************************************************************
1529 SUBROUTINE calculate_w_matrix_almo(matrix_w, almo_scf_env)
1530 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_w
1531 TYPE(almo_scf_env_type) :: almo_scf_env
1532
1533 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_w_matrix_almo'
1534
1535 INTEGER :: handle, ispin
1536 REAL(kind=dp) :: scaling
1537 TYPE(dbcsr_type) :: tmp_nn1, tmp_no1, tmp_oo1, tmp_oo2
1538
1539 CALL timeset(routinen, handle)
1540
1541 IF (almo_scf_env%nspins == 1) THEN
1542 scaling = 2.0_dp
1543 ELSE
1544 scaling = 1.0_dp
1545 END IF
1546
1547 DO ispin = 1, almo_scf_env%nspins
1548
1549 CALL dbcsr_create(tmp_nn1, template=almo_scf_env%matrix_s(1), &
1550 matrix_type=dbcsr_type_no_symmetry)
1551 CALL dbcsr_create(tmp_no1, template=almo_scf_env%matrix_t(ispin), &
1552 matrix_type=dbcsr_type_no_symmetry)
1553 CALL dbcsr_create(tmp_oo1, template=almo_scf_env%matrix_sigma_inv(ispin), &
1554 matrix_type=dbcsr_type_no_symmetry)
1555 CALL dbcsr_create(tmp_oo2, template=almo_scf_env%matrix_sigma_inv(ispin), &
1556 matrix_type=dbcsr_type_no_symmetry)
1557
1558 CALL dbcsr_copy(tmp_nn1, almo_scf_env%matrix_ks(ispin))
1559 ! 1. TMP_NO1=F.T
1560 CALL dbcsr_multiply("N", "N", scaling, tmp_nn1, almo_scf_env%matrix_t(ispin), &
1561 0.0_dp, tmp_no1, filter_eps=almo_scf_env%eps_filter)
1562 ! 2. TMP_OO1=T^(tr).TMP_NO1=T^(tr).(FT)
1563 CALL dbcsr_multiply("T", "N", 1.0_dp, almo_scf_env%matrix_t(ispin), tmp_no1, &
1564 0.0_dp, tmp_oo1, filter_eps=almo_scf_env%eps_filter)
1565 ! 3. TMP_OO2=TMP_OO1.siginv=(T^(tr)FT).siginv
1566 CALL dbcsr_multiply("N", "N", 1.0_dp, tmp_oo1, almo_scf_env%matrix_sigma_inv(ispin), &
1567 0.0_dp, tmp_oo2, filter_eps=almo_scf_env%eps_filter)
1568 ! 4. TMP_OO1=siginv.TMP_OO2=siginv.(T^(tr)FTsiginv)
1569 CALL dbcsr_multiply("N", "N", 1.0_dp, almo_scf_env%matrix_sigma_inv(ispin), tmp_oo2, &
1570 0.0_dp, tmp_oo1, filter_eps=almo_scf_env%eps_filter)
1571 ! 5. TMP_NO1=T.TMP_OO1.=T.(siginvT^(tr)FTsiginv)
1572 CALL dbcsr_multiply("N", "N", 1.0_dp, almo_scf_env%matrix_t(ispin), tmp_oo1, &
1573 0.0_dp, tmp_no1, filter_eps=almo_scf_env%eps_filter)
1574 ! 6. TMP_NN1=TMP_NO1.T^(tr)=(TsiginvT^(tr)FTsiginv).T^(tr)=RFR
1575 CALL dbcsr_multiply("N", "T", 1.0_dp, tmp_no1, almo_scf_env%matrix_t(ispin), &
1576 0.0_dp, tmp_nn1, filter_eps=almo_scf_env%eps_filter)
1577 CALL matrix_almo_to_qs(tmp_nn1, matrix_w(ispin)%matrix, almo_scf_env%mat_distr_aos)
1578
1579 CALL dbcsr_release(tmp_nn1)
1580 CALL dbcsr_release(tmp_no1)
1581 CALL dbcsr_release(tmp_oo1)
1582 CALL dbcsr_release(tmp_oo2)
1583
1584 END DO
1585
1586 CALL timestop(handle)
1587
1588 END SUBROUTINE calculate_w_matrix_almo
1589
1590END MODULE almo_scf_qs
1591
double distance_squared(double *A, double *B)
Definition grpp_utils.c:68
Interface between ALMO SCF and QS.
Definition almo_scf_qs.F:14
subroutine, public almo_scf_update_ks_energy(qs_env, energy, energy_singles_corr)
update qs_env total energy
subroutine, public construct_qs_mos(qs_env, almo_scf_env)
Create MOs in the QS env to be able to return ALMOs to QS.
subroutine, public almo_dm_to_almo_ks(qs_env, matrix_p, matrix_ks, energy_total, eps_filter, mat_distr_aos, smear, kts_sum)
uses the ALMO density matrix to compute ALMO KS matrix and the new energy
subroutine, public calculate_w_matrix_almo(matrix_w, almo_scf_env)
Compute matrix W (energy-weighted density matrix) that is needed for the evaluation of forces.
subroutine, public matrix_almo_create(matrix_new, matrix_qs, almo_scf_env, name_new, size_keys, symmetry_new, spin_key, init_domains)
create the ALMO matrix templates
subroutine, public init_almo_ks_matrix_via_qs(qs_env, matrix_ks, mat_distr_aos, eps_filter)
Initialization of the QS and ALMO KS matrix.
subroutine, public almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
return density matrix to the qs_env
subroutine, public matrix_qs_to_almo(matrix_qs, matrix_almo, mat_distr_aos)
convert between two types of matrices: QS style to ALMO style
subroutine, public almo_scf_construct_quencher(qs_env, almo_scf_env)
Creates the matrix that imposes absolute locality on MOs.
Types for all ALMO-based methods.
integer, parameter, public almo_mat_dim_occ
integer, parameter, public almo_mat_dim_virt_full
integer, parameter, public almo_mat_dim_aobasis
integer, parameter, public almo_mat_dim_virt
integer, parameter, public almo_mat_dim_virt_disc
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Handles all functions related to the CELL.
Definition cell_types.F:15
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_distribution_release(dist)
...
subroutine, public dbcsr_distribution_new(dist, template, group, pgrid, row_dist, col_dist, reuse_arrays)
...
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_get_block_p(matrix, row, col, block, found, row_size, col_size)
...
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)
...
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)
...
subroutine, public dbcsr_get_stored_coordinates(matrix, row, column, processor)
...
subroutine, public dbcsr_work_create(matrix, nblks_guess, sizedata_guess, n, work_mutable)
...
subroutine, public dbcsr_filter(matrix, eps)
...
subroutine, public dbcsr_finalize(matrix)
...
subroutine, public dbcsr_set(matrix, alpha)
...
subroutine, public dbcsr_release(matrix)
...
subroutine, public dbcsr_complete_redistribute(matrix, redist)
...
integer function, public dbcsr_get_num_blocks(matrix)
...
subroutine, public dbcsr_put_block(matrix, row, col, block, summation)
...
subroutine, public dbcsr_distribution_get(dist, row_dist, col_dist, nrows, ncols, has_threads, group, mynode, numnodes, nprows, npcols, myprow, mypcol, pgrid, subgroups_defined, prow_group, pcol_group)
...
DBCSR operations in CP2K.
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_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 ...
recursive integer function, public cp_logger_get_default_unit_nr(logger, local, skip_not_ionode)
asks the default unit number of the given logger. try to use cp_logger_get_unit_nr
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
unit conversion facility
Definition cp_units.F:30
real(kind=dp) function, public cp_unit_to_cp2k(value, unit_str, defaults, power)
converts to the internal cp2k units to the given unit
Definition cp_units.F:1210
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public almo_constraint_distance
integer, parameter, public do_bondparm_covalent
integer, parameter, public almo_mat_distr_molecular
integer, parameter, public almo_domain_layout_molecular
integer, parameter, public do_bondparm_vdw
integer, parameter, public almo_mat_distr_atomic
integer, parameter, public almo_constraint_block_diagonal
integer, parameter, public almo_constraint_ao_overlap
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Interface to the message passing library MPI.
Define the data structure for the molecule information.
subroutine, public get_molecule_set_info(molecule_set, atom_to_mol, mol_to_first_atom, mol_to_last_atom, mol_to_nelectrons, mol_to_nbasis, mol_to_charge, mol_to_multiplicity)
returns information about molecules in the set.
Define the data structure for the particle information.
Perform a QUICKSTEP wavefunction optimization (single point)
Definition qs_energy.F:14
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, 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.
subroutine, public set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, mimic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, rhoz_cneo_set, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, harris_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, wanniercentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Set the QUICKSTEP environment.
routines that build the Kohn-Sham matrix (i.e calculate the coulomb and xc parts
subroutine, public qs_ks_update_qs_env(qs_env, calculate_forces, just_energy, print_active)
updates the Kohn Sham matrix of the given qs_env (facility method)
subroutine, public qs_ks_did_change(ks_env, s_mstruct_changed, rho_changed, potential_changed, full_reset)
tells that some of the things relevant to the ks calculation did change. has to be called when change...
subroutine, public set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, exc_accint, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, kpoints, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
subroutine, public init_mo_set(mo_set, fm_pool, fm_ref, fm_struct, name, counter)
initializes an allocated mo_set. eigenvalues, mo_coeff, occupation_numbers are valid only after this ...
subroutine, public allocate_mo_set(mo_set, nao, nmo, nelectron, n_el_f, maxocc, flexible_electron_count)
Allocates a mo set and partially initializes it (nao,nmo,nelectron, and flexible_electron_count are v...
subroutine, public deallocate_mo_set(mo_set)
Deallocate a wavefunction data structure.
Define the neighbor list data types and the corresponding functionality.
subroutine, public neighbor_list_iterator_create(iterator_set, nl, search, nthread)
Neighbor list iterator functions.
subroutine, public neighbor_list_iterator_release(iterator_set)
...
integer function, public neighbor_list_iterate(iterator_set, mepos)
...
subroutine, public get_iterator_info(iterator_set, mepos, ikind, jkind, nkind, ilist, nlist, inode, nnode, iatom, jatom, r, cell)
...
methods of the rho structure (defined in qs_rho_types)
subroutine, public qs_rho_update_rho(rho_struct, qs_env, rho_xc_external, local_rho_set, task_list_external, task_list_external_soft, pw_env_external, para_env_external)
updates rho_r and rho_g to the rhorho_ao. if use_kinetic_energy_density also computes tau_r and tau_g...
superstucture that hold various representations of the density and keeps track of which ones are vali...
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
module that contains the definitions of the scf types
subroutine, public scf_env_create(scf_env)
allocates and initialize an scf_env
Type defining parameters related to the simulation cell.
Definition cell_types.F:60
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...
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
keeps the density in various representations, keeping track of which ones are valid.