(git:374b731)
Loading...
Searching...
No Matches
dm_ls_scf_qs.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Routines for a linear scaling quickstep SCF run based on the density
10!> matrix, with a focus on the interface between dm_ls_scf and qs
11!> \par History
12!> 2011.04 created [Joost VandeVondele]
13!> \author Joost VandeVondele
14! **************************************************************************************************
24 USE dbcsr_api, ONLY: &
25 dbcsr_complete_redistribute, dbcsr_copy, dbcsr_copy_into_existing, dbcsr_create, &
26 dbcsr_desymmetrize, dbcsr_distribution_get, dbcsr_distribution_hold, &
27 dbcsr_distribution_new, dbcsr_distribution_release, dbcsr_distribution_type, &
28 dbcsr_finalize, dbcsr_get_info, dbcsr_multiply, dbcsr_nblkrows_total, dbcsr_p_type, &
29 dbcsr_release, dbcsr_set, dbcsr_type, dbcsr_type_real_8
30 USE dm_ls_scf_types, ONLY: ls_cluster_atomic,&
31 ls_cluster_molecular,&
36 USE kinds, ONLY: default_string_length,&
37 dp
41 USE pw_env_types, ONLY: pw_env_get,&
43 USE pw_methods, ONLY: pw_zero
44 USE pw_pool_types, ONLY: pw_pool_p_type,&
46 USE pw_types, ONLY: pw_c1d_gs_type,&
59 USE qs_ks_types, ONLY: qs_ks_did_change,&
68 USE qs_rho_types, ONLY: qs_rho_get,&
72#include "./base/base_uses.f90"
73
74 IMPLICIT NONE
75
76 PRIVATE
77
78 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'dm_ls_scf_qs'
79
83
84CONTAINS
85
86! **************************************************************************************************
87!> \brief create a matrix for use (and as a template) in ls based on a qs template
88!> \param matrix_ls ...
89!> \param matrix_qs ...
90!> \param ls_mstruct ...
91!> \par History
92!> 2011.03 created [Joost VandeVondele]
93!> 2015.09 add support for PAO [Ole Schuett]
94!> \author Joost VandeVondele
95! **************************************************************************************************
96 SUBROUTINE matrix_ls_create(matrix_ls, matrix_qs, ls_mstruct)
97 TYPE(dbcsr_type) :: matrix_ls, matrix_qs
98 TYPE(ls_mstruct_type), INTENT(IN) :: ls_mstruct
99
100 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_ls_create'
101
102 CHARACTER(len=default_string_length) :: name
103 INTEGER :: handle, iatom, imol, jatom, &
104 ls_data_type, natom, nmol
105 INTEGER, ALLOCATABLE, DIMENSION(:), TARGET :: atom_to_cluster, atom_to_cluster_primus, &
106 clustered_blk_sizes, primus_of_mol
107 INTEGER, DIMENSION(:), POINTER :: clustered_col_dist, clustered_row_dist, &
108 ls_blk_sizes, ls_col_dist, ls_row_dist
109 TYPE(dbcsr_distribution_type) :: ls_dist, ls_dist_clustered
110
111 CALL timeset(routinen, handle)
112
113 ! Defaults -----------------------------------------------------------------------------------
114 CALL dbcsr_get_info(matrix_qs, col_blk_size=ls_blk_sizes, distribution=ls_dist)
115 CALL dbcsr_distribution_hold(ls_dist)
116 CALL dbcsr_distribution_get(ls_dist, row_dist=ls_row_dist, col_dist=ls_col_dist)
117 ls_data_type = dbcsr_type_real_8
118
119 ! PAO ----------------------------------------------------------------------------------------
120 IF (ls_mstruct%do_pao) THEN
121 CALL dbcsr_get_info(ls_mstruct%matrix_A, col_blk_size=ls_blk_sizes)
122 END IF
123
124 ! Clustering ---------------------------------------------------------------------------------
125 SELECT CASE (ls_mstruct%cluster_type)
126 CASE (ls_cluster_atomic)
127 ! do nothing
129 ! create format of the clustered matrix
130 natom = dbcsr_nblkrows_total(matrix_qs)
131 nmol = maxval(ls_mstruct%atom_to_molecule)
132 ALLOCATE (atom_to_cluster_primus(natom))
133 ALLOCATE (atom_to_cluster(natom))
134 ALLOCATE (primus_of_mol(nmol))
135 DO iatom = 1, natom
136 atom_to_cluster(iatom) = ls_mstruct%atom_to_molecule(iatom)
137 ! the first atom of the molecule is the primus
138 ! if the number of atoms per molecule is independent of system size, this is not a quadratic loop
139 ! it assumes that all atoms of the molecule are consecutive.
140 DO jatom = iatom, 1, -1
141 IF (ls_mstruct%atom_to_molecule(jatom) == atom_to_cluster(iatom)) THEN
142 atom_to_cluster_primus(iatom) = jatom
143 ELSE
144 EXIT
145 END IF
146 END DO
147 primus_of_mol(atom_to_cluster(iatom)) = atom_to_cluster_primus(iatom)
148 END DO
149
150 ! row
151 ALLOCATE (clustered_row_dist(nmol))
152 DO imol = 1, nmol
153 clustered_row_dist(imol) = ls_row_dist(primus_of_mol(imol))
154 END DO
155
156 ! col
157 ALLOCATE (clustered_col_dist(nmol))
158 DO imol = 1, nmol
159 clustered_col_dist(imol) = ls_col_dist(primus_of_mol(imol))
160 END DO
161
162 ALLOCATE (clustered_blk_sizes(nmol))
163 clustered_blk_sizes = 0
164 DO iatom = 1, natom
165 clustered_blk_sizes(atom_to_cluster(iatom)) = clustered_blk_sizes(atom_to_cluster(iatom)) + &
166 ls_blk_sizes(iatom)
167 END DO
168 ls_blk_sizes => clustered_blk_sizes ! redirect pointer
169
170 ! create new distribution
171 CALL dbcsr_distribution_new(ls_dist_clustered, &
172 template=ls_dist, &
173 row_dist=clustered_row_dist, &
174 col_dist=clustered_col_dist, &
175 reuse_arrays=.true.)
176 CALL dbcsr_distribution_release(ls_dist)
177 ls_dist = ls_dist_clustered
178
179 CASE DEFAULT
180 cpabort("Unknown LS cluster type")
181 END SELECT
182
183 ! Create actual matrix -----------------------------------------------------------------------
184 CALL dbcsr_get_info(matrix_qs, name=name)
185 CALL dbcsr_create(matrix_ls, &
186 name=name, &
187 dist=ls_dist, &
188 matrix_type="S", &
189 data_type=ls_data_type, &
190 row_blk_size=ls_blk_sizes, &
191 col_blk_size=ls_blk_sizes)
192 CALL dbcsr_distribution_release(ls_dist)
193 CALL dbcsr_finalize(matrix_ls)
194
195 CALL timestop(handle)
196
197 END SUBROUTINE matrix_ls_create
198
199! **************************************************************************************************
200!> \brief first link to QS, copy a QS matrix to LS matrix
201!> used to isolate QS style matrices from LS style
202!> will be useful for future features (e.g. precision, symmetry, blocking, ...)
203!> \param matrix_ls ...
204!> \param matrix_qs ...
205!> \param ls_mstruct ...
206!> \param covariant ...
207!> \par History
208!> 2010.10 created [Joost VandeVondele]
209!> 2015.09 add support for PAO [Ole Schuett]
210!> \author Joost VandeVondele
211! **************************************************************************************************
212 SUBROUTINE matrix_qs_to_ls(matrix_ls, matrix_qs, ls_mstruct, covariant)
213 TYPE(dbcsr_type) :: matrix_ls, matrix_qs
214 TYPE(ls_mstruct_type), INTENT(IN), TARGET :: ls_mstruct
215 LOGICAL, INTENT(IN) :: covariant
216
217 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_qs_to_ls'
218
219 INTEGER :: handle
220 INTEGER, DIMENSION(:), POINTER :: pao_blk_sizes
221 TYPE(dbcsr_type) :: matrix_pao, matrix_tmp
222 TYPE(dbcsr_type), POINTER :: matrix_trafo
223
224 CALL timeset(routinen, handle)
225
226 IF (.NOT. ls_mstruct%do_pao) THEN
227 CALL matrix_cluster(matrix_ls, matrix_qs, ls_mstruct)
228
229 ELSE ! using pao
230 CALL dbcsr_get_info(ls_mstruct%matrix_A, col_blk_size=pao_blk_sizes)
231 CALL dbcsr_create(matrix_pao, &
232 matrix_type="N", &
233 template=matrix_qs, &
234 row_blk_size=pao_blk_sizes, &
235 col_blk_size=pao_blk_sizes)
236
237 matrix_trafo => ls_mstruct%matrix_A ! contra-variant
238 IF (covariant) matrix_trafo => ls_mstruct%matrix_B ! co-variant
239 CALL dbcsr_create(matrix_tmp, template=matrix_trafo)
240
241 CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_qs, matrix_trafo, 0.0_dp, matrix_tmp)
242 CALL dbcsr_multiply("T", "N", 1.0_dp, matrix_trafo, matrix_tmp, 0.0_dp, matrix_pao)
243 CALL dbcsr_release(matrix_tmp)
244
245 CALL matrix_cluster(matrix_ls, matrix_pao, ls_mstruct)
246 CALL dbcsr_release(matrix_pao)
247 END IF
248
249 CALL timestop(handle)
250
251 END SUBROUTINE matrix_qs_to_ls
252
253! **************************************************************************************************
254!> \brief Performs molecular blocking and reduction to single precision if enabled
255!> \param matrix_out ...
256!> \param matrix_in ...
257!> \param ls_mstruct ...
258!> \author Ole Schuett
259! **************************************************************************************************
260 SUBROUTINE matrix_cluster(matrix_out, matrix_in, ls_mstruct)
261 TYPE(dbcsr_type) :: matrix_out, matrix_in
262 TYPE(ls_mstruct_type), INTENT(IN) :: ls_mstruct
263
264 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_cluster'
265
266 INTEGER :: handle
267 TYPE(dbcsr_type) :: matrix_in_nosym
268
269 CALL timeset(routinen, handle)
270
271 SELECT CASE (ls_mstruct%cluster_type)
272 CASE (ls_cluster_atomic)
273 CALL dbcsr_copy(matrix_out, matrix_in) ! takes care of an eventual data_type conversion
274
276 ! desymmetrize the qs matrix
277 CALL dbcsr_create(matrix_in_nosym, template=matrix_in, matrix_type="N")
278 CALL dbcsr_desymmetrize(matrix_in, matrix_in_nosym)
279
280 ! perform the magic complete redistribute copy
281 CALL dbcsr_complete_redistribute(matrix_in_nosym, matrix_out);
282 CALL dbcsr_release(matrix_in_nosym)
283
284 CASE DEFAULT
285 cpabort("Unknown LS cluster type")
286 END SELECT
287
288 CALL timestop(handle)
289
290 END SUBROUTINE matrix_cluster
291
292! **************************************************************************************************
293!> \brief second link to QS, copy a LS matrix to QS matrix
294!> used to isolate QS style matrices from LS style
295!> will be useful for future features (e.g. precision, symmetry, blocking, ...)
296!> \param matrix_qs ...
297!> \param matrix_ls ...
298!> \param ls_mstruct ...
299!> \param covariant ...
300!> \param keep_sparsity If set dbcsr_copy_into_existing will be used, by default set to .TRUE.
301!> \par History
302!> 2010.10 created [Joost VandeVondele]
303!> 2015.09 add support for PAO [Ole Schuett]
304!> \author Joost VandeVondele
305! **************************************************************************************************
306 SUBROUTINE matrix_ls_to_qs(matrix_qs, matrix_ls, ls_mstruct, covariant, keep_sparsity)
307 TYPE(dbcsr_type) :: matrix_qs, matrix_ls
308 TYPE(ls_mstruct_type), INTENT(IN), TARGET :: ls_mstruct
309 LOGICAL :: covariant
310 LOGICAL, OPTIONAL :: keep_sparsity
311
312 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_ls_to_qs'
313
314 INTEGER :: handle
315 INTEGER, DIMENSION(:), POINTER :: pao_blk_sizes
316 LOGICAL :: my_keep_sparsity
317 TYPE(dbcsr_type) :: matrix_declustered, matrix_tmp1, &
318 matrix_tmp2
319 TYPE(dbcsr_type), POINTER :: matrix_trafo
320
321 CALL timeset(routinen, handle)
322
323 my_keep_sparsity = .true.
324 IF (PRESENT(keep_sparsity)) &
325 my_keep_sparsity = keep_sparsity
326
327 IF (.NOT. ls_mstruct%do_pao) THEN
328 CALL dbcsr_create(matrix_declustered, template=matrix_qs)
329 CALL matrix_decluster(matrix_declustered, matrix_ls, ls_mstruct)
330 IF (my_keep_sparsity) THEN
331 CALL dbcsr_copy_into_existing(matrix_qs, matrix_declustered) ! preserve sparsity of matrix_qs
332 ELSE
333 CALL dbcsr_copy(matrix_qs, matrix_declustered) ! overwrite sparsity of matrix_qs
334 END IF
335 CALL dbcsr_release(matrix_declustered)
336
337 ELSE ! using pao
338 CALL dbcsr_get_info(ls_mstruct%matrix_A, col_blk_size=pao_blk_sizes)
339 CALL dbcsr_create(matrix_declustered, &
340 template=matrix_qs, &
341 row_blk_size=pao_blk_sizes, &
342 col_blk_size=pao_blk_sizes)
343
344 CALL matrix_decluster(matrix_declustered, matrix_ls, ls_mstruct)
345
346 matrix_trafo => ls_mstruct%matrix_B ! contra-variant
347 IF (covariant) matrix_trafo => ls_mstruct%matrix_A ! co-variant
348 CALL dbcsr_create(matrix_tmp1, template=matrix_trafo)
349 CALL dbcsr_create(matrix_tmp2, template=matrix_qs)
350 CALL dbcsr_multiply("N", "N", 1.0_dp, matrix_trafo, matrix_declustered, 0.0_dp, matrix_tmp1)
351 CALL dbcsr_multiply("N", "T", 1.0_dp, matrix_tmp1, matrix_trafo, 0.0_dp, matrix_tmp2)
352 IF (my_keep_sparsity) THEN
353 CALL dbcsr_copy_into_existing(matrix_qs, matrix_tmp2) ! preserve sparsity of matrix_qs
354 ELSE
355 CALL dbcsr_copy(matrix_qs, matrix_tmp2) ! overwrite sparsity of matrix_qs
356 END IF
357 CALL dbcsr_release(matrix_declustered)
358 CALL dbcsr_release(matrix_tmp1)
359 CALL dbcsr_release(matrix_tmp2)
360 END IF
361
362 CALL timestop(handle)
363
364 END SUBROUTINE matrix_ls_to_qs
365
366! **************************************************************************************************
367!> \brief Reverses molecular blocking and reduction to single precision if enabled
368!> \param matrix_out ...
369!> \param matrix_in ...
370!> \param ls_mstruct ...
371!> \author Ole Schuett
372! **************************************************************************************************
373 SUBROUTINE matrix_decluster(matrix_out, matrix_in, ls_mstruct)
374 TYPE(dbcsr_type) :: matrix_out, matrix_in
375 TYPE(ls_mstruct_type), INTENT(IN) :: ls_mstruct
376
377 CHARACTER(len=*), PARAMETER :: routinen = 'matrix_decluster'
378
379 INTEGER :: handle
380
381 CALL timeset(routinen, handle)
382
383 SELECT CASE (ls_mstruct%cluster_type)
384 CASE (ls_cluster_atomic)
385 CALL dbcsr_copy(matrix_out, matrix_in) ! takes care of an eventual data_type conversion
386
388 ! perform the magic complete redistribute copy
389 CALL dbcsr_complete_redistribute(matrix_in, matrix_out)
390
391 CASE DEFAULT
392 cpabort("Unknown LS cluster type")
393 END SELECT
394
395 CALL timestop(handle)
396
397 END SUBROUTINE matrix_decluster
398
399! **************************************************************************************************
400!> \brief further required initialization of QS.
401!> Might be factored-out since this seems common code with the other SCF.
402!> \param qs_env ...
403!> \par History
404!> 2010.10 created [Joost VandeVondele]
405!> \author Joost VandeVondele
406! **************************************************************************************************
407 SUBROUTINE ls_scf_init_qs(qs_env)
408 TYPE(qs_environment_type), POINTER :: qs_env
409
410 CHARACTER(len=*), PARAMETER :: routinen = 'ls_scf_init_qs'
411
412 INTEGER :: handle, ispin, nspin, unit_nr
413 TYPE(cp_logger_type), POINTER :: logger
414 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
415 TYPE(dft_control_type), POINTER :: dft_control
416 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
417 POINTER :: sab_orb
418 TYPE(qs_ks_env_type), POINTER :: ks_env
419
420 NULLIFY (sab_orb)
421 CALL timeset(routinen, handle)
422
423 ! get a useful output_unit
424 logger => cp_get_default_logger()
425 IF (logger%para_env%is_source()) THEN
426 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
427 ELSE
428 unit_nr = -1
429 END IF
430
431 ! get basic quantities from the qs_env
432 CALL get_qs_env(qs_env, dft_control=dft_control, &
433 matrix_s=matrix_s, &
434 matrix_ks=matrix_ks, &
435 ks_env=ks_env, &
436 sab_orb=sab_orb)
437
438 nspin = dft_control%nspins
439
440 ! we might have to create matrix_ks
441 IF (.NOT. ASSOCIATED(matrix_ks)) THEN
442 CALL dbcsr_allocate_matrix_set(matrix_ks, nspin)
443 DO ispin = 1, nspin
444 ALLOCATE (matrix_ks(ispin)%matrix)
445 CALL dbcsr_create(matrix_ks(ispin)%matrix, template=matrix_s(1)%matrix)
446 CALL cp_dbcsr_alloc_block_from_nbl(matrix_ks(ispin)%matrix, sab_orb)
447 CALL dbcsr_set(matrix_ks(ispin)%matrix, 0.0_dp)
448 END DO
449 CALL set_ks_env(ks_env, matrix_ks=matrix_ks)
450 END IF
451
452 CALL timestop(handle)
453
454 END SUBROUTINE ls_scf_init_qs
455
456! **************************************************************************************************
457!> \brief get an atomic initial guess
458!> \param qs_env ...
459!> \param energy ...
460!> \par History
461!> 2012.11 created [Joost VandeVondele]
462!> \author Joost VandeVondele
463! **************************************************************************************************
464 SUBROUTINE ls_scf_qs_atomic_guess(qs_env, energy)
465 TYPE(qs_environment_type), POINTER :: qs_env
466 REAL(kind=dp) :: energy
467
468 CHARACTER(len=*), PARAMETER :: routinen = 'ls_scf_qs_atomic_guess'
469
470 INTEGER :: handle, nspin, unit_nr
471 INTEGER, DIMENSION(2) :: nelectron_spin
472 LOGICAL :: has_unit_metric
473 TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
474 TYPE(cp_logger_type), POINTER :: logger
475 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s, rho_ao
476 TYPE(dft_control_type), POINTER :: dft_control
477 TYPE(mp_para_env_type), POINTER :: para_env
478 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
479 TYPE(qs_energy_type), POINTER :: qs_energy
480 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
481 TYPE(qs_ks_env_type), POINTER :: ks_env
482 TYPE(qs_rho_type), POINTER :: rho
483
484 CALL timeset(routinen, handle)
485 NULLIFY (rho, rho_ao)
486
487 ! get a useful output_unit
488 logger => cp_get_default_logger()
489 IF (logger%para_env%is_source()) THEN
490 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
491 ELSE
492 unit_nr = -1
493 END IF
494
495 ! get basic quantities from the qs_env
496 CALL get_qs_env(qs_env, dft_control=dft_control, &
497 matrix_s=matrix_s, &
498 matrix_ks=matrix_ks, &
499 ks_env=ks_env, &
500 energy=qs_energy, &
501 atomic_kind_set=atomic_kind_set, &
502 qs_kind_set=qs_kind_set, &
503 particle_set=particle_set, &
504 has_unit_metric=has_unit_metric, &
505 para_env=para_env, &
506 nelectron_spin=nelectron_spin, &
507 rho=rho)
508
509 CALL qs_rho_get(rho, rho_ao=rho_ao)
510
511 nspin = dft_control%nspins
512
513 ! create an initial atomic guess
514 IF (dft_control%qs_control%dftb .OR. dft_control%qs_control%semi_empirical .OR. &
515 dft_control%qs_control%xtb) THEN
516 CALL calculate_mopac_dm(rho_ao, matrix_s(1)%matrix, has_unit_metric, &
517 dft_control, particle_set, atomic_kind_set, qs_kind_set, &
518 nspin, nelectron_spin, para_env)
519 ELSE
520 CALL calculate_atomic_block_dm(rho_ao, matrix_s(1)%matrix, atomic_kind_set, qs_kind_set, &
521 nspin, nelectron_spin, unit_nr, para_env)
522 END IF
523
524 CALL qs_rho_update_rho(rho, qs_env=qs_env)
525 CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.true.)
526 CALL qs_ks_update_qs_env(qs_env, calculate_forces=.false., just_energy=.false.)
527
528 energy = qs_energy%total
529
530 CALL timestop(handle)
531
532 END SUBROUTINE ls_scf_qs_atomic_guess
533
534! **************************************************************************************************
535!> \brief use the density matrix in ls_scf_env to compute the new energy and KS matrix
536!> \param qs_env ...
537!> \param ls_scf_env ...
538!> \param energy_new ...
539!> \param iscf ...
540!> \par History
541!> 2011.04 created [Joost VandeVondele]
542!> 2015.02 added gspace density mixing [Patrick Seewald]
543!> \author Joost VandeVondele
544! **************************************************************************************************
545 SUBROUTINE ls_scf_dm_to_ks(qs_env, ls_scf_env, energy_new, iscf)
546 TYPE(qs_environment_type), POINTER :: qs_env
547 TYPE(ls_scf_env_type) :: ls_scf_env
548 REAL(kind=dp) :: energy_new
549 INTEGER, INTENT(IN) :: iscf
550
551 CHARACTER(len=*), PARAMETER :: routinen = 'ls_scf_dm_to_ks'
552
553 INTEGER :: handle, ispin, nspin, unit_nr
554 TYPE(cp_logger_type), POINTER :: logger
555 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: rho_ao
556 TYPE(mp_para_env_type), POINTER :: para_env
557 TYPE(qs_energy_type), POINTER :: energy
558 TYPE(qs_rho_type), POINTER :: rho
559
560 NULLIFY (energy, rho, rho_ao)
561 CALL timeset(routinen, handle)
562
563 logger => cp_get_default_logger()
564 IF (logger%para_env%is_source()) THEN
565 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
566 ELSE
567 unit_nr = -1
568 END IF
569
570 nspin = ls_scf_env%nspins
571 CALL get_qs_env(qs_env, para_env=para_env, energy=energy, rho=rho)
572 CALL qs_rho_get(rho, rho_ao=rho_ao)
573
574 ! set the new density matrix
575 DO ispin = 1, nspin
576 CALL matrix_ls_to_qs(rho_ao(ispin)%matrix, ls_scf_env%matrix_p(ispin), &
577 ls_scf_env%ls_mstruct, covariant=.false.)
578 END DO
579
580 ! compute the corresponding KS matrix and new energy, mix density if requested
581 CALL qs_rho_update_rho(rho, qs_env=qs_env)
582 IF (ls_scf_env%do_rho_mixing) THEN
583 IF (ls_scf_env%density_mixing_method .EQ. direct_mixing_nr) &
584 cpabort("Direct P mixing not implemented in linear scaling SCF. ")
585 IF (ls_scf_env%density_mixing_method >= gspace_mixing_nr) THEN
586 IF (iscf .GT. max(ls_scf_env%mixing_store%nskip_mixing, 1)) THEN
587 CALL gspace_mixing(qs_env, ls_scf_env%density_mixing_method, &
588 ls_scf_env%mixing_store, rho, para_env, &
589 iscf - 1)
590 IF (unit_nr > 0) THEN
591 WRITE (unit_nr, '(A57)') &
592 "*********************************************************"
593 WRITE (unit_nr, '(A13,F5.3,A20,A6,A7,I3)') &
594 " Using ALPHA=", ls_scf_env%mixing_store%alpha, &
595 " to mix rho: method=", ls_scf_env%mixing_store%iter_method, ", iscf=", iscf
596 WRITE (unit_nr, '(A8,F5.3,A6,F5.3,A8)') &
597 " rho_nw=", ls_scf_env%mixing_store%alpha, "*rho + ", &
598 1.0_dp - ls_scf_env%mixing_store%alpha, "*rho_old"
599 WRITE (unit_nr, '(A57)') &
600 "*********************************************************"
601 END IF
602 END IF
603 END IF
604 END IF
605
606 CALL qs_ks_did_change(qs_env%ks_env, rho_changed=.true.)
607 CALL qs_ks_update_qs_env(qs_env, calculate_forces=.false., &
608 just_energy=.false., print_active=.true.)
609 energy_new = energy%total
610
611 CALL timestop(handle)
612
613 END SUBROUTINE ls_scf_dm_to_ks
614
615! **************************************************************************************************
616!> \brief ...
617!> \param qs_env ...
618!> \param ls_scf_env ...
619!> \param matrix_p_ls ...
620!> \param unit_nr ...
621!> \param title ...
622!> \param stride ...
623! **************************************************************************************************
624 SUBROUTINE write_matrix_to_cube(qs_env, ls_scf_env, matrix_p_ls, unit_nr, title, stride)
625 TYPE(qs_environment_type), POINTER :: qs_env
626 TYPE(ls_scf_env_type) :: ls_scf_env
627 TYPE(dbcsr_type), INTENT(IN) :: matrix_p_ls
628 INTEGER, INTENT(IN) :: unit_nr
629 CHARACTER(LEN=*), INTENT(IN) :: title
630 INTEGER, DIMENSION(:), POINTER :: stride
631
632 CHARACTER(len=*), PARAMETER :: routinen = 'write_matrix_to_cube'
633
634 INTEGER :: handle
635 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks
636 TYPE(dbcsr_type), TARGET :: matrix_p_qs
637 TYPE(particle_list_type), POINTER :: particles
638 TYPE(pw_c1d_gs_type) :: wf_g
639 TYPE(pw_env_type), POINTER :: pw_env
640 TYPE(pw_pool_p_type), DIMENSION(:), POINTER :: pw_pools
641 TYPE(pw_pool_type), POINTER :: auxbas_pw_pool
642 TYPE(pw_r3d_rs_type) :: wf_r
643 TYPE(qs_ks_env_type), POINTER :: ks_env
644 TYPE(qs_subsys_type), POINTER :: subsys
645
646 CALL timeset(routinen, handle)
647
648 NULLIFY (ks_env, pw_env, auxbas_pw_pool, pw_pools, particles, subsys, matrix_ks)
649
650 CALL get_qs_env(qs_env, &
651 ks_env=ks_env, &
652 subsys=subsys, &
653 pw_env=pw_env, &
654 matrix_ks=matrix_ks)
655
656 CALL qs_subsys_get(subsys, particles=particles)
657
658 ! convert the density matrix (ls style) to QS style
659 CALL dbcsr_copy(matrix_p_qs, matrix_ks(1)%matrix)
660 CALL dbcsr_set(matrix_p_qs, 0.0_dp) !zero matrix creation
661 CALL matrix_ls_to_qs(matrix_p_qs, matrix_p_ls, ls_scf_env%ls_mstruct, covariant=.false.)
662
663 ! Print total electronic density
664 CALL pw_env_get(pw_env=pw_env, &
665 auxbas_pw_pool=auxbas_pw_pool, &
666 pw_pools=pw_pools)
667 CALL auxbas_pw_pool%create_pw(pw=wf_r)
668 CALL pw_zero(wf_r)
669 CALL auxbas_pw_pool%create_pw(pw=wf_g)
670 CALL pw_zero(wf_g)
671 CALL calculate_rho_elec(matrix_p=matrix_p_qs, &
672 rho=wf_r, &
673 rho_gspace=wf_g, &
674 ks_env=ks_env)
675
676 ! write this to a cube
677 CALL cp_pw_to_cube(wf_r, unit_nr=unit_nr, title=title, &
678 particles=particles, stride=stride)
679
680 !free memory
681 CALL auxbas_pw_pool%give_back_pw(wf_r)
682 CALL auxbas_pw_pool%give_back_pw(wf_g)
683 CALL dbcsr_release(matrix_p_qs)
684
685 CALL timestop(handle)
686
687 END SUBROUTINE write_matrix_to_cube
688
689! **************************************************************************************************
690!> \brief Initialize g-space density mixing
691!> \param qs_env ...
692!> \param ls_scf_env ...
693! **************************************************************************************************
694 SUBROUTINE rho_mixing_ls_init(qs_env, ls_scf_env)
695 TYPE(qs_environment_type), POINTER :: qs_env
696 TYPE(ls_scf_env_type) :: ls_scf_env
697
698 CHARACTER(len=*), PARAMETER :: routinen = 'rho_mixing_ls_init'
699
700 INTEGER :: handle
701 TYPE(dft_control_type), POINTER :: dft_control
702 TYPE(qs_rho_type), POINTER :: rho
703 TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom
704
705 CALL timeset(routinen, handle)
706
707 CALL get_qs_env(qs_env, dft_control=dft_control, rho=rho)
708
709 CALL mixing_allocate(qs_env, ls_scf_env%density_mixing_method, nspins=ls_scf_env%nspins, &
710 mixing_store=ls_scf_env%mixing_store)
711 IF (ls_scf_env%density_mixing_method >= gspace_mixing_nr) THEN
712 IF (dft_control%qs_control%gapw) THEN
713 CALL get_qs_env(qs_env, rho_atom_set=rho_atom)
714 CALL mixing_init(ls_scf_env%density_mixing_method, rho, ls_scf_env%mixing_store, &
715 ls_scf_env%para_env, rho_atom=rho_atom)
716 ELSEIF (dft_control%qs_control%dftb .OR. dft_control%qs_control%xtb) THEN
717 CALL charge_mixing_init(ls_scf_env%mixing_store)
718 ELSEIF (dft_control%qs_control%semi_empirical) THEN
719 cpabort('SE Code not possible')
720 ELSE
721 CALL mixing_init(ls_scf_env%density_mixing_method, rho, ls_scf_env%mixing_store, &
722 ls_scf_env%para_env)
723 END IF
724 END IF
725 CALL timestop(handle)
726 END SUBROUTINE rho_mixing_ls_init
727
728END MODULE dm_ls_scf_qs
Define the atomic kind types and their sub types.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
DBCSR operations in CP2K.
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
A wrapper around pw_to_cube() which accepts particle_list_type.
subroutine, public cp_pw_to_cube(pw, unit_nr, title, particles, stride, zero_tails, silent, mpi_io)
...
Routines for a linear scaling quickstep SCF run based on the density matrix, with a focus on the inte...
subroutine, public matrix_ls_to_qs(matrix_qs, matrix_ls, ls_mstruct, covariant, keep_sparsity)
second link to QS, copy a LS matrix to QS matrix used to isolate QS style matrices from LS style will...
subroutine, public matrix_decluster(matrix_out, matrix_in, ls_mstruct)
Reverses molecular blocking and reduction to single precision if enabled.
subroutine, public ls_scf_dm_to_ks(qs_env, ls_scf_env, energy_new, iscf)
use the density matrix in ls_scf_env to compute the new energy and KS matrix
subroutine, public write_matrix_to_cube(qs_env, ls_scf_env, matrix_p_ls, unit_nr, title, stride)
...
subroutine, public rho_mixing_ls_init(qs_env, ls_scf_env)
Initialize g-space density mixing.
subroutine, public matrix_ls_create(matrix_ls, matrix_qs, ls_mstruct)
create a matrix for use (and as a template) in ls based on a qs template
subroutine, public matrix_qs_to_ls(matrix_ls, matrix_qs, ls_mstruct, covariant)
first link to QS, copy a QS matrix to LS matrix used to isolate QS style matrices from LS style will ...
subroutine, public ls_scf_init_qs(qs_env)
further required initialization of QS. Might be factored-out since this seems common code with the ot...
subroutine, public ls_scf_qs_atomic_guess(qs_env, energy)
get an atomic initial guess
Types needed for a linear scaling quickstep SCF run based on the density matrix.
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public ls_cluster_molecular
integer, parameter, public ls_cluster_atomic
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
integer, parameter, public default_string_length
Definition kinds.F:57
Interface to the message passing library MPI.
represent a simple array based list of the given type
Define the data structure for the particle information.
container for various plainwaves related things
subroutine, public pw_env_get(pw_env, pw_pools, cube_info, gridlevel_info, auxbas_pw_pool, auxbas_grid, auxbas_rs_desc, auxbas_rs_grid, rs_descs, rs_grids, xc_pw_pool, vdw_pw_pool, poisson_env, interp_section)
returns the various attributes of the pw env
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Routine to return block diagonal density matrix. Blocks correspond to the atomic densities.
subroutine, public calculate_atomic_block_dm(pmatrix, matrix_s, atomic_kind_set, qs_kind_set, nspin, nelectron_spin, ounit, para_env)
returns a block diagonal density matrix. Blocks correspond to the atomic densities.
Calculate the plane wave density by collocating the primitive Gaussian functions (pgf).
subroutine, public calculate_rho_elec(matrix_p, matrix_p_kp, rho, rho_gspace, total_rho, ks_env, soft_valid, compute_tau, compute_grad, basis_type, der_type, idir, task_list_external, pw_env_external)
computes the density corresponding to a given density matrix on the grid
module that contains the definitions of the scf types
integer, parameter, public direct_mixing_nr
integer, parameter, public gspace_mixing_nr
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, 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_nonbond, sab_almo, sab_kp, sab_kp_nosym, 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, 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, ecoul_1c, rho0_s_rs, rho0_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, 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, rhs)
Get the QUICKSTEP environment.
subroutine, public gspace_mixing(qs_env, mixing_method, mixing_store, rho, para_env, iter_count)
Driver for the g-space mixing, calls the proper routine given the requested method.
Routines to somehow generate an initial guess.
subroutine, public calculate_mopac_dm(pmat, matrix_s, has_unit_metric, dft_control, particle_set, atomic_kind_set, qs_kind_set, nspin, nelectron_spin, para_env)
returns a block diagonal density matrix. Blocks correspond to the mopac initial guess.
Define the quickstep kind type and their sub types.
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 set_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, 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, 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_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, task_list, task_list_soft, subsys, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env)
...
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...
elemental subroutine, public charge_mixing_init(mixing_store)
initialiation needed when charge mixing is used
subroutine, public mixing_init(mixing_method, rho, mixing_store, para_env, rho_atom)
initialiation needed when gspace mixing is used
subroutine, public mixing_allocate(qs_env, mixing_method, p_mix_new, p_delta, nspins, mixing_store)
allocation needed when density mixing is used
Define the neighbor list data types and the corresponding functionality.
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...
types that represent a quickstep subsys
subroutine, public qs_subsys_get(subsys, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, cp_subsys, nelectron_total, nelectron_spin)
...
Provides all information about an atomic kind.
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment
contained for different pw related things
to create arrays of pools
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
Provides all information about a quickstep kind.
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.