(git:a660c7f)
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 ! start linear-scaling replacement:
292 ! works only for molecular blocks AND molecular distributions
293 DO row = 1, nblkrows_tot
294 tr = .false.
295 iblock_row = row
296 iblock_col = row
297 CALL dbcsr_get_stored_coordinates(matrix_new, iblock_row, iblock_col, hold)
298
299 IF (hold == mynode) THEN
300
301 active = .true.
302
303 one_dim_is_mo = .false.
304 DO dimen = 1, 2 ! 1 - row, 2 - column dimension
305 IF (size_keys(dimen) == almo_mat_dim_occ) one_dim_is_mo = .true.
306 END DO
307 IF (one_dim_is_mo) THEN
308 IF (almo_scf_env%nocc_of_domain(row, spin_key) == 0) active = .false.
309 END IF
310
311 one_dim_is_mo = .false.
312 DO dimen = 1, 2
313 IF (size_keys(dimen) == almo_mat_dim_virt) one_dim_is_mo = .true.
314 END DO
315 IF (one_dim_is_mo) THEN
316 IF (almo_scf_env%nvirt_of_domain(row, spin_key) == 0) active = .false.
317 END IF
318
319 one_dim_is_mo = .false.
320 DO dimen = 1, 2
321 IF (size_keys(dimen) == almo_mat_dim_virt_disc) one_dim_is_mo = .true.
322 END DO
323 IF (one_dim_is_mo) THEN
324 IF (almo_scf_env%nvirt_disc_of_domain(row, spin_key) == 0) active = .false.
325 END IF
326
327 one_dim_is_mo = .false.
328 DO dimen = 1, 2
329 IF (size_keys(dimen) == almo_mat_dim_virt_full) one_dim_is_mo = .true.
330 END DO
331 IF (one_dim_is_mo) THEN
332 IF (almo_scf_env%nvirt_full_of_domain(row, spin_key) == 0) active = .false.
333 END IF
334
335 IF (active) THEN
336 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
337 new_block(:, :) = 1.0_dp
338 CALL dbcsr_put_block(matrix_new, iblock_row, iblock_col, new_block)
339 DEALLOCATE (new_block)
340 END IF
341
342 END IF ! mynode
343 END DO
344 ! end lnear-scaling replacement
345
346 END IF ! init_domains
347
348 CALL dbcsr_finalize(matrix_new)
349
350 CALL timestop(handle)
351
352 END SUBROUTINE matrix_almo_create
353
354! **************************************************************************************************
355!> \brief convert between two types of matrices: QS style to ALMO style
356!> \param matrix_qs ...
357!> \param matrix_almo ...
358!> \param mat_distr_aos ...
359!> \par History
360!> 2011.06 created [Rustam Z Khaliullin]
361!> \author Rustam Z Khaliullin
362! **************************************************************************************************
363 SUBROUTINE matrix_qs_to_almo(matrix_qs, matrix_almo, mat_distr_aos)
364
365 TYPE(dbcsr_type) :: matrix_qs, matrix_almo
366 INTEGER :: mat_distr_aos
367
368 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_qs_to_almo'
369
370 INTEGER :: handle
371 TYPE(dbcsr_type) :: matrix_qs_nosym
372
373 CALL timeset(routinen, handle)
374 !RZK-warning if it's not a N(AO)xN(AO) matrix then stop
375
376 SELECT CASE (mat_distr_aos)
378 ! automatic data_type conversion
379 CALL dbcsr_copy(matrix_almo, matrix_qs)
381 ! desymmetrize the qs matrix
382 CALL dbcsr_create(matrix_qs_nosym, template=matrix_qs, matrix_type=dbcsr_type_no_symmetry)
383 CALL dbcsr_desymmetrize(matrix_qs, matrix_qs_nosym)
384
385 ! perform the magic complete_redistribute
386 ! before calling complete_redistribute set all blocks to zero
387 ! otherwise the non-zero elements of the redistributed matrix,
388 ! which are in zero-blocks of the original matrix, will remain
389 ! in the final redistributed matrix. this is a bug in
390 ! complete_redistribute. RZK-warning it should be later corrected by calling
391 ! dbcsr_set to 0.0 from within complete_redistribute
392 CALL dbcsr_set(matrix_almo, 0.0_dp)
393 CALL dbcsr_complete_redistribute(matrix_qs_nosym, matrix_almo)
394 CALL dbcsr_release(matrix_qs_nosym)
395
396 CASE DEFAULT
397 cpabort("Unknown mat_distr_aos for matrix_qs_to_almo")
398 END SELECT
399
400 CALL timestop(handle)
401
402 END SUBROUTINE matrix_qs_to_almo
403
404! **************************************************************************************************
405!> \brief convert between two types of matrices: ALMO style to QS style
406!> \param matrix_almo ...
407!> \param matrix_qs ...
408!> \param mat_distr_aos ...
409!> \par History
410!> 2011.06 created [Rustam Z Khaliullin]
411!> \author Rustam Z Khaliullin
412! **************************************************************************************************
413 SUBROUTINE matrix_almo_to_qs(matrix_almo, matrix_qs, mat_distr_aos)
414 TYPE(dbcsr_type) :: matrix_almo, matrix_qs
415 INTEGER, INTENT(IN) :: mat_distr_aos
416
417 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_almo_to_qs'
418
419 INTEGER :: handle
420 TYPE(dbcsr_type) :: matrix_almo_redist
421
422 CALL timeset(routinen, handle)
423 ! RZK-warning if it's not a N(AO)xN(AO) matrix then stop
424
425 SELECT CASE (mat_distr_aos)
427 CALL dbcsr_copy(matrix_qs, matrix_almo, keep_sparsity=.true.)
429 CALL dbcsr_create(matrix_almo_redist, template=matrix_qs)
430 CALL dbcsr_complete_redistribute(matrix_almo, matrix_almo_redist)
431 CALL dbcsr_set(matrix_qs, 0.0_dp)
432 CALL dbcsr_copy(matrix_qs, matrix_almo_redist, keep_sparsity=.true.)
433 CALL dbcsr_release(matrix_almo_redist)
434 CASE DEFAULT
435 cpabort("Unknown mat_distr_aos for matrix_almo_to_qs")
436 END SELECT
437
438 CALL timestop(handle)
439
440 END SUBROUTINE matrix_almo_to_qs
441
442! **************************************************************************************************
443!> \brief Initialization of the QS and ALMO KS matrix
444!> \param qs_env ...
445!> \param matrix_ks ...
446!> \param mat_distr_aos ...
447!> \param eps_filter ...
448!> \par History
449!> 2011.05 created [Rustam Z Khaliullin]
450!> \author Rustam Z Khaliullin
451! **************************************************************************************************
452 SUBROUTINE init_almo_ks_matrix_via_qs(qs_env, matrix_ks, mat_distr_aos, eps_filter)
453
454 TYPE(qs_environment_type), POINTER :: qs_env
455 TYPE(dbcsr_type), DIMENSION(:) :: matrix_ks
456 INTEGER :: mat_distr_aos
457 REAL(kind=dp) :: eps_filter
458
459 CHARACTER(len=*), PARAMETER :: routinen = 'init_almo_ks_matrix_via_qs'
460
461 INTEGER :: handle, ispin, nspin
462 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_qs_ks, matrix_qs_s
463 TYPE(dft_control_type), POINTER :: dft_control
464 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
465 POINTER :: sab_orb
466 TYPE(qs_ks_env_type), POINTER :: ks_env
467
468 CALL timeset(routinen, handle)
469
470 NULLIFY (sab_orb)
471
472 ! get basic quantities from the qs_env
473 CALL get_qs_env(qs_env, &
474 dft_control=dft_control, &
475 matrix_s=matrix_qs_s, &
476 matrix_ks=matrix_qs_ks, &
477 ks_env=ks_env, &
478 sab_orb=sab_orb)
479
480 nspin = dft_control%nspins
481
482 ! create matrix_ks in the QS env if necessary
483 IF (.NOT. ASSOCIATED(matrix_qs_ks)) THEN
484 CALL dbcsr_allocate_matrix_set(matrix_qs_ks, nspin)
485 DO ispin = 1, nspin
486 ALLOCATE (matrix_qs_ks(ispin)%matrix)
487 CALL dbcsr_create(matrix_qs_ks(ispin)%matrix, &
488 template=matrix_qs_s(1)%matrix)
489 CALL cp_dbcsr_alloc_block_from_nbl(matrix_qs_ks(ispin)%matrix, sab_orb)
490 CALL dbcsr_set(matrix_qs_ks(ispin)%matrix, 0.0_dp)
491 END DO
492 CALL set_ks_env(ks_env, matrix_ks=matrix_qs_ks)
493 END IF
494
495 ! copy to ALMO
496 DO ispin = 1, nspin
497 CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, matrix_ks(ispin), mat_distr_aos)
498 CALL dbcsr_filter(matrix_ks(ispin), eps_filter)
499 END DO
500
501 CALL timestop(handle)
502
503 END SUBROUTINE init_almo_ks_matrix_via_qs
504
505! **************************************************************************************************
506!> \brief Create MOs in the QS env to be able to return ALMOs to QS
507!> \param qs_env ...
508!> \param almo_scf_env ...
509!> \par History
510!> 2016.12 created [Yifei Shi]
511!> \author Yifei Shi
512! **************************************************************************************************
513 SUBROUTINE construct_qs_mos(qs_env, almo_scf_env)
514
515 TYPE(qs_environment_type), POINTER :: qs_env
516 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env
517
518 CHARACTER(len=*), PARAMETER :: routinen = 'construct_qs_mos'
519
520 INTEGER :: handle, ispin, ncol_fm, nrow_fm
521 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
522 TYPE(cp_fm_type) :: mo_fm_copy
523 TYPE(dft_control_type), POINTER :: dft_control
524 TYPE(mo_set_type), DIMENSION(:), POINTER :: mos
525 TYPE(qs_scf_env_type), POINTER :: scf_env
526
527 CALL timeset(routinen, handle)
528
529 ! create and init scf_env (this is necessary to return MOs to qs)
530 NULLIFY (mos, fm_struct_tmp, scf_env)
531 ALLOCATE (scf_env)
532 CALL scf_env_create(scf_env)
533
534 !CALL qs_scf_env_initialize(qs_env, scf_env)
535 CALL set_qs_env(qs_env, scf_env=scf_env)
536 CALL get_qs_env(qs_env, dft_control=dft_control, mos=mos)
537
538 CALL dbcsr_get_info(almo_scf_env%matrix_t(1), nfullrows_total=nrow_fm, nfullcols_total=ncol_fm)
539
540 ! allocate and init mo_set
541 DO ispin = 1, almo_scf_env%nspins
542 CALL dbcsr_get_info(almo_scf_env%matrix_t(ispin), nfullrows_total=nrow_fm, nfullcols_total=ncol_fm)
543
544 ! Currently only fm version of mo_set is usable.
545 ! First transform the matrix_t to fm version
546 ! Empty the containers to prevent memory leaks
547 CALL deallocate_mo_set(mos(ispin))
548
549 IF (almo_scf_env%nspins == 1) THEN
550 CALL allocate_mo_set(mo_set=mos(ispin), &
551 nao=nrow_fm, &
552 nmo=ncol_fm, &
553 nelectron=almo_scf_env%nelectrons_total, &
554 n_el_f=real(almo_scf_env%nelectrons_total, dp), &
555 maxocc=2.0_dp, &
556 flexible_electron_count=dft_control%relax_multiplicity)
557 ELSE IF (almo_scf_env%nspins == 2) THEN
558 CALL allocate_mo_set(mo_set=mos(ispin), &
559 nao=nrow_fm, &
560 nmo=ncol_fm, &
561 nelectron=sum(almo_scf_env%nocc_of_domain(:, ispin)), &
562 n_el_f=real(sum(almo_scf_env%nocc_of_domain(:, ispin)), dp), &
563 maxocc=1.0_dp, &
564 flexible_electron_count=dft_control%relax_multiplicity)
565 END IF
566
567 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=nrow_fm, ncol_global=ncol_fm, &
568 context=almo_scf_env%blacs_env, &
569 para_env=almo_scf_env%para_env)
570
571 CALL cp_fm_create(mo_fm_copy, fm_struct_tmp, name="t_orthogonal_converted_to_fm")
572 CALL cp_fm_struct_release(fm_struct_tmp)
573 !CALL copy_dbcsr_to_fm(almo_scf_env%matrix_t(ispin), mo_fm_copy)
574
575 CALL init_mo_set(mos(ispin), fm_ref=mo_fm_copy, name='fm_mo')
576
577 CALL cp_fm_release(mo_fm_copy)
578
579 END DO
580
581 CALL timestop(handle)
582
583 END SUBROUTINE construct_qs_mos
584
585! **************************************************************************************************
586!> \brief return density matrix to the qs_env
587!> \param qs_env ...
588!> \param matrix_p ...
589!> \param mat_distr_aos ...
590!> \par History
591!> 2011.05 created [Rustam Z Khaliullin]
592!> \author Rustam Z Khaliullin
593! **************************************************************************************************
594 SUBROUTINE almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
595 TYPE(qs_environment_type), POINTER :: qs_env
596 TYPE(dbcsr_type), DIMENSION(:) :: matrix_p
597 INTEGER, INTENT(IN) :: mat_distr_aos
598
599 CHARACTER(len=*), PARAMETER :: routinen = 'almo_dm_to_qs_env'
600
601 INTEGER :: handle, ispin, nspins
602 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
603 TYPE(qs_rho_type), POINTER :: rho
604
605 CALL timeset(routinen, handle)
606
607 NULLIFY (rho, rho_ao)
608 nspins = SIZE(matrix_p)
609 CALL get_qs_env(qs_env, rho=rho)
610 CALL qs_rho_get(rho, rho_ao=rho_ao)
611
612 ! set the new density matrix
613 DO ispin = 1, nspins
614 CALL matrix_almo_to_qs(matrix_p(ispin), &
615 rho_ao(ispin)%matrix, &
616 mat_distr_aos)
617 END DO
618 CALL qs_rho_update_rho(rho, qs_env=qs_env)
619 CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.true.)
620
621 CALL timestop(handle)
622
623 END SUBROUTINE almo_dm_to_qs_env
624
625! **************************************************************************************************
626!> \brief uses the ALMO density matrix
627!> to compute KS matrix (inside QS environment) and the new energy
628!> \param qs_env ...
629!> \param matrix_p ...
630!> \param energy_total ...
631!> \param mat_distr_aos ...
632!> \param smear ...
633!> \param kTS_sum ...
634!> \par History
635!> 2011.05 created [Rustam Z Khaliullin]
636!> 2018.09 smearing support [Ruben Staub]
637!> \author Rustam Z Khaliullin
638! **************************************************************************************************
639 SUBROUTINE almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos, smear, kTS_sum)
640 TYPE(qs_environment_type), POINTER :: qs_env
641 TYPE(dbcsr_type), DIMENSION(:) :: matrix_p
642 REAL(kind=dp) :: energy_total
643 INTEGER, INTENT(IN) :: mat_distr_aos
644 LOGICAL, INTENT(IN), OPTIONAL :: smear
645 REAL(kind=dp), INTENT(IN), OPTIONAL :: kts_sum
646
647 CHARACTER(len=*), PARAMETER :: routinen = 'almo_dm_to_qs_ks'
648
649 INTEGER :: handle
650 LOGICAL :: smearing
651 REAL(kind=dp) :: entropic_term
652 TYPE(qs_energy_type), POINTER :: energy
653
654 CALL timeset(routinen, handle)
655
656 IF (PRESENT(smear)) THEN
657 smearing = smear
658 ELSE
659 smearing = .false.
660 END IF
661
662 IF (PRESENT(kts_sum)) THEN
663 entropic_term = kts_sum
664 ELSE
665 entropic_term = 0.0_dp
666 END IF
667
668 NULLIFY (energy)
669 CALL get_qs_env(qs_env, energy=energy)
670 CALL almo_dm_to_qs_env(qs_env, matrix_p, mat_distr_aos)
671 CALL qs_ks_update_qs_env(qs_env, calculate_forces=.false., just_energy=.false., &
672 print_active=.true.)
673
674 !! Add electronic entropy contribution if smearing is requested
675 !! Previous QS entropy is replaced by the sum of the entropy for each spin
676 IF (smearing) THEN
677 energy%total = energy%total - energy%kTS + entropic_term
678 END IF
679
680 energy_total = energy%total
681
682 CALL timestop(handle)
683
684 END SUBROUTINE almo_dm_to_qs_ks
685
686! **************************************************************************************************
687!> \brief uses the ALMO density matrix
688!> to compute ALMO KS matrix and the new energy
689!> \param qs_env ...
690!> \param matrix_p ...
691!> \param matrix_ks ...
692!> \param energy_total ...
693!> \param eps_filter ...
694!> \param mat_distr_aos ...
695!> \param smear ...
696!> \param kTS_sum ...
697!> \par History
698!> 2011.05 created [Rustam Z Khaliullin]
699!> 2018.09 smearing support [Ruben Staub]
700!> \author Rustam Z Khaliullin
701! **************************************************************************************************
702 SUBROUTINE almo_dm_to_almo_ks(qs_env, matrix_p, matrix_ks, energy_total, eps_filter, &
703 mat_distr_aos, smear, kTS_sum)
704
705 TYPE(qs_environment_type), POINTER :: qs_env
706 TYPE(dbcsr_type), DIMENSION(:) :: matrix_p, matrix_ks
707 REAL(kind=dp) :: energy_total, eps_filter
708 INTEGER, INTENT(IN) :: mat_distr_aos
709 LOGICAL, INTENT(IN), OPTIONAL :: smear
710 REAL(kind=dp), INTENT(IN), OPTIONAL :: kts_sum
711
712 CHARACTER(len=*), PARAMETER :: routinen = 'almo_dm_to_almo_ks'
713
714 INTEGER :: handle, ispin, nspins
715 LOGICAL :: smearing
716 REAL(kind=dp) :: entropic_term
717 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_qs_ks
718
719 CALL timeset(routinen, handle)
720
721 IF (PRESENT(smear)) THEN
722 smearing = smear
723 ELSE
724 smearing = .false.
725 END IF
726
727 IF (PRESENT(kts_sum)) THEN
728 entropic_term = kts_sum
729 ELSE
730 entropic_term = 0.0_dp
731 END IF
732
733 ! update KS matrix in the QS env
734 CALL almo_dm_to_qs_ks(qs_env, matrix_p, energy_total, mat_distr_aos, &
735 smear=smearing, &
736 kts_sum=entropic_term)
737
738 nspins = SIZE(matrix_ks)
739
740 ! get KS matrix from the QS env and convert to the ALMO format
741 CALL get_qs_env(qs_env, matrix_ks=matrix_qs_ks)
742 DO ispin = 1, nspins
743 CALL matrix_qs_to_almo(matrix_qs_ks(ispin)%matrix, matrix_ks(ispin), mat_distr_aos)
744 CALL dbcsr_filter(matrix_ks(ispin), eps_filter)
745 END DO
746
747 CALL timestop(handle)
748
749 END SUBROUTINE almo_dm_to_almo_ks
750
751! **************************************************************************************************
752!> \brief update qs_env total energy
753!> \param qs_env ...
754!> \param energy ...
755!> \param energy_singles_corr ...
756!> \par History
757!> 2013.03 created [Rustam Z Khaliullin]
758!> \author Rustam Z Khaliullin
759! **************************************************************************************************
760 SUBROUTINE almo_scf_update_ks_energy(qs_env, energy, energy_singles_corr)
761
762 TYPE(qs_environment_type), POINTER :: qs_env
763 REAL(kind=dp), INTENT(IN), OPTIONAL :: energy, energy_singles_corr
764
765 TYPE(qs_energy_type), POINTER :: qs_energy
766
767 CALL get_qs_env(qs_env, energy=qs_energy)
768
769 IF (PRESENT(energy_singles_corr)) THEN
770 qs_energy%singles_corr = energy_singles_corr
771 ELSE
772 qs_energy%singles_corr = 0.0_dp
773 END IF
774
775 IF (PRESENT(energy)) THEN
776 qs_energy%total = energy
777 END IF
778
779 qs_energy%total = qs_energy%total + qs_energy%singles_corr
780
781 END SUBROUTINE almo_scf_update_ks_energy
782
783! **************************************************************************************************
784!> \brief Creates the matrix that imposes absolute locality on MOs
785!> \param qs_env ...
786!> \param almo_scf_env ...
787!> \par History
788!> 2011.11 created [Rustam Z. Khaliullin]
789!> \author Rustam Z. Khaliullin
790! **************************************************************************************************
791 SUBROUTINE almo_scf_construct_quencher(qs_env, almo_scf_env)
792
793 TYPE(qs_environment_type), POINTER :: qs_env
794 TYPE(almo_scf_env_type), INTENT(INOUT) :: almo_scf_env
795
796 CHARACTER(len=*), PARAMETER :: routinen = 'almo_scf_construct_quencher'
797
798 CHARACTER :: sym
799 INTEGER :: col, contact_atom_1, contact_atom_2, domain_col, domain_map_local_entries, &
800 domain_row, global_entries, global_list_length, grid1, groupid, handle, hold, iatom, &
801 iatom2, iblock_col, iblock_row, idomain, idomain2, ientry, igrid, ineig, ineighbor, &
802 inode, inode2, ipair, ispin, jatom, jatom2, jdomain2, local_list_length, &
803 max_domain_neighbors, max_neig, mynode, nblkcols_tot, nblkrows_tot, nblks, ndomains, &
804 neig_temp, nnode2, nnodes, row, unit_nr
805 INTEGER, ALLOCATABLE, DIMENSION(:) :: current_number_neighbors, domain_entries_cpu, &
806 domain_map_global, domain_map_local, first_atom_of_molecule, global_list, &
807 last_atom_of_molecule, list_length_cpu, list_offset_cpu, local_list, offset_for_cpu
808 INTEGER, ALLOCATABLE, DIMENSION(:, :) :: domain_grid, domain_neighbor_list, &
809 domain_neighbor_list_excessive
810 INTEGER, DIMENSION(:), POINTER :: col_blk_size, row_blk_size
811 LOGICAL :: already_listed, block_active, &
812 delayed_increment, found, &
813 max_neig_fails, tr
814 REAL(kind=dp) :: contact1_radius, contact2_radius, &
815 distance, distance_squared, overlap, &
816 r0, r1, s0, s1, trial_distance_squared
817 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: new_block
818 REAL(kind=dp), DIMENSION(3) :: rab
819 REAL(kind=dp), DIMENSION(:, :), POINTER :: p_old_block
820 TYPE(cell_type), POINTER :: cell
821 TYPE(cp_logger_type), POINTER :: logger
822 TYPE(dbcsr_distribution_type) :: dist
823 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
824 TYPE(dbcsr_type) :: matrix_s_sym
825 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
826 TYPE(mp_comm_type) :: group
828 DIMENSION(:), POINTER :: nl_iterator, nl_iterator2
829 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
830 POINTER :: sab_almo
831 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
832
833 CALL timeset(routinen, handle)
834
835 ! get a useful output_unit
836 logger => cp_get_default_logger()
837 IF (logger%para_env%is_source()) THEN
838 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
839 ELSE
840 unit_nr = -1
841 END IF
842
843 ndomains = almo_scf_env%ndomains
844
845 CALL get_qs_env(qs_env=qs_env, &
846 particle_set=particle_set, &
847 molecule_set=molecule_set, &
848 cell=cell, &
849 matrix_s=matrix_s, &
850 sab_almo=sab_almo)
851
852 ! if we are dealing with molecules get info about them
853 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular .OR. &
854 almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
855 ALLOCATE (first_atom_of_molecule(almo_scf_env%nmolecules))
856 ALLOCATE (last_atom_of_molecule(almo_scf_env%nmolecules))
857 CALL get_molecule_set_info(molecule_set, &
858 mol_to_first_atom=first_atom_of_molecule, &
859 mol_to_last_atom=last_atom_of_molecule)
860 END IF
861
862 ! create a symmetrized copy of the ao overlap
863 CALL dbcsr_create(matrix_s_sym, &
864 template=almo_scf_env%matrix_s(1), &
865 matrix_type=dbcsr_type_no_symmetry)
866 CALL dbcsr_get_info(almo_scf_env%matrix_s(1), &
867 matrix_type=sym)
868 IF (sym == dbcsr_type_no_symmetry) THEN
869 CALL dbcsr_copy(matrix_s_sym, almo_scf_env%matrix_s(1))
870 ELSE
871 CALL dbcsr_desymmetrize(almo_scf_env%matrix_s(1), &
872 matrix_s_sym)
873 END IF
874
875 ALLOCATE (almo_scf_env%quench_t(almo_scf_env%nspins))
876 ALLOCATE (almo_scf_env%domain_map(almo_scf_env%nspins))
877
878 DO ispin = 1, almo_scf_env%nspins
879
880 ! create the sparsity template for the occupied orbitals
881 CALL matrix_almo_create(matrix_new=almo_scf_env%quench_t(ispin), &
882 matrix_qs=matrix_s(1)%matrix, &
883 almo_scf_env=almo_scf_env, &
884 name_new="T_QUENCHER", &
886 symmetry_new=dbcsr_type_no_symmetry, &
887 spin_key=ispin, &
888 init_domains=.false.)
889
890 ! initialize distance quencher
891 CALL dbcsr_work_create(almo_scf_env%quench_t(ispin), work_mutable=.true.)
892 CALL dbcsr_get_info(almo_scf_env%quench_t(ispin), distribution=dist, &
893 nblkrows_total=nblkrows_tot, nblkcols_total=nblkcols_tot)
894 CALL dbcsr_distribution_get(dist, numnodes=nnodes, group=groupid, mynode=mynode)
895 CALL group%set_handle(groupid)
896
897 ! create global atom neighbor list from the local lists
898 ! first, calculate number of local pairs
899 local_list_length = 0
900 CALL neighbor_list_iterator_create(nl_iterator, sab_almo)
901 DO WHILE (neighbor_list_iterate(nl_iterator) == 0)
902 ! nnode - total number of neighbors for iatom
903 ! inode - current neighbor count
904 CALL get_iterator_info(nl_iterator, &
905 iatom=iatom2, jatom=jatom2, inode=inode2, nnode=nnode2)
906 IF (inode2 == 1) THEN
907 local_list_length = local_list_length + nnode2
908 END IF
909 END DO
910 CALL neighbor_list_iterator_release(nl_iterator)
911
912 ! second, extract the local list to an array
913 ALLOCATE (local_list(2*local_list_length))
914 local_list(:) = 0
915 local_list_length = 0
916 CALL neighbor_list_iterator_create(nl_iterator2, sab_almo)
917 DO WHILE (neighbor_list_iterate(nl_iterator2) == 0)
918 CALL get_iterator_info(nl_iterator2, &
919 iatom=iatom2, jatom=jatom2)
920 local_list(2*local_list_length + 1) = iatom2
921 local_list(2*local_list_length + 2) = jatom2
922 local_list_length = local_list_length + 1
923 END DO ! end loop over pairs of atoms
924 CALL neighbor_list_iterator_release(nl_iterator2)
925
926 ! third, communicate local length to the other nodes
927 ALLOCATE (list_length_cpu(nnodes), list_offset_cpu(nnodes))
928 CALL group%allgather(2*local_list_length, list_length_cpu)
929
930 ! fourth, create a global list
931 list_offset_cpu(1) = 0
932 DO inode = 2, nnodes
933 list_offset_cpu(inode) = list_offset_cpu(inode - 1) + &
934 list_length_cpu(inode - 1)
935 END DO
936 global_list_length = list_offset_cpu(nnodes) + list_length_cpu(nnodes)
937
938 ! fifth, communicate all list data
939 ALLOCATE (global_list(global_list_length))
940 CALL group%allgatherv(local_list, global_list, &
941 list_length_cpu, list_offset_cpu)
942 DEALLOCATE (list_length_cpu, list_offset_cpu)
943 DEALLOCATE (local_list)
944
945 ! calculate maximum number of atoms surrounding the domain
946 ALLOCATE (current_number_neighbors(almo_scf_env%ndomains))
947 current_number_neighbors(:) = 0
948 global_list_length = global_list_length/2
949 DO ipair = 1, global_list_length
950 iatom2 = global_list(2*(ipair - 1) + 1)
951 jatom2 = global_list(2*(ipair - 1) + 2)
952 idomain2 = almo_scf_env%domain_index_of_atom(iatom2)
953 jdomain2 = almo_scf_env%domain_index_of_atom(jatom2)
954 ! add to the list
955 current_number_neighbors(idomain2) = current_number_neighbors(idomain2) + 1
956 ! add j,i with i,j
957 IF (idomain2 /= jdomain2) THEN
958 current_number_neighbors(jdomain2) = current_number_neighbors(jdomain2) + 1
959 END IF
960 END DO
961 max_domain_neighbors = maxval(current_number_neighbors)
962
963 ! use the global atom neighbor list to create a global domain neighbor list
964 ALLOCATE (domain_neighbor_list_excessive(ndomains, max_domain_neighbors))
965 current_number_neighbors(:) = 1
966 DO ipair = 1, ndomains
967 domain_neighbor_list_excessive(ipair, 1) = ipair
968 END DO
969 DO ipair = 1, global_list_length
970 iatom2 = global_list(2*(ipair - 1) + 1)
971 jatom2 = global_list(2*(ipair - 1) + 2)
972 idomain2 = almo_scf_env%domain_index_of_atom(iatom2)
973 jdomain2 = almo_scf_env%domain_index_of_atom(jatom2)
974 already_listed = .false.
975 DO ineighbor = 1, current_number_neighbors(idomain2)
976 IF (domain_neighbor_list_excessive(idomain2, ineighbor) == jdomain2) THEN
977 already_listed = .true.
978 EXIT
979 END IF
980 END DO
981 IF (.NOT. already_listed) THEN
982 ! add to the list
983 current_number_neighbors(idomain2) = current_number_neighbors(idomain2) + 1
984 domain_neighbor_list_excessive(idomain2, current_number_neighbors(idomain2)) = jdomain2
985 ! add j,i with i,j
986 IF (idomain2 /= jdomain2) THEN
987 current_number_neighbors(jdomain2) = current_number_neighbors(jdomain2) + 1
988 domain_neighbor_list_excessive(jdomain2, current_number_neighbors(jdomain2)) = idomain2
989 END IF
990 END IF
991 END DO ! end loop over pairs of atoms
992 DEALLOCATE (global_list)
993
994 max_domain_neighbors = maxval(current_number_neighbors)
995 ALLOCATE (domain_neighbor_list(ndomains, max_domain_neighbors))
996 domain_neighbor_list(:, :) = 0
997 domain_neighbor_list(:, :) = domain_neighbor_list_excessive(:, 1:max_domain_neighbors)
998 DEALLOCATE (domain_neighbor_list_excessive)
999
1000 ALLOCATE (almo_scf_env%domain_map(ispin)%index1(ndomains))
1001 ALLOCATE (almo_scf_env%domain_map(ispin)%pairs(max_domain_neighbors*ndomains, 2))
1002 almo_scf_env%domain_map(ispin)%pairs(:, :) = 0
1003 almo_scf_env%domain_map(ispin)%index1(:) = 0
1004 domain_map_local_entries = 0
1005
1006 ! RZK-warning intermediate [0,1] quencher values are ill-defined
1007 ! for molecules (not continuous and conceptually inadequate)
1008
1009 CALL dbcsr_get_info(almo_scf_env%quench_t(ispin), &
1010 row_blk_size=row_blk_size, col_blk_size=col_blk_size)
1011 ! O(N) loop over domain pairs
1012 DO row = 1, nblkrows_tot
1013 DO col = 1, current_number_neighbors(row)
1014 tr = .false.
1015 iblock_row = row
1016 iblock_col = domain_neighbor_list(row, col)
1017 CALL dbcsr_get_stored_coordinates(almo_scf_env%quench_t(ispin), &
1018 iblock_row, iblock_col, hold)
1019
1020 IF (hold == mynode) THEN
1021
1022 ! Translate indices of distribution blocks to indices of domain blocks
1023 ! Rows are AOs
1024 domain_row = almo_scf_env%domain_index_of_ao_block(iblock_row)
1025 ! Columns are electrons (i.e. MOs)
1026 domain_col = almo_scf_env%domain_index_of_mo_block(iblock_col)
1027
1028 SELECT CASE (almo_scf_env%constraint_type)
1030
1031 block_active = .false.
1032 ! type of electron groups
1033 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
1034
1035 ! type of ao domains
1036 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1037
1038 ! ao domains are molecular / electron groups are molecular
1039 IF (domain_row == domain_col) THEN
1040 block_active = .true.
1041 END IF
1042
1043 ELSE ! ao domains are atomic
1044
1045 ! ao domains are atomic / electron groups are molecular
1046 cpabort("Illegal: atomic domains and molecular groups")
1047
1048 END IF
1049
1050 ELSE ! electron groups are atomic
1051
1052 ! type of ao domains
1053 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1054
1055 ! ao domains are molecular / electron groups are atomic
1056 cpabort("Illegal: molecular domains and atomic groups")
1057
1058 ELSE
1059
1060 ! ao domains are atomic / electron groups are atomic
1061 IF (domain_row == domain_col) THEN
1062 block_active = .true.
1063 END IF
1064
1065 END IF
1066
1067 END IF ! end type of electron groups
1068
1069 IF (block_active) THEN
1070
1071 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1072 new_block(:, :) = 1.0_dp
1073 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1074 DEALLOCATE (new_block)
1075
1076 IF (domain_map_local_entries >= max_domain_neighbors*almo_scf_env%ndomains) THEN
1077 cpabort("weird... max_domain_neighbors is exceeded")
1078 END IF
1079 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1080 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1081 domain_map_local_entries = domain_map_local_entries + 1
1082
1083 END IF
1084
1086
1087 ! type of electron groups
1088 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
1089
1090 ! type of ao domains
1091 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1092
1093 ! ao domains are molecular / electron groups are molecular
1094
1095 ! compute the maximum overlap between the atoms of the two molecules
1096 CALL dbcsr_get_block_p(matrix_s_sym, iblock_row, iblock_col, p_old_block, found)
1097 IF (found) THEN
1098 overlap = maxval(abs(p_old_block))
1099 ELSE
1100 overlap = 0.0_dp
1101 END IF
1102
1103 ELSE ! ao domains are atomic
1104
1105 ! ao domains are atomic / electron groups are molecular
1106 ! overlap_between_atom_and_molecule(atom=domain_row,molecule=domain_col)
1107 cpabort("atomic domains and molecular groups - NYI")
1108
1109 END IF
1110
1111 ELSE ! electron groups are atomic
1112
1113 ! type of ao domains
1114 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1115
1116 ! ao domains are molecular / electron groups are atomic
1117 ! overlap_between_atom_and_molecule(atom=domain_col,molecule=domain_row)
1118 cpabort("molecular domains and atomic groups - NYI")
1119
1120 ELSE
1121
1122 ! ao domains are atomic / electron groups are atomic
1123 ! compute max overlap between atoms: domain_row and domain_col
1124 CALL dbcsr_get_block_p(matrix_s_sym, iblock_row, iblock_col, p_old_block, found)
1125 IF (found) THEN
1126 overlap = maxval(abs(p_old_block))
1127 ELSE
1128 overlap = 0.0_dp
1129 END IF
1130
1131 END IF
1132
1133 END IF ! end type of electron groups
1134
1135 s0 = -log10(abs(almo_scf_env%quencher_s0))
1136 s1 = -log10(abs(almo_scf_env%quencher_s1))
1137 IF (overlap == 0.0_dp) THEN
1138 overlap = -log10(abs(almo_scf_env%eps_filter)) + 100.0_dp
1139 ELSE
1140 overlap = -log10(overlap)
1141 END IF
1142 IF (s0 < 0.0_dp) THEN
1143 cpabort("S0 is less than zero")
1144 END IF
1145 IF (s1 <= 0.0_dp) THEN
1146 cpabort("S1 is less than or equal to zero")
1147 END IF
1148 IF (s0 >= s1) THEN
1149 cpabort("S0 is greater than or equal to S1")
1150 END IF
1151
1152 ! Fill in non-zero blocks if AOs are close to the electron center
1153 IF (overlap < s1) THEN
1154 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1155 IF (overlap <= s0) THEN
1156 new_block(:, :) = 1.0_dp
1157 ELSE
1158 new_block(:, :) = 1.0_dp/(1.0_dp + exp(-(s0 - s1)/(s0 - overlap) - (s0 - s1)/(overlap - s1)))
1159 END IF
1160
1161 IF (abs(new_block(1, 1)) > abs(almo_scf_env%eps_filter)) THEN
1162 IF (domain_map_local_entries >= max_domain_neighbors*almo_scf_env%ndomains) THEN
1163 cpabort("weird... max_domain_neighbors is exceeded")
1164 END IF
1165 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1166 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1167 domain_map_local_entries = domain_map_local_entries + 1
1168 END IF
1169
1170 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1171 DEALLOCATE (new_block)
1172
1173 END IF
1174
1176
1177 ! type of electron groups
1178 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular) THEN
1179
1180 ! type of ao domains
1181 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1182
1183 ! ao domains are molecular / electron groups are molecular
1184
1185 ! compute distance between molecules: domain_row and domain_col
1186 ! distance between molecules is defined as the smallest
1187 ! distance among all atom pairs
1188 IF (domain_row == domain_col) THEN
1189 distance = 0.0_dp
1190 contact_atom_1 = first_atom_of_molecule(domain_row)
1191 contact_atom_2 = first_atom_of_molecule(domain_col)
1192 ELSE
1193 distance_squared = 1.0e+100_dp
1194 contact_atom_1 = -1
1195 contact_atom_2 = -1
1196 DO iatom = first_atom_of_molecule(domain_row), last_atom_of_molecule(domain_row)
1197 DO jatom = first_atom_of_molecule(domain_col), last_atom_of_molecule(domain_col)
1198 rab(:) = pbc(particle_set(iatom)%r(:), particle_set(jatom)%r(:), cell)
1199 trial_distance_squared = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
1200 IF (trial_distance_squared < distance_squared) THEN
1201 distance_squared = trial_distance_squared
1202 contact_atom_1 = iatom
1203 contact_atom_2 = jatom
1204 END IF
1205 END DO ! jatom
1206 END DO ! iatom
1207 cpassert(contact_atom_1 > 0)
1208 distance = sqrt(distance_squared)
1209 END IF
1210
1211 ELSE ! ao domains are atomic
1212
1213 ! ao domains are atomic / electron groups are molecular
1214 !distance_between_atom_and_molecule(atom=domain_row,molecule=domain_col)
1215 cpabort("atomic domains and molecular groups - NYI")
1216
1217 END IF
1218
1219 ELSE ! electron groups are atomic
1220
1221 ! type of ao domains
1222 IF (almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1223
1224 ! ao domains are molecular / electron groups are atomic
1225 !distance_between_atom_and_molecule(atom=domain_col,molecule=domain_row)
1226 cpabort("molecular domains and atomic groups - NYI")
1227
1228 ELSE
1229
1230 ! ao domains are atomic / electron groups are atomic
1231 ! compute distance between atoms: domain_row and domain_col
1232 rab(:) = pbc(particle_set(domain_row)%r(:), particle_set(domain_col)%r(:), cell)
1233 distance = sqrt(rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3))
1234 contact_atom_1 = domain_row
1235 contact_atom_2 = domain_col
1236
1237 END IF
1238
1239 END IF ! end type of electron groups
1240
1241 ! get atomic radii to compute distance cutoff threshold
1242 IF (almo_scf_env%quencher_radius_type == do_bondparm_covalent) THEN
1243 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_1)%atomic_kind, &
1244 rcov=contact1_radius)
1245 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_2)%atomic_kind, &
1246 rcov=contact2_radius)
1247 ELSE IF (almo_scf_env%quencher_radius_type == do_bondparm_vdw) THEN
1248 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_1)%atomic_kind, &
1249 rvdw=contact1_radius)
1250 CALL get_atomic_kind(atomic_kind=particle_set(contact_atom_2)%atomic_kind, &
1251 rvdw=contact2_radius)
1252 ELSE
1253 cpabort("Illegal quencher_radius_type")
1254 END IF
1255 contact1_radius = cp_unit_to_cp2k(contact1_radius, "angstrom")
1256 contact2_radius = cp_unit_to_cp2k(contact2_radius, "angstrom")
1257
1258 !RZK-warning the procedure is faulty for molecules:
1259 ! the closest contacts should be found using
1260 ! the element specific radii
1261
1262 ! compute inner and outer cutoff radii
1263 r0 = almo_scf_env%quencher_r0_factor*(contact1_radius + contact2_radius)
1264 !+almo_scf_env%quencher_r0_shift
1265 r1 = almo_scf_env%quencher_r1_factor*(contact1_radius + contact2_radius)
1266 !+almo_scf_env%quencher_r1_shift
1267
1268 IF (r0 < 0.0_dp) THEN
1269 cpabort("R0 is less than zero")
1270 END IF
1271 IF (r1 <= 0.0_dp) THEN
1272 cpabort("R1 is less than or equal to zero")
1273 END IF
1274 IF (r0 > r1) THEN
1275 cpabort("R0 is greater than or equal to R1")
1276 END IF
1277
1278 ! Fill in non-zero blocks if AOs are close to the electron center
1279 IF (distance < r1) THEN
1280 ALLOCATE (new_block(row_blk_size(iblock_row), col_blk_size(iblock_col)))
1281 IF (distance <= r0) THEN
1282 new_block(:, :) = 1.0_dp
1283 ELSE
1284 ! remove the intermediate values from the quencher temporarily
1285 cpabort("distance > r0 not yet validated") ! Unexplained in https://github.com/cp2k/cp2k/pull/5345
1286 new_block(:, :) = 1.0_dp/(1.0_dp + exp((r1 - r0)/(r0 - distance) + (r1 - r0)/(r1 - distance)))
1287 END IF
1288
1289 IF (abs(new_block(1, 1)) > abs(almo_scf_env%eps_filter)) THEN
1290 IF (domain_map_local_entries >= max_domain_neighbors*almo_scf_env%ndomains) THEN
1291 cpabort("weird... max_domain_neighbors is exceeded")
1292 END IF
1293 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 1) = iblock_row
1294 almo_scf_env%domain_map(ispin)%pairs(domain_map_local_entries + 1, 2) = iblock_col
1295 domain_map_local_entries = domain_map_local_entries + 1
1296 END IF
1297
1298 CALL dbcsr_put_block(almo_scf_env%quench_t(ispin), iblock_row, iblock_col, new_block)
1299 DEALLOCATE (new_block)
1300 END IF
1301
1302 CASE DEFAULT
1303 cpabort("Illegal constraint type")
1304 END SELECT
1305
1306 END IF ! mynode
1307
1308 END DO
1309 END DO ! end O(N) loop over pairs
1310
1311 DEALLOCATE (domain_neighbor_list)
1312 DEALLOCATE (current_number_neighbors)
1313
1314 CALL dbcsr_finalize(almo_scf_env%quench_t(ispin))
1315
1316 CALL dbcsr_filter(almo_scf_env%quench_t(ispin), &
1317 almo_scf_env%eps_filter)
1318
1319 ! check that both domain_map and quench_t have the same number of entries
1320 nblks = dbcsr_get_num_blocks(almo_scf_env%quench_t(ispin))
1321 IF (nblks /= domain_map_local_entries) THEN
1322 cpabort("number of blocks is wrong")
1323 END IF
1324
1325 ! first, communicate map sizes on the other nodes
1326 ALLOCATE (domain_entries_cpu(nnodes), offset_for_cpu(nnodes))
1327 CALL group%allgather(2*domain_map_local_entries, domain_entries_cpu)
1328
1329 ! second, create
1330 offset_for_cpu(1) = 0
1331 DO inode = 2, nnodes
1332 offset_for_cpu(inode) = offset_for_cpu(inode - 1) + &
1333 domain_entries_cpu(inode - 1)
1334 END DO
1335 global_entries = offset_for_cpu(nnodes) + domain_entries_cpu(nnodes)
1336
1337 ! communicate all entries
1338 ALLOCATE (domain_map_global(global_entries))
1339 ALLOCATE (domain_map_local(2*domain_map_local_entries))
1340 DO ientry = 1, domain_map_local_entries
1341 domain_map_local(2*(ientry - 1) + 1) = almo_scf_env%domain_map(ispin)%pairs(ientry, 1)
1342 domain_map_local(2*ientry) = almo_scf_env%domain_map(ispin)%pairs(ientry, 2)
1343 END DO
1344 CALL group%allgatherv(domain_map_local, domain_map_global, &
1345 domain_entries_cpu, offset_for_cpu)
1346 DEALLOCATE (domain_entries_cpu, offset_for_cpu)
1347 DEALLOCATE (domain_map_local)
1348
1349 DEALLOCATE (almo_scf_env%domain_map(ispin)%index1)
1350 DEALLOCATE (almo_scf_env%domain_map(ispin)%pairs)
1351 ALLOCATE (almo_scf_env%domain_map(ispin)%index1(ndomains))
1352 ALLOCATE (almo_scf_env%domain_map(ispin)%pairs(global_entries/2, 2))
1353 almo_scf_env%domain_map(ispin)%pairs(:, :) = 0
1354 almo_scf_env%domain_map(ispin)%index1(:) = 0
1355
1356 ! unpack the received data into a local variable
1357 ! since we do not know the maximum global number of neighbors
1358 ! try one. if fails increase the maximum number and try again
1359 ! until it succeeds
1360 max_neig = max_domain_neighbors
1361 max_neig_fails = .true.
1362 max_neig_loop: DO WHILE (max_neig_fails)
1363 ALLOCATE (domain_grid(almo_scf_env%ndomains, 0:max_neig))
1364 domain_grid(:, :) = 0
1365 ! init the number of collected neighbors
1366 domain_grid(:, 0) = 1
1367 ! loop over the records
1368 global_entries = global_entries/2
1369 DO ientry = 1, global_entries
1370 ! get the center
1371 grid1 = domain_map_global(2*ientry)
1372 ! get the neighbor
1373 ineig = domain_map_global(2*(ientry - 1) + 1)
1374 ! check boundaries
1375 IF (domain_grid(grid1, 0) > max_neig) THEN
1376 ! this neighbor will overstep the boundaries
1377 ! stop the trial and increase the max number of neighbors
1378 DEALLOCATE (domain_grid)
1379 max_neig = max_neig*2
1380 cycle max_neig_loop
1381 END IF
1382 ! for the current center loop over the collected neighbors
1383 ! to insert the current record in a numerical order
1384 delayed_increment = .false.
1385 DO igrid = 1, domain_grid(grid1, 0)
1386 ! compare the current neighbor with that already in the 'book'
1387 IF (ineig < domain_grid(grid1, igrid)) THEN
1388 ! if this one is smaller then insert it here and pick up the one
1389 ! from the book to continue inserting
1390 neig_temp = ineig
1391 ineig = domain_grid(grid1, igrid)
1392 domain_grid(grid1, igrid) = neig_temp
1393 ELSE
1394 IF (domain_grid(grid1, igrid) == 0) THEN
1395 ! got the empty slot now - insert the record
1396 domain_grid(grid1, igrid) = ineig
1397 ! increase the record counter but do it outside the loop
1398 delayed_increment = .true.
1399 END IF
1400 END IF
1401 END DO
1402 IF (delayed_increment) THEN
1403 domain_grid(grid1, 0) = domain_grid(grid1, 0) + 1
1404 ELSE
1405 ! should not be here - all records must be inserted
1406 cpabort("all records must be inserted")
1407 END IF
1408 END DO
1409 max_neig_fails = .false.
1410 END DO max_neig_loop
1411 DEALLOCATE (domain_map_global)
1412
1413 ientry = 1
1414 DO idomain = 1, almo_scf_env%ndomains
1415 DO ineig = 1, domain_grid(idomain, 0) - 1
1416 almo_scf_env%domain_map(ispin)%pairs(ientry, 1) = domain_grid(idomain, ineig)
1417 almo_scf_env%domain_map(ispin)%pairs(ientry, 2) = idomain
1418 ientry = ientry + 1
1419 END DO
1420 almo_scf_env%domain_map(ispin)%index1(idomain) = ientry
1421 END DO
1422 DEALLOCATE (domain_grid)
1423
1424 END DO ! ispin
1425 IF (almo_scf_env%nspins == 2) THEN
1426 CALL dbcsr_copy(almo_scf_env%quench_t(2), &
1427 almo_scf_env%quench_t(1))
1428 almo_scf_env%domain_map(2)%pairs(:, :) = &
1429 almo_scf_env%domain_map(1)%pairs(:, :)
1430 almo_scf_env%domain_map(2)%index1(:) = &
1431 almo_scf_env%domain_map(1)%index1(:)
1432 END IF
1433
1434 CALL dbcsr_release(matrix_s_sym)
1435
1436 IF (almo_scf_env%domain_layout_mos == almo_domain_layout_molecular .OR. &
1437 almo_scf_env%domain_layout_aos == almo_domain_layout_molecular) THEN
1438 DEALLOCATE (first_atom_of_molecule)
1439 DEALLOCATE (last_atom_of_molecule)
1440 END IF
1441
1442 CALL timestop(handle)
1443
1444 END SUBROUTINE almo_scf_construct_quencher
1445
1446! *****************************************************************************
1447!> \brief Compute matrix W (energy-weighted density matrix) that is needed
1448!> for the evaluation of forces
1449!> \param matrix_w ...
1450!> \param almo_scf_env ...
1451!> \par History
1452!> 2015.03 created [Rustam Z. Khaliullin]
1453!> \author Rustam Z. Khaliullin
1454! **************************************************************************************************
1455 SUBROUTINE calculate_w_matrix_almo(matrix_w, almo_scf_env)
1456 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_w
1457 TYPE(almo_scf_env_type) :: almo_scf_env
1458
1459 CHARACTER(len=*), PARAMETER :: routinen = 'calculate_w_matrix_almo'
1460
1461 INTEGER :: handle, ispin
1462 REAL(kind=dp) :: scaling
1463 TYPE(dbcsr_type) :: tmp_nn1, tmp_no1, tmp_oo1, tmp_oo2
1464
1465 CALL timeset(routinen, handle)
1466
1467 IF (almo_scf_env%nspins == 1) THEN
1468 scaling = 2.0_dp
1469 ELSE
1470 scaling = 1.0_dp
1471 END IF
1472
1473 DO ispin = 1, almo_scf_env%nspins
1474
1475 CALL dbcsr_create(tmp_nn1, template=almo_scf_env%matrix_s(1), &
1476 matrix_type=dbcsr_type_no_symmetry)
1477 CALL dbcsr_create(tmp_no1, template=almo_scf_env%matrix_t(ispin), &
1478 matrix_type=dbcsr_type_no_symmetry)
1479 CALL dbcsr_create(tmp_oo1, template=almo_scf_env%matrix_sigma_inv(ispin), &
1480 matrix_type=dbcsr_type_no_symmetry)
1481 CALL dbcsr_create(tmp_oo2, template=almo_scf_env%matrix_sigma_inv(ispin), &
1482 matrix_type=dbcsr_type_no_symmetry)
1483
1484 CALL dbcsr_copy(tmp_nn1, almo_scf_env%matrix_ks(ispin))
1485 ! 1. TMP_NO1=F.T
1486 CALL dbcsr_multiply("N", "N", scaling, tmp_nn1, almo_scf_env%matrix_t(ispin), &
1487 0.0_dp, tmp_no1, filter_eps=almo_scf_env%eps_filter)
1488 ! 2. TMP_OO1=T^(tr).TMP_NO1=T^(tr).(FT)
1489 CALL dbcsr_multiply("T", "N", 1.0_dp, almo_scf_env%matrix_t(ispin), tmp_no1, &
1490 0.0_dp, tmp_oo1, filter_eps=almo_scf_env%eps_filter)
1491 ! 3. TMP_OO2=TMP_OO1.siginv=(T^(tr)FT).siginv
1492 CALL dbcsr_multiply("N", "N", 1.0_dp, tmp_oo1, almo_scf_env%matrix_sigma_inv(ispin), &
1493 0.0_dp, tmp_oo2, filter_eps=almo_scf_env%eps_filter)
1494 ! 4. TMP_OO1=siginv.TMP_OO2=siginv.(T^(tr)FTsiginv)
1495 CALL dbcsr_multiply("N", "N", 1.0_dp, almo_scf_env%matrix_sigma_inv(ispin), tmp_oo2, &
1496 0.0_dp, tmp_oo1, filter_eps=almo_scf_env%eps_filter)
1497 ! 5. TMP_NO1=T.TMP_OO1.=T.(siginvT^(tr)FTsiginv)
1498 CALL dbcsr_multiply("N", "N", 1.0_dp, almo_scf_env%matrix_t(ispin), tmp_oo1, &
1499 0.0_dp, tmp_no1, filter_eps=almo_scf_env%eps_filter)
1500 ! 6. TMP_NN1=TMP_NO1.T^(tr)=(TsiginvT^(tr)FTsiginv).T^(tr)=RFR
1501 CALL dbcsr_multiply("N", "T", 1.0_dp, tmp_no1, almo_scf_env%matrix_t(ispin), &
1502 0.0_dp, tmp_nn1, filter_eps=almo_scf_env%eps_filter)
1503 CALL matrix_almo_to_qs(tmp_nn1, matrix_w(ispin)%matrix, almo_scf_env%mat_distr_aos)
1504
1505 CALL dbcsr_release(tmp_nn1)
1506 CALL dbcsr_release(tmp_no1)
1507 CALL dbcsr_release(tmp_oo1)
1508 CALL dbcsr_release(tmp_oo2)
1509
1510 END DO
1511
1512 CALL timestop(handle)
1513
1514 END SUBROUTINE calculate_w_matrix_almo
1515
1516END MODULE almo_scf_qs
1517
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:1222
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.