(git:9754b87)
Loading...
Searching...
No Matches
qs_scf_diagonalization.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Different diagonalization schemes that can be used
10!> for the iterative solution of the eigenvalue problem
11!> \par History
12!> started from routines previously located in the qs_scf module
13!> 05.2009
14! **************************************************************************************************
20 USE cp_cfm_diag, ONLY: cp_cfm_geeig,&
22 USE cp_cfm_types, ONLY: cp_cfm_create,&
28 USE cp_dbcsr_api, ONLY: &
30 dbcsr_set, dbcsr_type, dbcsr_type_antisymmetric, dbcsr_type_no_symmetry, &
31 dbcsr_type_symmetric
50 USE cp_fm_types, ONLY: &
58 USE input_constants, ONLY: &
63 USE kinds, ONLY: dp
68 USE kpoint_types, ONLY: get_kpoint_info,&
71 USE machine, ONLY: m_flush,&
89 USE qs_ks_types, ONLY: qs_ks_did_change,&
91 USE qs_matrix_pools, ONLY: mpools_get,&
99 USE qs_mo_types, ONLY: get_mo_set,&
105 USE qs_rho_types, ONLY: qs_rho_get,&
117 USE qs_scf_types, ONLY: qs_scf_env_type,&
120#include "./base/base_uses.f90"
121
122 IMPLICIT NONE
123
124 PRIVATE
125
126 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_scf_diagonalization'
127
131
132CONTAINS
133
134! **************************************************************************************************
135!> \brief the inner loop of scf, specific to diagonalization with S matrix
136!> basically, in goes the ks matrix out goes a new p matrix
137!> \param scf_env ...
138!> \param mos ...
139!> \param matrix_ks ...
140!> \param matrix_s ...
141!> \param scf_control ...
142!> \param scf_section ...
143!> \param diis_step ...
144!> \par History
145!> 03.2006 created [Joost VandeVondele]
146! **************************************************************************************************
147
148 SUBROUTINE general_eigenproblem(scf_env, mos, matrix_ks, &
149 matrix_s, scf_control, scf_section, &
150 diis_step)
151
152 TYPE(qs_scf_env_type), POINTER :: scf_env
153 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mos
154 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
155 TYPE(scf_control_type), POINTER :: scf_control
156 TYPE(section_vals_type), POINTER :: scf_section
157 LOGICAL, INTENT(INOUT) :: diis_step
158
159 INTEGER :: ispin, nspin
160 LOGICAL :: do_level_shift, owns_ortho, use_jacobi
161 REAL(kind=dp) :: diis_error, eps_diis
162 TYPE(cp_fm_type), POINTER :: ortho
163 TYPE(dbcsr_type), POINTER :: ortho_dbcsr
164
165 nspin = SIZE(matrix_ks)
166 NULLIFY (ortho, ortho_dbcsr)
167
168 DO ispin = 1, nspin
169 CALL copy_dbcsr_to_fm(matrix_ks(ispin)%matrix, &
170 scf_env%scf_work1(ispin))
171 END DO
172
173 eps_diis = scf_control%eps_diis
174
175 IF (scf_env%iter_count > 1 .AND. .NOT. scf_env%skip_diis) THEN
176 CALL qs_diis_b_step(scf_env%scf_diis_buffer, mos, scf_env%scf_work1, &
177 scf_env%scf_work2, scf_env%iter_delta, diis_error, diis_step, &
178 eps_diis, scf_control%nmixing, &
179 s_matrix=matrix_s, &
180 scf_section=scf_section)
181 ELSE
182 diis_step = .false.
183 END IF
184
185 do_level_shift = ((scf_control%level_shift /= 0.0_dp) .AND. &
186 ((scf_control%density_guess == core_guess) .OR. &
187 (scf_env%iter_count > 1)))
188
189 IF ((scf_env%iter_count > 1) .AND. &
190 (scf_env%iter_delta < scf_control%diagonalization%eps_jacobi)) THEN
191 use_jacobi = .true.
192 ELSE
193 use_jacobi = .false.
194 END IF
195
196 IF (diis_step) THEN
197 scf_env%iter_param = diis_error
198 IF (use_jacobi) THEN
199 scf_env%iter_method = "DIIS/Jacobi"
200 ELSE
201 scf_env%iter_method = "DIIS/Diag."
202 END IF
203 ELSE
204 IF (scf_env%mixing_method == 0) THEN
205 scf_env%iter_method = "NoMix/Diag."
206 ELSE IF (scf_env%mixing_method == 1) THEN
207 scf_env%iter_param = scf_env%p_mix_alpha
208 IF (use_jacobi) THEN
209 scf_env%iter_method = "P_Mix/Jacobi"
210 ELSE
211 scf_env%iter_method = "P_Mix/Diag."
212 END IF
213 ELSEIF (scf_env%mixing_method > 1) THEN
214 scf_env%iter_param = scf_env%mixing_store%alpha
215 IF (use_jacobi) THEN
216 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Jacobi"
217 ELSE
218 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Diag."
219 END IF
220 END IF
221 END IF
222
223 IF (scf_env%cholesky_method == cholesky_dbcsr) THEN
224 ortho_dbcsr => scf_env%ortho_dbcsr
225 DO ispin = 1, nspin
226 CALL eigensolver_dbcsr(matrix_ks=matrix_ks(ispin)%matrix, matrix_ks_fm=scf_env%scf_work1(ispin), &
227 mo_set=mos(ispin), &
228 ortho_dbcsr=ortho_dbcsr, &
229 ksbuf1=scf_env%buf1_dbcsr, ksbuf2=scf_env%buf2_dbcsr)
230 END DO
231
232 ELSE IF (scf_env%cholesky_method > cholesky_off) THEN
233 IF (scf_env%cholesky_method == cholesky_inverse) THEN
234 ortho => scf_env%ortho_m1
235 ELSE
236 ortho => scf_env%ortho
237 END IF
238
239 owns_ortho = .false.
240 IF (.NOT. ASSOCIATED(ortho)) THEN
241 ALLOCATE (ortho)
242 owns_ortho = .true.
243 END IF
244
245 DO ispin = 1, nspin
246 IF (do_level_shift) THEN
247 CALL eigensolver(matrix_ks_fm=scf_env%scf_work1(ispin), &
248 mo_set=mos(ispin), &
249 ortho=ortho, &
250 work=scf_env%scf_work2, &
251 cholesky_method=scf_env%cholesky_method, &
252 do_level_shift=do_level_shift, &
253 level_shift=scf_control%level_shift, &
254 matrix_u_fm=scf_env%ortho, &
255 use_jacobi=use_jacobi)
256 ELSE
257 CALL eigensolver(matrix_ks_fm=scf_env%scf_work1(ispin), &
258 mo_set=mos(ispin), &
259 ortho=ortho, &
260 work=scf_env%scf_work2, &
261 cholesky_method=scf_env%cholesky_method, &
262 do_level_shift=do_level_shift, &
263 level_shift=scf_control%level_shift, &
264 use_jacobi=use_jacobi)
265 END IF
266 END DO
267
268 IF (owns_ortho) DEALLOCATE (ortho)
269 ELSE
270 ortho => scf_env%ortho
271
272 owns_ortho = .false.
273 IF (.NOT. ASSOCIATED(ortho)) THEN
274 ALLOCATE (ortho)
275 owns_ortho = .true.
276 END IF
277
278 IF (do_level_shift) THEN
279 DO ispin = 1, nspin
280 IF (ASSOCIATED(scf_env%scf_work1_red) .AND. ASSOCIATED(scf_env%scf_work2_red) &
281 .AND. ASSOCIATED(scf_env%ortho_red) .AND. ASSOCIATED(scf_env%ortho_m1_red)) THEN
282 CALL eigensolver_symm(matrix_ks_fm=scf_env%scf_work1(ispin), &
283 mo_set=mos(ispin), &
284 ortho=ortho, &
285 work=scf_env%scf_work2, &
286 do_level_shift=do_level_shift, &
287 level_shift=scf_control%level_shift, &
288 matrix_u_fm=scf_env%ortho_m1, &
289 use_jacobi=use_jacobi, &
290 jacobi_threshold=scf_control%diagonalization%jacobi_threshold, &
291 matrix_ks_fm_red=scf_env%scf_work1_red(ispin), &
292 ortho_red=scf_env%ortho_red, &
293 work_red=scf_env%scf_work2_red, &
294 matrix_u_fm_red=scf_env%ortho_m1_red)
295 ELSE
296 CALL eigensolver_symm(matrix_ks_fm=scf_env%scf_work1(ispin), &
297 mo_set=mos(ispin), &
298 ortho=ortho, &
299 work=scf_env%scf_work2, &
300 do_level_shift=do_level_shift, &
301 level_shift=scf_control%level_shift, &
302 matrix_u_fm=scf_env%ortho_m1, &
303 use_jacobi=use_jacobi, &
304 jacobi_threshold=scf_control%diagonalization%jacobi_threshold)
305 END IF
306 END DO
307 ELSE
308 DO ispin = 1, nspin
309 IF (ASSOCIATED(scf_env%scf_work1_red) .AND. ASSOCIATED(scf_env%scf_work2_red) &
310 .AND. ASSOCIATED(scf_env%ortho_red)) THEN
311 CALL eigensolver_symm(matrix_ks_fm=scf_env%scf_work1(ispin), &
312 mo_set=mos(ispin), &
313 ortho=ortho, &
314 work=scf_env%scf_work2, &
315 do_level_shift=do_level_shift, &
316 level_shift=scf_control%level_shift, &
317 use_jacobi=use_jacobi, &
318 jacobi_threshold=scf_control%diagonalization%jacobi_threshold, &
319 matrix_ks_fm_red=scf_env%scf_work1_red(ispin), &
320 ortho_red=scf_env%ortho_red, &
321 work_red=scf_env%scf_work2_red)
322 ELSE
323 CALL eigensolver_symm(matrix_ks_fm=scf_env%scf_work1(ispin), &
324 mo_set=mos(ispin), &
325 ortho=ortho, &
326 work=scf_env%scf_work2, &
327 do_level_shift=do_level_shift, &
328 level_shift=scf_control%level_shift, &
329 use_jacobi=use_jacobi, &
330 jacobi_threshold=scf_control%diagonalization%jacobi_threshold)
331 END IF
332 END DO
333 END IF
334
335 IF (owns_ortho) DEALLOCATE (ortho)
336 END IF
337
338 END SUBROUTINE general_eigenproblem
339
340! **************************************************************************************************
341!> \brief ...
342!> \param scf_env ...
343!> \param mos ...
344!> \param matrix_ks ...
345!> \param matrix_s ...
346!> \param scf_control ...
347!> \param scf_section ...
348!> \param diis_step ...
349! **************************************************************************************************
350 SUBROUTINE do_general_diag(scf_env, mos, matrix_ks, &
351 matrix_s, scf_control, scf_section, &
352 diis_step)
353
354 TYPE(qs_scf_env_type), POINTER :: scf_env
355 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
356 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
357 TYPE(scf_control_type), POINTER :: scf_control
358 TYPE(section_vals_type), POINTER :: scf_section
359 LOGICAL, INTENT(INOUT) :: diis_step
360
361 INTEGER :: ispin, nspin
362 REAL(kind=dp) :: total_zeff_corr
363
364 nspin = SIZE(matrix_ks)
365
366 CALL general_eigenproblem(scf_env, mos, matrix_ks, &
367 matrix_s, scf_control, scf_section, diis_step)
368
369 total_zeff_corr = 0.0_dp
370 total_zeff_corr = scf_env%sum_zeff_corr
371
372 IF (abs(total_zeff_corr) > 0.0_dp) THEN
373 CALL set_mo_occupation(mo_array=mos, &
374 smear=scf_control%smear, tot_zeff_corr=total_zeff_corr)
375 ELSE
376 CALL set_mo_occupation(mo_array=mos, &
377 smear=scf_control%smear)
378 END IF
379
380 DO ispin = 1, nspin
381 CALL calculate_density_matrix(mos(ispin), &
382 scf_env%p_mix_new(ispin, 1)%matrix)
383 END DO
384
385 END SUBROUTINE do_general_diag
386
387! **************************************************************************************************
388!> \brief Kpoint diagonalization routine
389!> Transforms matrices to kpoint, distributes kpoint groups, performs
390!> general diagonalization (no storgae of overlap decomposition), stores
391!> MOs, calculates occupation numbers, calculates density matrices
392!> in kpoint representation, transforms density matrices to real space
393!> \param matrix_ks Kohn-sham matrices (RS indices, global)
394!> \param matrix_s Overlap matrices (RS indices, global)
395!> \param kpoints Kpoint environment
396!> \param scf_env SCF environment
397!> \param scf_control SCF control variables
398!> \param update_p ...
399!> \param diis_step ...
400!> \param diis_error ...
401!> \param qs_env ...
402!> \par History
403!> 08.2014 created [JGH]
404! **************************************************************************************************
405 SUBROUTINE do_general_diag_kp(matrix_ks, matrix_s, kpoints, scf_env, scf_control, update_p, &
406 diis_step, diis_error, qs_env)
407
408 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_ks, matrix_s
409 TYPE(kpoint_type), POINTER :: kpoints
410 TYPE(qs_scf_env_type), POINTER :: scf_env
411 TYPE(scf_control_type), POINTER :: scf_control
412 LOGICAL, INTENT(IN) :: update_p
413 LOGICAL, INTENT(INOUT) :: diis_step
414 REAL(dp), INTENT(INOUT), OPTIONAL :: diis_error
415 TYPE(qs_environment_type), OPTIONAL, POINTER :: qs_env
416
417 CHARACTER(len=*), PARAMETER :: routinen = 'do_general_diag_kp'
418 COMPLEX(KIND=dp), PARAMETER :: cone = cmplx(1.0_dp, 0.0_dp, kind=dp), &
419 czero = cmplx(0.0_dp, 0.0_dp, kind=dp), ione = cmplx(0.0_dp, 1.0_dp, kind=dp)
420
421 COMPLEX(KIND=dp), ALLOCATABLE, DIMENSION(:) :: coeffs
422 INTEGER :: handle, ib, igroup, ik, ikp, indx, &
423 ispin, jb, kplocal, nb, nkp, &
424 nkp_groups, nspin
425 INTEGER, DIMENSION(2) :: kp_range
426 INTEGER, DIMENSION(:, :), POINTER :: kp_dist
427 INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
428 LOGICAL :: do_diis, my_kpgrp, use_real_wfn
429 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues
430 REAL(kind=dp), DIMENSION(:, :), POINTER :: xkp
431 TYPE(copy_info_type), ALLOCATABLE, DIMENSION(:, :) :: info
432 TYPE(cp_cfm_type) :: cksmat, cmos, csmat, cwork
433 TYPE(cp_fm_pool_p_type), DIMENSION(:), POINTER :: ao_ao_fm_pools
434 TYPE(cp_fm_struct_type), POINTER :: matrix_struct, mo_struct
435 TYPE(cp_fm_type) :: fmdummy, fmlocal, rksmat, rsmat
436 TYPE(cp_fm_type), DIMENSION(:), POINTER :: fmwork
437 TYPE(cp_fm_type), POINTER :: imos, mo_coeff, rmos
438 TYPE(dbcsr_type), POINTER :: cmatrix, rmatrix, tmpmat
439 TYPE(kpoint_env_type), POINTER :: kp
440 TYPE(mp_para_env_type), POINTER :: para_env, para_env_global
441 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
442 POINTER :: sab_nl
443 TYPE(qs_matrix_pools_type), POINTER :: mpools
444 TYPE(section_vals_type), POINTER :: scf_section
445
446 CALL timeset(routinen, handle)
447
448 NULLIFY (sab_nl)
449 CALL get_kpoint_info(kpoints, nkp=nkp, xkp=xkp, use_real_wfn=use_real_wfn, kp_range=kp_range, &
450 nkp_groups=nkp_groups, kp_dist=kp_dist, sab_nl=sab_nl, &
451 cell_to_index=cell_to_index)
452 cpassert(ASSOCIATED(sab_nl))
453 kplocal = kp_range(2) - kp_range(1) + 1
454
455 !Whether we use DIIS for k-points
456 do_diis = .false.
457 IF (scf_env%iter_count > 1 .AND. .NOT. scf_env%skip_diis .AND. .NOT. use_real_wfn &
458 .AND. PRESENT(diis_error) .AND. PRESENT(qs_env)) do_diis = .true.
459
460 ! allocate some work matrices
461 ALLOCATE (rmatrix, cmatrix, tmpmat)
462 CALL dbcsr_create(rmatrix, template=matrix_ks(1, 1)%matrix, &
463 matrix_type=dbcsr_type_symmetric)
464 CALL dbcsr_create(cmatrix, template=matrix_ks(1, 1)%matrix, &
465 matrix_type=dbcsr_type_antisymmetric)
466 CALL dbcsr_create(tmpmat, template=matrix_ks(1, 1)%matrix, &
467 matrix_type=dbcsr_type_no_symmetry)
468 CALL cp_dbcsr_alloc_block_from_nbl(rmatrix, sab_nl)
469 CALL cp_dbcsr_alloc_block_from_nbl(cmatrix, sab_nl)
470
471 fmwork => scf_env%scf_work1
472
473 ! fm pools to be used within a kpoint group
474 CALL get_kpoint_info(kpoints, mpools=mpools)
475 CALL mpools_get(mpools, ao_ao_fm_pools=ao_ao_fm_pools)
476
477 CALL fm_pool_create_fm(ao_ao_fm_pools(1)%pool, fmlocal)
478 CALL cp_fm_get_info(fmlocal, matrix_struct=matrix_struct)
479
480 IF (use_real_wfn) THEN
481 CALL cp_fm_create(rksmat, matrix_struct)
482 CALL cp_fm_create(rsmat, matrix_struct)
483 ELSE
484 CALL cp_cfm_create(cksmat, matrix_struct)
485 CALL cp_cfm_create(csmat, matrix_struct)
486 CALL cp_cfm_create(cwork, matrix_struct)
487 kp => kpoints%kp_env(1)%kpoint_env
488 CALL get_mo_set(kp%mos(1, 1), mo_coeff=mo_coeff)
489 CALL cp_fm_get_info(mo_coeff, matrix_struct=mo_struct)
490 CALL cp_cfm_create(cmos, mo_struct)
491 END IF
492
493 para_env => kpoints%blacs_env_all%para_env
494 nspin = SIZE(matrix_ks, 1)
495 ALLOCATE (info(kplocal*nspin*nkp_groups, 4))
496
497 ! Setup and start all the communication
498 indx = 0
499 DO ikp = 1, kplocal
500 DO ispin = 1, nspin
501 DO igroup = 1, nkp_groups
502 ! number of current kpoint
503 ik = kp_dist(1, igroup) + ikp - 1
504 my_kpgrp = (ik >= kpoints%kp_range(1) .AND. ik <= kpoints%kp_range(2))
505 indx = indx + 1
506 IF (use_real_wfn) THEN
507 ! FT of matrices KS and S, then transfer to FM type
508 CALL dbcsr_set(rmatrix, 0.0_dp)
509 CALL rskp_transform(rmatrix=rmatrix, rsmat=matrix_ks, ispin=ispin, &
510 xkp=xkp(1:3, ik), cell_to_index=cell_to_index, sab_nl=sab_nl)
511 CALL dbcsr_desymmetrize(rmatrix, tmpmat)
512 CALL copy_dbcsr_to_fm(tmpmat, fmwork(1))
513 ! s matrix is not spin dependent
514 CALL dbcsr_set(rmatrix, 0.0_dp)
515 CALL rskp_transform(rmatrix=rmatrix, rsmat=matrix_s, ispin=1, &
516 xkp=xkp(1:3, ik), cell_to_index=cell_to_index, sab_nl=sab_nl)
517 CALL dbcsr_desymmetrize(rmatrix, tmpmat)
518 CALL copy_dbcsr_to_fm(tmpmat, fmwork(3))
519 ELSE
520 ! FT of matrices KS and S, then transfer to FM type
521 CALL dbcsr_set(rmatrix, 0.0_dp)
522 CALL dbcsr_set(cmatrix, 0.0_dp)
523 CALL rskp_transform(rmatrix=rmatrix, cmatrix=cmatrix, rsmat=matrix_ks, ispin=ispin, &
524 xkp=xkp(1:3, ik), cell_to_index=cell_to_index, sab_nl=sab_nl)
525 CALL dbcsr_desymmetrize(rmatrix, tmpmat)
526 CALL copy_dbcsr_to_fm(tmpmat, fmwork(1))
527 CALL dbcsr_desymmetrize(cmatrix, tmpmat)
528 CALL copy_dbcsr_to_fm(tmpmat, fmwork(2))
529 ! s matrix is not spin dependent, double the work
530 CALL dbcsr_set(rmatrix, 0.0_dp)
531 CALL dbcsr_set(cmatrix, 0.0_dp)
532 CALL rskp_transform(rmatrix=rmatrix, cmatrix=cmatrix, rsmat=matrix_s, ispin=1, &
533 xkp=xkp(1:3, ik), cell_to_index=cell_to_index, sab_nl=sab_nl)
534 CALL dbcsr_desymmetrize(rmatrix, tmpmat)
535 CALL copy_dbcsr_to_fm(tmpmat, fmwork(3))
536 CALL dbcsr_desymmetrize(cmatrix, tmpmat)
537 CALL copy_dbcsr_to_fm(tmpmat, fmwork(4))
538 END IF
539 ! transfer to kpoint group
540 ! redistribution of matrices, new blacs environment
541 ! fmwork -> fmlocal -> rksmat/cksmat
542 ! fmwork -> fmlocal -> rsmat/csmat
543 IF (use_real_wfn) THEN
544 IF (my_kpgrp) THEN
545 CALL cp_fm_start_copy_general(fmwork(1), rksmat, para_env, info(indx, 1))
546 CALL cp_fm_start_copy_general(fmwork(3), rsmat, para_env, info(indx, 2))
547 ELSE
548 CALL cp_fm_start_copy_general(fmwork(1), fmdummy, para_env, info(indx, 1))
549 CALL cp_fm_start_copy_general(fmwork(3), fmdummy, para_env, info(indx, 2))
550 END IF
551 ELSE
552 IF (my_kpgrp) THEN
553 CALL cp_fm_start_copy_general(fmwork(1), fmlocal, para_env, info(indx, 1))
554 CALL cp_fm_start_copy_general(fmwork(2), fmlocal, para_env, info(indx, 2))
555 CALL cp_fm_start_copy_general(fmwork(3), fmlocal, para_env, info(indx, 3))
556 CALL cp_fm_start_copy_general(fmwork(4), fmlocal, para_env, info(indx, 4))
557 ELSE
558 CALL cp_fm_start_copy_general(fmwork(1), fmdummy, para_env, info(indx, 1))
559 CALL cp_fm_start_copy_general(fmwork(2), fmdummy, para_env, info(indx, 2))
560 CALL cp_fm_start_copy_general(fmwork(3), fmdummy, para_env, info(indx, 3))
561 CALL cp_fm_start_copy_general(fmwork(4), fmdummy, para_env, info(indx, 4))
562 END IF
563 END IF
564 END DO
565 END DO
566 END DO
567
568 ! Finish communication then diagonalise in each group
569 IF (do_diis) THEN
570 CALL get_qs_env(qs_env, para_env=para_env_global)
571 scf_section => section_vals_get_subs_vals(qs_env%input, "DFT%SCF")
572 CALL qs_diis_b_info_kp(kpoints%scf_diis_buffer, ib, nb)
573 indx = 0
574 DO ikp = 1, kplocal
575 DO ispin = 1, nspin
576 DO igroup = 1, nkp_groups
577 ! number of current kpoint
578 ik = kp_dist(1, igroup) + ikp - 1
579 my_kpgrp = (ik >= kpoints%kp_range(1) .AND. ik <= kpoints%kp_range(2))
580 indx = indx + 1
581 IF (my_kpgrp) THEN
582 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 1))
583 CALL cp_cfm_scale_and_add_fm(czero, cksmat, cone, fmlocal)
584 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 2))
585 CALL cp_cfm_scale_and_add_fm(cone, cksmat, ione, fmlocal)
586 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 3))
587 CALL cp_cfm_scale_and_add_fm(czero, csmat, cone, fmlocal)
588 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 4))
589 CALL cp_cfm_scale_and_add_fm(cone, csmat, ione, fmlocal)
590 END IF
591 END DO !igroup
592
593 kp => kpoints%kp_env(ikp)%kpoint_env
594 CALL qs_diis_b_calc_err_kp(kpoints%scf_diis_buffer, ib, kp%mos, cksmat, csmat, &
595 ispin, ikp, kplocal, scf_section)
596
597 END DO !ispin
598 END DO !ikp
599
600 ALLOCATE (coeffs(nb))
601 CALL qs_diis_b_step_kp(kpoints%scf_diis_buffer, coeffs, ib, nb, scf_env%iter_delta, diis_error, &
602 diis_step, scf_control%eps_diis, nspin, nkp, kplocal, scf_control%nmixing, &
603 scf_section, para_env_global)
604
605 !build the ks matrices and idagonalize
606 DO ikp = 1, kplocal
607 DO ispin = 1, nspin
608 kp => kpoints%kp_env(ikp)%kpoint_env
609 CALL cp_cfm_to_cfm(kpoints%scf_diis_buffer%smat(ikp), csmat)
610
611 CALL cp_cfm_scale(czero, cksmat)
612 DO jb = 1, nb
613 CALL cp_cfm_scale_and_add(cone, cksmat, coeffs(jb), kpoints%scf_diis_buffer%param(jb, ispin, ikp))
614 END DO
615
616 CALL get_mo_set(kp%mos(1, ispin), mo_coeff=rmos, eigenvalues=eigenvalues)
617 CALL get_mo_set(kp%mos(2, ispin), mo_coeff=imos)
618 IF (scf_env%cholesky_method == cholesky_off) THEN
619 CALL cp_cfm_geeig_canon(cksmat, csmat, cmos, eigenvalues, cwork, &
620 scf_control%eps_eigval)
621 ELSE
622 CALL cp_cfm_geeig(cksmat, csmat, cmos, eigenvalues, cwork)
623 END IF
624 ! copy eigenvalues to imag set (keep them in sync)
625 kp%mos(2, ispin)%eigenvalues = eigenvalues
626 ! split real and imaginary part of mos
627 CALL cp_cfm_to_fm(cmos, rmos, imos)
628 END DO
629 END DO
630
631 ELSE !no DIIS
632 diis_step = .false.
633 indx = 0
634 DO ikp = 1, kplocal
635 DO ispin = 1, nspin
636 DO igroup = 1, nkp_groups
637 ! number of current kpoint
638 ik = kp_dist(1, igroup) + ikp - 1
639 my_kpgrp = (ik >= kpoints%kp_range(1) .AND. ik <= kpoints%kp_range(2))
640 indx = indx + 1
641 IF (my_kpgrp) THEN
642 IF (use_real_wfn) THEN
643 CALL cp_fm_finish_copy_general(rksmat, info(indx, 1))
644 CALL cp_fm_finish_copy_general(rsmat, info(indx, 2))
645 ELSE
646 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 1))
647 CALL cp_cfm_scale_and_add_fm(czero, cksmat, cone, fmlocal)
648 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 2))
649 CALL cp_cfm_scale_and_add_fm(cone, cksmat, ione, fmlocal)
650 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 3))
651 CALL cp_cfm_scale_and_add_fm(czero, csmat, cone, fmlocal)
652 CALL cp_fm_finish_copy_general(fmlocal, info(indx, 4))
653 CALL cp_cfm_scale_and_add_fm(cone, csmat, ione, fmlocal)
654 END IF
655 END IF
656 END DO
657
658 ! Each kpoint group has now information on a kpoint to be diagonalized
659 ! General eigensolver Hermite or Symmetric
660 kp => kpoints%kp_env(ikp)%kpoint_env
661 IF (use_real_wfn) THEN
662 CALL get_mo_set(kp%mos(1, ispin), mo_coeff=mo_coeff, eigenvalues=eigenvalues)
663 IF (scf_env%cholesky_method == cholesky_off) THEN
664 CALL cp_fm_geeig_canon(rksmat, rsmat, mo_coeff, eigenvalues, fmlocal, &
665 scf_control%eps_eigval)
666 ELSE
667 CALL cp_fm_geeig(rksmat, rsmat, mo_coeff, eigenvalues, fmlocal)
668 END IF
669 ELSE
670 CALL get_mo_set(kp%mos(1, ispin), mo_coeff=rmos, eigenvalues=eigenvalues)
671 CALL get_mo_set(kp%mos(2, ispin), mo_coeff=imos)
672 IF (scf_env%cholesky_method == cholesky_off) THEN
673 CALL cp_cfm_geeig_canon(cksmat, csmat, cmos, eigenvalues, cwork, &
674 scf_control%eps_eigval)
675 ELSE
676 CALL cp_cfm_geeig(cksmat, csmat, cmos, eigenvalues, cwork)
677 END IF
678 ! copy eigenvalues to imag set (keep them in sync)
679 kp%mos(2, ispin)%eigenvalues = eigenvalues
680 ! split real and imaginary part of mos
681 CALL cp_cfm_to_fm(cmos, rmos, imos)
682 END IF
683 END DO
684 END DO
685 END IF
686
687 ! Clean up communication
688 indx = 0
689 DO ikp = 1, kplocal
690 DO ispin = 1, nspin
691 DO igroup = 1, nkp_groups
692 ! number of current kpoint
693 ik = kp_dist(1, igroup) + ikp - 1
694 my_kpgrp = (ik >= kpoints%kp_range(1) .AND. ik <= kpoints%kp_range(2))
695 indx = indx + 1
696 IF (use_real_wfn) THEN
697 CALL cp_fm_cleanup_copy_general(info(indx, 1))
698 CALL cp_fm_cleanup_copy_general(info(indx, 2))
699 ELSE
700 CALL cp_fm_cleanup_copy_general(info(indx, 1))
701 CALL cp_fm_cleanup_copy_general(info(indx, 2))
702 CALL cp_fm_cleanup_copy_general(info(indx, 3))
703 CALL cp_fm_cleanup_copy_general(info(indx, 4))
704 END IF
705 END DO
706 END DO
707 END DO
708
709 ! All done
710 DEALLOCATE (info)
711
712 IF (update_p) THEN
713 ! MO occupations
714 CALL kpoint_set_mo_occupation(kpoints, scf_control%smear)
715
716 ! density matrices
717 CALL kpoint_density_matrices(kpoints)
718 ! density matrices in real space
719 CALL kpoint_density_transform(kpoints, scf_env%p_mix_new, .false., &
720 matrix_s(1, 1)%matrix, sab_nl, fmwork)
721 END IF
722
723 CALL dbcsr_deallocate_matrix(rmatrix)
724 CALL dbcsr_deallocate_matrix(cmatrix)
725 CALL dbcsr_deallocate_matrix(tmpmat)
726
727 IF (use_real_wfn) THEN
728 CALL cp_fm_release(rksmat)
729 CALL cp_fm_release(rsmat)
730 ELSE
731 CALL cp_cfm_release(cksmat)
732 CALL cp_cfm_release(csmat)
733 CALL cp_cfm_release(cwork)
734 CALL cp_cfm_release(cmos)
735 END IF
736 CALL fm_pool_give_back_fm(ao_ao_fm_pools(1)%pool, fmlocal)
737
738 CALL timestop(handle)
739
740 END SUBROUTINE do_general_diag_kp
741
742! **************************************************************************************************
743!> \brief inner loop within MOS subspace, to refine occupation and density,
744!> before next diagonalization of the Hamiltonian
745!> \param qs_env ...
746!> \param scf_env ...
747!> \param subspace_env ...
748!> \param mos ...
749!> \param rho ...
750!> \param ks_env ...
751!> \param scf_section ...
752!> \param scf_control ...
753!> \par History
754!> 09.2009 created [MI]
755!> \note it is assumed that when diagonalization is used, also some mixing procedure is active
756! **************************************************************************************************
757 SUBROUTINE do_scf_diag_subspace(qs_env, scf_env, subspace_env, mos, rho, &
758 ks_env, scf_section, scf_control)
759
760 TYPE(qs_environment_type), POINTER :: qs_env
761 TYPE(qs_scf_env_type), POINTER :: scf_env
762 TYPE(subspace_env_type), POINTER :: subspace_env
763 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
764 TYPE(qs_rho_type), POINTER :: rho
765 TYPE(qs_ks_env_type), POINTER :: ks_env
766 TYPE(section_vals_type), POINTER :: scf_section
767 TYPE(scf_control_type), POINTER :: scf_control
768
769 CHARACTER(LEN=*), PARAMETER :: routinen = 'do_scf_diag_subspace'
770 REAL(kind=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
771
772 INTEGER :: handle, i, iloop, ispin, nao, nmo, &
773 nspin, output_unit
774 LOGICAL :: converged
775 REAL(dp) :: ene_diff, ene_old, iter_delta, max_val, &
776 sum_band, sum_val, t1, t2
777 REAL(kind=dp), DIMENSION(:), POINTER :: mo_eigenvalues, mo_occupations
778 TYPE(cp_1d_r_p_type), ALLOCATABLE, DIMENSION(:) :: eval_first, occ_first
779 TYPE(cp_fm_type) :: work
780 TYPE(cp_fm_type), POINTER :: c0, chc, evec, mo_coeff
781 TYPE(cp_logger_type), POINTER :: logger
782 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s, rho_ao
783 TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: rho_ao_kp
784 TYPE(dft_control_type), POINTER :: dft_control
785 TYPE(mp_para_env_type), POINTER :: para_env
786 TYPE(qs_energy_type), POINTER :: energy
787 TYPE(rho_atom_type), DIMENSION(:), POINTER :: rho_atom
788
789 CALL timeset(routinen, handle)
790 NULLIFY (c0, chc, energy, evec, matrix_ks, mo_coeff, mo_eigenvalues, &
791 mo_occupations, dft_control, rho_ao, rho_ao_kp)
792
793 logger => cp_get_default_logger()
794 output_unit = cp_print_key_unit_nr(logger, scf_section, "PRINT%DIAG_SUB_SCF", &
795 extension=".scfLog")
796
797 !Extra loop keeping mos unchanged and refining the subspace occupation
798 nspin = SIZE(mos)
799 CALL qs_rho_get(rho, rho_ao=rho_ao, rho_ao_kp=rho_ao_kp)
800
801 ALLOCATE (eval_first(nspin))
802 ALLOCATE (occ_first(nspin))
803 DO ispin = 1, nspin
804 CALL get_mo_set(mo_set=mos(ispin), &
805 nmo=nmo, &
806 eigenvalues=mo_eigenvalues, &
807 occupation_numbers=mo_occupations)
808 ALLOCATE (eval_first(ispin)%array(nmo))
809 ALLOCATE (occ_first(ispin)%array(nmo))
810 eval_first(ispin)%array(1:nmo) = mo_eigenvalues(1:nmo)
811 occ_first(ispin)%array(1:nmo) = mo_occupations(1:nmo)
812 END DO
813
814 DO ispin = 1, nspin
815 ! does not yet handle k-points
816 CALL dbcsr_copy(subspace_env%p_matrix_store(ispin)%matrix, rho_ao(ispin)%matrix)
817 CALL dbcsr_copy(rho_ao(ispin)%matrix, scf_env%p_mix_new(ispin, 1)%matrix)
818 END DO
819
820 subspace_env%p_matrix_mix => scf_env%p_mix_new
821
822 NULLIFY (matrix_ks, energy, para_env, matrix_s)
823 CALL get_qs_env(qs_env, &
824 matrix_ks=matrix_ks, &
825 energy=energy, &
826 matrix_s=matrix_s, &
827 para_env=para_env, &
828 dft_control=dft_control)
829
830 ! mixing storage allocation
831 IF (subspace_env%mixing_method >= gspace_mixing_nr) THEN
832 CALL mixing_allocate(qs_env, subspace_env%mixing_method, scf_env%p_mix_new, &
833 scf_env%p_delta, nspin, subspace_env%mixing_store)
834 IF (dft_control%qs_control%gapw) THEN
835 CALL get_qs_env(qs_env=qs_env, rho_atom_set=rho_atom)
836 CALL mixing_init(subspace_env%mixing_method, rho, subspace_env%mixing_store, &
837 para_env, rho_atom=rho_atom)
838 ELSEIF (dft_control%qs_control%dftb .OR. dft_control%qs_control%xtb) THEN
839 CALL charge_mixing_init(subspace_env%mixing_store)
840 ELSEIF (dft_control%qs_control%semi_empirical) THEN
841 cpabort('SE Code not possible')
842 ELSE
843 CALL mixing_init(subspace_env%mixing_method, rho, subspace_env%mixing_store, para_env)
844 END IF
845 END IF
846
847 ene_old = 0.0_dp
848 ene_diff = 0.0_dp
849 IF (output_unit > 0) THEN
850 WRITE (output_unit, "(/T19,A)") '<<<<<<<<< SUBSPACE ROTATION <<<<<<<<<<'
851 WRITE (output_unit, "(T4,A,T13,A,T21,A,T38,A,T51,A,T65,A/,T4,A)") &
852 "In-step", "Time", "Convergence", "Band ene.", "Total ene.", "Energy diff.", repeat("-", 74)
853 END IF
854
855 ! recalculate density matrix here
856
857 ! update of density
858 CALL qs_rho_update_rho(rho, qs_env=qs_env)
859
860 DO iloop = 1, subspace_env%max_iter
861 t1 = m_walltime()
862 converged = .false.
863 ene_old = energy%total
864
865 CALL qs_ks_did_change(ks_env, rho_changed=.true.)
866 CALL qs_ks_update_qs_env(qs_env, calculate_forces=.false., &
867 just_energy=.false., print_active=.false.)
868
869 max_val = 0.0_dp
870 sum_val = 0.0_dp
871 sum_band = 0.0_dp
872 DO ispin = 1, SIZE(matrix_ks)
873 CALL get_mo_set(mo_set=mos(ispin), &
874 nao=nao, &
875 nmo=nmo, &
876 eigenvalues=mo_eigenvalues, &
877 occupation_numbers=mo_occupations, &
878 mo_coeff=mo_coeff)
879
880 !compute C'HC
881 chc => subspace_env%chc_mat(ispin)
882 evec => subspace_env%c_vec(ispin)
883 c0 => subspace_env%c0(ispin)
884 CALL cp_fm_to_fm(mo_coeff, c0)
885 CALL cp_fm_create(work, c0%matrix_struct)
886 CALL cp_dbcsr_sm_fm_multiply(matrix_ks(ispin)%matrix, c0, work, nmo)
887 CALL parallel_gemm('T', 'N', nmo, nmo, nao, rone, c0, work, rzero, chc)
888 CALL cp_fm_release(work)
889 !diagonalize C'HC
890 CALL choose_eigv_solver(chc, evec, mo_eigenvalues)
891
892 !rotate the mos by the eigenvectors of C'HC
893 CALL parallel_gemm('N', 'N', nao, nmo, nmo, rone, c0, evec, rzero, mo_coeff)
894
895 CALL set_mo_occupation(mo_set=mos(ispin), &
896 smear=scf_control%smear)
897
898 ! does not yet handle k-points
899 CALL calculate_density_matrix(mos(ispin), &
900 subspace_env%p_matrix_mix(ispin, 1)%matrix)
901
902 DO i = 1, nmo
903 sum_band = sum_band + mo_eigenvalues(i)*mo_occupations(i)
904 END DO
905
906 !check for self consistency
907 END DO
908
909 IF (subspace_env%mixing_method == direct_mixing_nr) THEN
910 CALL scf_env_density_mixing(subspace_env%p_matrix_mix, &
911 scf_env%mixing_store, rho_ao_kp, para_env, iter_delta, iloop)
912 ELSE
913 CALL self_consistency_check(rho_ao_kp, scf_env%p_delta, para_env, &
914 subspace_env%p_matrix_mix, delta=iter_delta)
915 END IF
916
917 DO ispin = 1, nspin
918 ! does not yet handle k-points
919 CALL dbcsr_copy(rho_ao(ispin)%matrix, subspace_env%p_matrix_mix(ispin, 1)%matrix)
920 END DO
921 ! update of density
922 CALL qs_rho_update_rho(rho, qs_env=qs_env)
923 ! Mixing in reciprocal space
924 IF (subspace_env%mixing_method >= gspace_mixing_nr) THEN
925 CALL gspace_mixing(qs_env, scf_env%mixing_method, subspace_env%mixing_store, &
926 rho, para_env, scf_env%iter_count)
927 END IF
928
929 ene_diff = energy%total - ene_old
930 converged = (abs(ene_diff) < subspace_env%eps_ene .AND. &
931 iter_delta < subspace_env%eps_adapt*scf_env%iter_delta)
932 t2 = m_walltime()
933 IF (output_unit > 0) THEN
934 WRITE (output_unit, "(T4,I5,T11,F8.3,T18,E14.4,T34,F12.5,T46,F16.8,T62,E14.4)") &
935 iloop, t2 - t1, iter_delta, sum_band, energy%total, ene_diff
936 CALL m_flush(output_unit)
937 END IF
938 IF (converged) THEN
939 IF (output_unit > 0) WRITE (output_unit, "(T10,A,I6,A,/)") &
940 " Reached convergence in ", iloop, " iterations "
941 EXIT
942 END IF
943
944 END DO ! iloop
945
946 NULLIFY (subspace_env%p_matrix_mix)
947 DO ispin = 1, nspin
948 ! does not yet handle k-points
949 CALL dbcsr_copy(scf_env%p_mix_new(ispin, 1)%matrix, rho_ao(ispin)%matrix)
950 CALL dbcsr_copy(rho_ao(ispin)%matrix, subspace_env%p_matrix_store(ispin)%matrix)
951
952 DEALLOCATE (eval_first(ispin)%array, occ_first(ispin)%array)
953 END DO
954 DEALLOCATE (eval_first, occ_first)
955
956 CALL timestop(handle)
957
958 END SUBROUTINE do_scf_diag_subspace
959
960! **************************************************************************************************
961!> \brief ...
962!> \param subspace_env ...
963!> \param qs_env ...
964!> \param mos ...
965! **************************************************************************************************
966 SUBROUTINE diag_subspace_allocate(subspace_env, qs_env, mos)
967
968 TYPE(subspace_env_type), POINTER :: subspace_env
969 TYPE(qs_environment_type), POINTER :: qs_env
970 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mos
971
972 CHARACTER(LEN=*), PARAMETER :: routinen = 'diag_subspace_allocate'
973
974 INTEGER :: handle, i, ispin, nmo, nspin
975 TYPE(cp_fm_struct_type), POINTER :: fm_struct_tmp
976 TYPE(cp_fm_type), POINTER :: mo_coeff
977 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_s
978 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
979 POINTER :: sab_orb
980
981 CALL timeset(routinen, handle)
982
983 NULLIFY (sab_orb, matrix_s)
984 CALL get_qs_env(qs_env=qs_env, sab_orb=sab_orb, &
985 matrix_s=matrix_s)
986
987 nspin = SIZE(mos)
988! *** allocate p_atrix_store ***
989 IF (.NOT. ASSOCIATED(subspace_env%p_matrix_store)) THEN
990 CALL dbcsr_allocate_matrix_set(subspace_env%p_matrix_store, nspin)
991
992 DO i = 1, nspin
993 ALLOCATE (subspace_env%p_matrix_store(i)%matrix)
994 CALL dbcsr_create(matrix=subspace_env%p_matrix_store(i)%matrix, template=matrix_s(1)%matrix, &
995 name="DENSITY_STORE", matrix_type=dbcsr_type_symmetric)
996 CALL cp_dbcsr_alloc_block_from_nbl(subspace_env%p_matrix_store(i)%matrix, &
997 sab_orb)
998 CALL dbcsr_set(subspace_env%p_matrix_store(i)%matrix, 0.0_dp)
999 END DO
1000
1001 END IF
1002
1003 ALLOCATE (subspace_env%chc_mat(nspin))
1004 ALLOCATE (subspace_env%c_vec(nspin))
1005 ALLOCATE (subspace_env%c0(nspin))
1006
1007 DO ispin = 1, nspin
1008 CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff, nmo=nmo)
1009 CALL cp_fm_create(subspace_env%c0(ispin), mo_coeff%matrix_struct)
1010 NULLIFY (fm_struct_tmp)
1011 CALL cp_fm_struct_create(fm_struct_tmp, nrow_global=nmo, ncol_global=nmo, &
1012 para_env=mo_coeff%matrix_struct%para_env, &
1013 context=mo_coeff%matrix_struct%context)
1014 CALL cp_fm_create(subspace_env%chc_mat(ispin), fm_struct_tmp, "chc")
1015 CALL cp_fm_create(subspace_env%c_vec(ispin), fm_struct_tmp, "vec")
1016 CALL cp_fm_struct_release(fm_struct_tmp)
1017 END DO
1018
1019 CALL timestop(handle)
1020
1021 END SUBROUTINE diag_subspace_allocate
1022
1023! **************************************************************************************************
1024!> \brief the inner loop of scf, specific to diagonalization without S matrix
1025!> basically, in goes the ks matrix out goes a new p matrix
1026!> \param scf_env ...
1027!> \param mos ...
1028!> \param matrix_ks ...
1029!> \param scf_control ...
1030!> \param scf_section ...
1031!> \param diis_step ...
1032!> \par History
1033!> 03.2006 created [Joost VandeVondele]
1034! **************************************************************************************************
1035 SUBROUTINE do_special_diag(scf_env, mos, matrix_ks, scf_control, &
1036 scf_section, diis_step)
1037
1038 TYPE(qs_scf_env_type), POINTER :: scf_env
1039 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
1040 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks
1041 TYPE(scf_control_type), POINTER :: scf_control
1042 TYPE(section_vals_type), POINTER :: scf_section
1043 LOGICAL, INTENT(INOUT) :: diis_step
1044
1045 INTEGER :: ispin, nspin
1046 LOGICAL :: do_level_shift, use_jacobi
1047 REAL(kind=dp) :: diis_error
1048
1049 nspin = SIZE(matrix_ks)
1050
1051 DO ispin = 1, nspin
1052 CALL copy_dbcsr_to_fm(matrix_ks(ispin)%matrix, scf_env%scf_work1(ispin))
1053 END DO
1054 IF (scf_env%iter_count > 1 .AND. .NOT. scf_env%skip_diis) THEN
1055 CALL qs_diis_b_step(scf_env%scf_diis_buffer, mos, scf_env%scf_work1, &
1056 scf_env%scf_work2, scf_env%iter_delta, diis_error, diis_step, &
1057 scf_control%eps_diis, scf_control%nmixing, &
1058 scf_section=scf_section)
1059 ELSE
1060 diis_step = .false.
1061 END IF
1062
1063 IF ((scf_env%iter_count > 1) .AND. (scf_env%iter_delta < scf_control%diagonalization%eps_jacobi)) THEN
1064 use_jacobi = .true.
1065 ELSE
1066 use_jacobi = .false.
1067 END IF
1068
1069 do_level_shift = ((scf_control%level_shift /= 0.0_dp) .AND. &
1070 ((scf_control%density_guess == core_guess) .OR. (scf_env%iter_count > 1)))
1071 IF (diis_step) THEN
1072 scf_env%iter_param = diis_error
1073 IF (use_jacobi) THEN
1074 scf_env%iter_method = "DIIS/Jacobi"
1075 ELSE
1076 scf_env%iter_method = "DIIS/Diag."
1077 END IF
1078 ELSE
1079 IF (scf_env%mixing_method == 1) THEN
1080 scf_env%iter_param = scf_env%p_mix_alpha
1081 IF (use_jacobi) THEN
1082 scf_env%iter_method = "P_Mix/Jacobi"
1083 ELSE
1084 scf_env%iter_method = "P_Mix/Diag."
1085 END IF
1086 ELSEIF (scf_env%mixing_method > 1) THEN
1087 scf_env%iter_param = scf_env%mixing_store%alpha
1088 IF (use_jacobi) THEN
1089 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Jacobi"
1090 ELSE
1091 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Diag."
1092 END IF
1093 END IF
1094 END IF
1095 scf_env%iter_delta = 0.0_dp
1096
1097 DO ispin = 1, nspin
1098 CALL eigensolver_simple(matrix_ks=scf_env%scf_work1(ispin), &
1099 mo_set=mos(ispin), &
1100 work=scf_env%scf_work2, &
1101 do_level_shift=do_level_shift, &
1102 level_shift=scf_control%level_shift, &
1103 use_jacobi=use_jacobi, &
1104 jacobi_threshold=scf_control%diagonalization%jacobi_threshold)
1105 END DO
1106
1107 CALL set_mo_occupation(mo_array=mos, &
1108 smear=scf_control%smear)
1109
1110 DO ispin = 1, nspin
1111 ! does not yet handle k-points
1112 CALL calculate_density_matrix(mos(ispin), &
1113 scf_env%p_mix_new(ispin, 1)%matrix)
1114 END DO
1115
1116 END SUBROUTINE do_special_diag
1117
1118! **************************************************************************************************
1119!> \brief the inner loop of scf, specific to iterative diagonalization using OT
1120!> with S matrix; basically, in goes the ks matrix out goes a new p matrix
1121!> \param scf_env ...
1122!> \param mos ...
1123!> \param matrix_ks ...
1124!> \param matrix_s ...
1125!> \param scf_control ...
1126!> \param scf_section ...
1127!> \param diis_step ...
1128!> \par History
1129!> 10.2008 created [JGH]
1130! **************************************************************************************************
1131 SUBROUTINE do_ot_diag(scf_env, mos, matrix_ks, matrix_s, &
1132 scf_control, scf_section, diis_step)
1133
1134 TYPE(qs_scf_env_type), POINTER :: scf_env
1135 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
1136 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
1137 TYPE(scf_control_type), POINTER :: scf_control
1138 TYPE(section_vals_type), POINTER :: scf_section
1139 LOGICAL, INTENT(INOUT) :: diis_step
1140
1141 INTEGER :: homo, ispin, nmo, nspin
1142 REAL(kind=dp) :: diis_error, eps_iter
1143 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues
1144 TYPE(cp_fm_type), POINTER :: mo_coeff
1145
1146 NULLIFY (eigenvalues)
1147
1148 nspin = SIZE(matrix_ks)
1149
1150 DO ispin = 1, nspin
1151 CALL copy_dbcsr_to_fm(matrix_ks(ispin)%matrix, &
1152 scf_env%scf_work1(ispin))
1153 END DO
1154
1155 IF ((scf_env%iter_count > 1) .AND. (.NOT. scf_env%skip_diis)) THEN
1156 CALL qs_diis_b_step(scf_env%scf_diis_buffer, mos, scf_env%scf_work1, &
1157 scf_env%scf_work2, scf_env%iter_delta, diis_error, diis_step, &
1158 scf_control%eps_diis, scf_control%nmixing, &
1159 s_matrix=matrix_s, &
1160 scf_section=scf_section)
1161 ELSE
1162 diis_step = .false.
1163 END IF
1164
1165 eps_iter = scf_control%diagonalization%eps_iter
1166 IF (diis_step) THEN
1167 scf_env%iter_param = diis_error
1168 scf_env%iter_method = "DIIS/OTdiag"
1169 DO ispin = 1, nspin
1170 CALL copy_fm_to_dbcsr(scf_env%scf_work1(ispin), &
1171 matrix_ks(ispin)%matrix, keep_sparsity=.true.)
1172 END DO
1173 eps_iter = max(eps_iter, scf_control%diagonalization%eps_adapt*diis_error)
1174 ELSE
1175 IF (scf_env%mixing_method == 1) THEN
1176 scf_env%iter_param = scf_env%p_mix_alpha
1177 scf_env%iter_method = "P_Mix/OTdiag."
1178 ELSEIF (scf_env%mixing_method > 1) THEN
1179 scf_env%iter_param = scf_env%mixing_store%alpha
1180 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/OTdiag."
1181 END IF
1182 END IF
1183
1184 scf_env%iter_delta = 0.0_dp
1185
1186 DO ispin = 1, nspin
1187 CALL get_mo_set(mos(ispin), &
1188 mo_coeff=mo_coeff, &
1189 eigenvalues=eigenvalues, &
1190 nmo=nmo, &
1191 homo=homo)
1192 CALL ot_eigensolver(matrix_h=matrix_ks(ispin)%matrix, &
1193 matrix_s=matrix_s(1)%matrix, &
1194 matrix_c_fm=mo_coeff, &
1195 preconditioner=scf_env%ot_preconditioner(1)%preconditioner, &
1196 eps_gradient=eps_iter, &
1197 iter_max=scf_control%diagonalization%max_iter, &
1198 silent=.true., &
1199 ot_settings=scf_control%diagonalization%ot_settings)
1200 CALL calculate_subspace_eigenvalues(mo_coeff, matrix_ks(ispin)%matrix, &
1201 evals_arg=eigenvalues, &
1202 do_rotation=.true.)
1203 CALL copy_fm_to_dbcsr(mos(ispin)%mo_coeff, &
1204 mos(ispin)%mo_coeff_b)
1205 !fm->dbcsr
1206 END DO
1207
1208 CALL set_mo_occupation(mo_array=mos, &
1209 smear=scf_control%smear)
1210
1211 DO ispin = 1, nspin
1212 ! does not yet handle k-points
1213 CALL calculate_density_matrix(mos(ispin), &
1214 scf_env%p_mix_new(ispin, 1)%matrix)
1215 END DO
1216
1217 END SUBROUTINE do_ot_diag
1218
1219! **************************************************************************************************
1220!> \brief Solve a set restricted open Kohn-Sham (ROKS) equations based on the
1221!> alpha and beta Kohn-Sham matrices from unrestricted Kohn-Sham.
1222!> \param scf_env ...
1223!> \param mos ...
1224!> \param matrix_ks ...
1225!> \param matrix_s ...
1226!> \param scf_control ...
1227!> \param scf_section ...
1228!> \param diis_step ...
1229!> \param orthogonal_basis ...
1230!> \par History
1231!> 04.2006 created [MK]
1232!> Revised (01.05.06,MK)
1233!> \note
1234!> this is only a high-spin ROKS.
1235! **************************************************************************************************
1236 SUBROUTINE do_roks_diag(scf_env, mos, matrix_ks, matrix_s, &
1237 scf_control, scf_section, diis_step, &
1238 orthogonal_basis)
1239
1240 ! Literature: - C. C. J. Roothaan, Rev. Mod. Phys. 32, 179 (1960)
1241 ! - M. F. Guest and V. R. Saunders, Mol. Phys. 28(3), 819 (1974)
1242 ! - M. Filatov and S. Shaik, Chem. Phys. Lett. 288, 689 (1998)
1243
1244 TYPE(qs_scf_env_type), POINTER :: scf_env
1245 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mos
1246 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
1247 TYPE(scf_control_type), POINTER :: scf_control
1248 TYPE(section_vals_type), POINTER :: scf_section
1249 LOGICAL, INTENT(INOUT) :: diis_step
1250 LOGICAL, INTENT(IN) :: orthogonal_basis
1251
1252 CHARACTER(LEN=*), PARAMETER :: routinen = 'do_roks_diag'
1253
1254 INTEGER :: handle, homoa, homob, imo, nalpha, nao, &
1255 nbeta, nmo
1256 REAL(kind=dp) :: diis_error, level_shift_loc
1257 REAL(kind=dp), DIMENSION(:), POINTER :: eiga, eigb, occa, occb
1258 TYPE(cp_fm_type), POINTER :: ksa, ksb, mo2ao, moa, mob, ortho, work
1259
1260! -------------------------------------------------------------------------
1261
1262 CALL timeset(routinen, handle)
1263
1264 IF (scf_env%cholesky_method == cholesky_inverse) THEN
1265 ortho => scf_env%ortho_m1
1266 ELSE
1267 ortho => scf_env%ortho
1268 END IF
1269 work => scf_env%scf_work2
1270
1271 ksa => scf_env%scf_work1(1)
1272 ksb => scf_env%scf_work1(2)
1273
1274 CALL copy_dbcsr_to_fm(matrix_ks(1)%matrix, ksa)
1275 CALL copy_dbcsr_to_fm(matrix_ks(2)%matrix, ksb)
1276
1277 ! Get MO information
1278
1279 CALL get_mo_set(mo_set=mos(1), &
1280 nao=nao, &
1281 nmo=nmo, &
1282 nelectron=nalpha, &
1283 homo=homoa, &
1284 eigenvalues=eiga, &
1285 occupation_numbers=occa, &
1286 mo_coeff=moa)
1287
1288 CALL get_mo_set(mo_set=mos(2), &
1289 nelectron=nbeta, &
1290 homo=homob, &
1291 eigenvalues=eigb, &
1292 occupation_numbers=occb, &
1293 mo_coeff=mob)
1294
1295 ! Define the amount of level-shifting
1296
1297 IF ((scf_control%level_shift /= 0.0_dp) .AND. &
1298 ((scf_control%density_guess == core_guess) .OR. &
1299 (scf_control%density_guess == restart_guess) .OR. &
1300 (scf_env%iter_count > 1))) THEN
1301 level_shift_loc = scf_control%level_shift
1302 ELSE
1303 level_shift_loc = 0.0_dp
1304 END IF
1305
1306 IF ((scf_env%iter_count > 1) .OR. &
1307 (scf_control%density_guess == core_guess) .OR. &
1308 (scf_control%density_guess == restart_guess)) THEN
1309
1310 ! Transform the spin unrestricted alpha and beta Kohn-Sham matrices
1311 ! from AO basis to MO basis: K(MO) = C(T)*K(AO)*C
1312
1313 CALL cp_fm_symm("L", "U", nao, nao, 1.0_dp, ksa, moa, 0.0_dp, work)
1314 CALL parallel_gemm("T", "N", nao, nao, nao, 1.0_dp, moa, work, 0.0_dp, ksa)
1315
1316 CALL cp_fm_symm("L", "U", nao, nao, 1.0_dp, ksb, moa, 0.0_dp, work)
1317 CALL parallel_gemm("T", "N", nao, nao, nao, 1.0_dp, moa, work, 0.0_dp, ksb)
1318
1319 ! Combine the spin unrestricted alpha and beta Kohn-Sham matrices
1320 ! in the MO basis
1321
1322 IF (scf_control%roks_scheme == general_roks) THEN
1323 CALL combine_ks_matrices(ksa, ksb, occa, occb, scf_control%roks_f, &
1324 nalpha, nbeta)
1325 ELSE IF (scf_control%roks_scheme == high_spin_roks) THEN
1326 CALL combine_ks_matrices(ksa, ksb, occa, occb, scf_control%roks_parameter)
1327 ELSE
1328 cpabort("Unknown ROKS scheme requested")
1329 END IF
1330
1331 ! Back-transform the restricted open Kohn-Sham matrix from MO basis
1332 ! to AO basis
1333
1334 IF (orthogonal_basis) THEN
1335 ! Q = C
1336 mo2ao => moa
1337 ELSE
1338 ! Q = S*C
1339 mo2ao => mob
1340!MK CALL copy_sm_to_fm(matrix_s(1)%matrix,work)
1341!MK CALL cp_fm_symm("L", "U",nao, nao, 1.0_dp, work, moa, 0.0_dp, mo2ao)
1342 CALL cp_dbcsr_sm_fm_multiply(matrix_s(1)%matrix, moa, mo2ao, nao)
1343 END IF
1344
1345 ! K(AO) = Q*K(MO)*Q(T)
1346
1347 CALL parallel_gemm("N", "T", nao, nao, nao, 1.0_dp, ksa, mo2ao, 0.0_dp, work)
1348 CALL parallel_gemm("N", "N", nao, nao, nao, 1.0_dp, mo2ao, work, 0.0_dp, ksa)
1349
1350 ELSE
1351
1352 ! No transformation matrix available, yet. The closed shell part,
1353 ! i.e. the beta Kohn-Sham matrix in AO basis, is taken.
1354 ! There might be better choices, anyhow.
1355
1356 CALL cp_fm_to_fm(ksb, ksa)
1357
1358 END IF
1359
1360 ! Update DIIS buffer and possibly perform DIIS extrapolation step
1361
1362 IF (scf_env%iter_count > 1) THEN
1363 IF (orthogonal_basis) THEN
1364 CALL qs_diis_b_step(diis_buffer=scf_env%scf_diis_buffer, &
1365 mo_array=mos, &
1366 kc=scf_env%scf_work1, &
1367 sc=work, &
1368 delta=scf_env%iter_delta, &
1369 error_max=diis_error, &
1370 diis_step=diis_step, &
1371 eps_diis=scf_control%eps_diis, &
1372 scf_section=scf_section, &
1373 roks=.true.)
1374 cpassert(scf_env%iter_delta == scf_env%iter_delta)
1375 ELSE
1376 CALL qs_diis_b_step(diis_buffer=scf_env%scf_diis_buffer, &
1377 mo_array=mos, &
1378 kc=scf_env%scf_work1, &
1379 sc=work, &
1380 delta=scf_env%iter_delta, &
1381 error_max=diis_error, &
1382 diis_step=diis_step, &
1383 eps_diis=scf_control%eps_diis, &
1384 scf_section=scf_section, &
1385 s_matrix=matrix_s, &
1386 roks=.true.)
1387 END IF
1388 END IF
1389
1390 IF (diis_step) THEN
1391 scf_env%iter_param = diis_error
1392 scf_env%iter_method = "DIIS/Diag."
1393 ELSE
1394 IF (scf_env%mixing_method == 1) THEN
1395 scf_env%iter_param = scf_env%p_mix_alpha
1396 scf_env%iter_method = "P_Mix/Diag."
1397 ELSEIF (scf_env%mixing_method > 1) THEN
1398 scf_env%iter_param = scf_env%mixing_store%alpha
1399 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Diag."
1400 END IF
1401 END IF
1402
1403 scf_env%iter_delta = 0.0_dp
1404
1405 IF (level_shift_loc /= 0.0_dp) THEN
1406
1407 ! Transform the current Kohn-Sham matrix from AO to MO basis
1408 ! for level-shifting using the current MO set
1409
1410 CALL cp_fm_symm("L", "U", nao, nao, 1.0_dp, ksa, moa, 0.0_dp, work)
1411 CALL parallel_gemm("T", "N", nao, nao, nao, 1.0_dp, moa, work, 0.0_dp, ksa)
1412
1413 ! Apply level-shifting using 50:50 split of the shift (could be relaxed)
1414
1415 DO imo = homob + 1, homoa
1416 CALL cp_fm_add_to_element(ksa, imo, imo, 0.5_dp*level_shift_loc)
1417 END DO
1418 DO imo = homoa + 1, nmo
1419 CALL cp_fm_add_to_element(ksa, imo, imo, level_shift_loc)
1420 END DO
1421
1422 ELSE IF (.NOT. orthogonal_basis) THEN
1423
1424 ! Transform the current Kohn-Sham matrix to an orthogonal basis
1425 SELECT CASE (scf_env%cholesky_method)
1426 CASE (cholesky_reduce)
1427 CALL cp_fm_cholesky_reduce(ksa, ortho)
1428 CASE (cholesky_restore)
1429 CALL cp_fm_uplo_to_full(ksa, work)
1430 CALL cp_fm_cholesky_restore(ksa, nao, ortho, work, &
1431 "SOLVE", pos="RIGHT")
1432 CALL cp_fm_cholesky_restore(work, nao, ortho, ksa, &
1433 "SOLVE", pos="LEFT", transa="T")
1434 CASE (cholesky_inverse)
1435 CALL cp_fm_uplo_to_full(ksa, work)
1436 CALL cp_fm_cholesky_restore(ksa, nao, ortho, work, &
1437 "MULTIPLY", pos="RIGHT")
1438 CALL cp_fm_cholesky_restore(work, nao, ortho, ksa, &
1439 "MULTIPLY", pos="LEFT", transa="T")
1440 CASE (cholesky_off)
1441 CALL cp_fm_symm("L", "U", nao, nao, 1.0_dp, ksa, ortho, 0.0_dp, work)
1442 CALL parallel_gemm("N", "N", nao, nao, nao, 1.0_dp, ortho, work, 0.0_dp, ksa)
1443 END SELECT
1444
1445 END IF
1446
1447 ! Diagonalization of the ROKS operator matrix
1448
1449 CALL choose_eigv_solver(ksa, work, eiga)
1450
1451 ! Back-transformation of the orthonormal eigenvectors if needed
1452
1453 IF (level_shift_loc /= 0.0_dp) THEN
1454 ! Use old MO set for back-transformation if level-shifting was applied
1455 CALL cp_fm_to_fm(moa, ortho)
1456 CALL parallel_gemm("N", "N", nao, nmo, nao, 1.0_dp, ortho, work, 0.0_dp, moa)
1457 ELSE
1458 IF (orthogonal_basis) THEN
1459 CALL cp_fm_to_fm(work, moa)
1460 ELSE
1461 SELECT CASE (scf_env%cholesky_method)
1463 CALL cp_fm_cholesky_restore(work, nmo, ortho, moa, "SOLVE")
1464 CASE (cholesky_inverse)
1465 CALL cp_fm_cholesky_restore(work, nmo, ortho, moa, "MULTIPLY")
1466 CASE (cholesky_off)
1467 CALL parallel_gemm("N", "N", nao, nmo, nao, 1.0_dp, ortho, work, 0.0_dp, moa)
1468 END SELECT
1469 END IF
1470 END IF
1471
1472 ! Correct MO eigenvalues, if level-shifting was applied
1473
1474 IF (level_shift_loc /= 0.0_dp) THEN
1475 DO imo = homob + 1, homoa
1476 eiga(imo) = eiga(imo) - 0.5_dp*level_shift_loc
1477 END DO
1478 DO imo = homoa + 1, nmo
1479 eiga(imo) = eiga(imo) - level_shift_loc
1480 END DO
1481 END IF
1482
1483 ! Update also the beta MO set
1484
1485 eigb(:) = eiga(:)
1486 CALL cp_fm_to_fm(moa, mob)
1487
1488 ! Calculate the new alpha and beta density matrix
1489
1490 ! does not yet handle k-points
1491 CALL calculate_density_matrix(mos(1), scf_env%p_mix_new(1, 1)%matrix)
1492 CALL calculate_density_matrix(mos(2), scf_env%p_mix_new(2, 1)%matrix)
1493
1494 CALL timestop(handle)
1495
1496 END SUBROUTINE do_roks_diag
1497
1498! **************************************************************************************************
1499!> \brief iterative diagonalization using the block Krylov-space approach
1500!> \param scf_env ...
1501!> \param mos ...
1502!> \param matrix_ks ...
1503!> \param scf_control ...
1504!> \param scf_section ...
1505!> \param check_moconv_only ...
1506!> \param
1507!> \par History
1508!> 05.2009 created [MI]
1509! **************************************************************************************************
1510
1511 SUBROUTINE do_block_krylov_diag(scf_env, mos, matrix_ks, &
1512 scf_control, scf_section, check_moconv_only)
1513
1514 TYPE(qs_scf_env_type), POINTER :: scf_env
1515 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
1516 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks
1517 TYPE(scf_control_type), POINTER :: scf_control
1518 TYPE(section_vals_type), POINTER :: scf_section
1519 LOGICAL, INTENT(IN), OPTIONAL :: check_moconv_only
1520
1521 CHARACTER(LEN=*), PARAMETER :: routinen = 'do_block_krylov_diag'
1522 REAL(kind=dp), PARAMETER :: rone = 1.0_dp, rzero = 0.0_dp
1523
1524 INTEGER :: handle, homo, ispin, iter, nao, nmo, &
1525 output_unit
1526 LOGICAL :: converged, my_check_moconv_only
1527 REAL(dp) :: eps_iter, t1, t2
1528 REAL(kind=dp), DIMENSION(:), POINTER :: mo_eigenvalues
1529 TYPE(cp_fm_type), POINTER :: c0, c1, chc, evec, ks, mo_coeff, ortho, &
1530 work
1531 TYPE(cp_logger_type), POINTER :: logger
1532
1533 logger => cp_get_default_logger()
1534 CALL timeset(routinen, handle)
1535
1536 output_unit = cp_print_key_unit_nr(logger, scf_section, "PRINT%LANCZOS", &
1537 extension=".scfLog")
1538
1539 my_check_moconv_only = .false.
1540 IF (PRESENT(check_moconv_only)) my_check_moconv_only = check_moconv_only
1541
1542 NULLIFY (mo_coeff, ortho, work, ks)
1543 NULLIFY (mo_eigenvalues)
1544 NULLIFY (c0, c1)
1545
1546 IF (scf_env%cholesky_method == cholesky_inverse) THEN
1547 ortho => scf_env%ortho_m1
1548 ELSE
1549 ortho => scf_env%ortho
1550 END IF
1551 work => scf_env%scf_work2
1552
1553 DO ispin = 1, SIZE(matrix_ks)
1554 CALL copy_dbcsr_to_fm(matrix_ks(ispin)%matrix, &
1555 scf_env%scf_work1(ispin))
1556 END DO
1557
1558 IF (scf_env%mixing_method == 1) THEN
1559 scf_env%iter_param = scf_env%p_mix_alpha
1560 scf_env%iter_method = "P_Mix/Lanczos"
1561 ELSE
1562! scf_env%iter_param = scf_env%mixing_store%alpha
1563 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Lanc."
1564 END IF
1565
1566 DO ispin = 1, SIZE(matrix_ks)
1567
1568 ks => scf_env%scf_work1(ispin)
1569 CALL cp_fm_uplo_to_full(ks, work)
1570
1571 CALL get_mo_set(mo_set=mos(ispin), &
1572 nao=nao, &
1573 nmo=nmo, &
1574 homo=homo, &
1575 eigenvalues=mo_eigenvalues, &
1576 mo_coeff=mo_coeff)
1577
1578 NULLIFY (c0, c1)
1579 c0 => scf_env%krylov_space%mo_conv(ispin)
1580 c1 => scf_env%krylov_space%mo_refine(ispin)
1581 SELECT CASE (scf_env%cholesky_method)
1582 CASE (cholesky_reduce)
1583 CALL cp_fm_cholesky_reduce(ks, ortho)
1584 CALL cp_fm_uplo_to_full(ks, work)
1585 CALL cp_fm_cholesky_restore(mo_coeff, nmo, ortho, c0, "MULTIPLY")
1586 CASE (cholesky_restore)
1587 CALL cp_fm_cholesky_restore(ks, nao, ortho, work, &
1588 "SOLVE", pos="RIGHT")
1589 CALL cp_fm_cholesky_restore(work, nao, ortho, ks, &
1590 "SOLVE", pos="LEFT", transa="T")
1591 CALL cp_fm_cholesky_restore(mo_coeff, nmo, ortho, c0, "MULTIPLY")
1592 CASE (cholesky_inverse)
1593 CALL cp_fm_cholesky_restore(ks, nao, ortho, work, &
1594 "MULTIPLY", pos="RIGHT")
1595 CALL cp_fm_cholesky_restore(work, nao, ortho, ks, &
1596 "MULTIPLY", pos="LEFT", transa="T")
1597 CALL cp_fm_cholesky_restore(mo_coeff, nmo, ortho, c0, "SOLVE")
1598 END SELECT
1599
1600 scf_env%krylov_space%nmo_nc = nmo
1601 scf_env%krylov_space%nmo_conv = 0
1602
1603 t1 = m_walltime()
1604 IF (output_unit > 0) THEN
1605 WRITE (output_unit, "(/T15,A)") '<<<<<<<<< LANCZOS REFINEMENT <<<<<<<<<<'
1606 WRITE (output_unit, "(T8,A,T15,A,T23,A,T36,A,T49,A,T60,A,/,T8,A)") &
1607 " Spin ", " Cycle ", &
1608 " conv. MOS ", " B2MAX ", " B2MIN ", " Time", repeat("-", 60)
1609 END IF
1610 eps_iter = max(scf_env%krylov_space%eps_conv, scf_env%krylov_space%eps_adapt*scf_env%iter_delta)
1611 iter = 0
1612 converged = .false.
1613 !Check convergence of MOS
1614 IF (my_check_moconv_only) THEN
1615
1616 CALL lanczos_refinement(scf_env%krylov_space, ks, c0, c1, mo_eigenvalues, &
1617 nao, eps_iter, ispin, check_moconv_only=my_check_moconv_only)
1618 t2 = m_walltime()
1619 IF (output_unit > 0) &
1620 WRITE (output_unit, '(T8,I3,T16,I5,T24,I6,T33,E12.4,2x,E12.4,T60,F8.3)') &
1621 ispin, iter, scf_env%krylov_space%nmo_conv, &
1622 scf_env%krylov_space%max_res_norm, scf_env%krylov_space%min_res_norm, t2 - t1
1623
1624 cycle
1625 ELSE
1626 !Block Lanczos refinement
1627 DO iter = 1, scf_env%krylov_space%max_iter
1628 CALL lanczos_refinement_2v(scf_env%krylov_space, ks, c0, c1, mo_eigenvalues, &
1629 nao, eps_iter, ispin)
1630 t2 = m_walltime()
1631 IF (output_unit > 0) THEN
1632 WRITE (output_unit, '(T8,I3,T16,I5,T24,I6,T33,E12.4,2x,E12.4,T60,F8.3)') &
1633 ispin, iter, scf_env%krylov_space%nmo_conv, &
1634 scf_env%krylov_space%max_res_norm, scf_env%krylov_space%min_res_norm, t2 - t1
1635 END IF
1636 t1 = m_walltime()
1637 IF (scf_env%krylov_space%max_res_norm < eps_iter) THEN
1638 converged = .true.
1639 IF (output_unit > 0) WRITE (output_unit, *) &
1640 " Reached convergence in ", iter, " iterations "
1641 EXIT
1642 END IF
1643 END DO
1644
1645 IF (.NOT. converged .AND. output_unit > 0) THEN
1646 WRITE (output_unit, "(T4, A)") " WARNING Lanczos refinement could "// &
1647 "not converge all the mos:"
1648 WRITE (output_unit, "(T40,A,T70,I10)") " number of not converged mos ", &
1649 scf_env%krylov_space%nmo_nc
1650 WRITE (output_unit, "(T40,A,T70,E10.2)") " max norm of the residual ", &
1651 scf_env%krylov_space%max_res_norm
1652
1653 END IF
1654
1655 ! For the moment skip the re-orthogonalization
1656 IF (.false.) THEN
1657 !Re-orthogonalization
1658 NULLIFY (chc, evec)
1659 chc => scf_env%krylov_space%chc_mat(ispin)
1660 evec => scf_env%krylov_space%c_vec(ispin)
1661 CALL parallel_gemm('N', 'N', nao, nmo, nao, rone, ks, c0, rzero, work)
1662 CALL parallel_gemm('T', 'N', nmo, nmo, nao, rone, c0, work, rzero, chc)
1663 !Diagonalize (C^t)HC
1664 CALL choose_eigv_solver(chc, evec, mo_eigenvalues)
1665 !Rotate the C vectors
1666 CALL parallel_gemm('N', 'N', nao, nmo, nmo, rone, c0, evec, rzero, c1)
1667 c0 => scf_env%krylov_space%mo_refine(ispin)
1668 END IF
1669
1670 IF (scf_env%cholesky_method == cholesky_inverse) THEN
1671 CALL cp_fm_cholesky_restore(c0, nmo, ortho, mo_coeff, "MULTIPLY")
1672 ELSE
1673 CALL cp_fm_cholesky_restore(c0, nmo, ortho, mo_coeff, "SOLVE")
1674 END IF
1675
1676 CALL set_mo_occupation(mo_set=mos(ispin), &
1677 smear=scf_control%smear)
1678
1679 ! does not yet handle k-points
1680 CALL calculate_density_matrix(mos(ispin), &
1681 scf_env%p_mix_new(ispin, 1)%matrix)
1682 END IF
1683 END DO ! ispin
1684
1685 IF (output_unit > 0) THEN
1686 WRITE (output_unit, "(T15,A/)") '<<<<<<<<< END LANCZOS REFINEMENT <<<<<<<<<<'
1687 END IF
1688
1689 CALL cp_print_key_finished_output(output_unit, logger, scf_section, &
1690 "PRINT%LANCZOS")
1691
1692 CALL timestop(handle)
1693
1694 END SUBROUTINE do_block_krylov_diag
1695
1696! **************************************************************************************************
1697!> \brief iterative diagonalization using the block davidson space approach
1698!> \param qs_env ...
1699!> \param scf_env ...
1700!> \param mos ...
1701!> \param matrix_ks ...
1702!> \param matrix_s ...
1703!> \param scf_control ...
1704!> \param scf_section ...
1705!> \param check_moconv_only ...
1706!> \param
1707!> \par History
1708!> 05.2011 created [MI]
1709! **************************************************************************************************
1710
1711 SUBROUTINE do_block_davidson_diag(qs_env, scf_env, mos, matrix_ks, matrix_s, &
1712 scf_control, scf_section, check_moconv_only)
1713
1714 TYPE(qs_environment_type), POINTER :: qs_env
1715 TYPE(qs_scf_env_type), POINTER :: scf_env
1716 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
1717 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: matrix_ks, matrix_s
1718 TYPE(scf_control_type), POINTER :: scf_control
1719 TYPE(section_vals_type), POINTER :: scf_section
1720 LOGICAL, INTENT(IN), OPTIONAL :: check_moconv_only
1721
1722 CHARACTER(LEN=*), PARAMETER :: routinen = 'do_block_davidson_diag'
1723
1724 INTEGER :: handle, ispin, nspins, output_unit
1725 LOGICAL :: do_prec, my_check_moconv_only
1726 TYPE(cp_logger_type), POINTER :: logger
1727
1728 logger => cp_get_default_logger()
1729 CALL timeset(routinen, handle)
1730
1731 output_unit = cp_print_key_unit_nr(logger, scf_section, "PRINT%DAVIDSON", &
1732 extension=".scfLog")
1733
1734 IF (output_unit > 0) &
1735 WRITE (output_unit, "(/T15,A)") '<<<<<<<<< DAVIDSON ITERATIONS <<<<<<<<<<'
1736
1737 IF (scf_env%mixing_method == 1) THEN
1738 scf_env%iter_param = scf_env%p_mix_alpha
1739 scf_env%iter_method = "P_Mix/Dav."
1740 ELSE
1741 scf_env%iter_param = scf_env%mixing_store%alpha
1742 scf_env%iter_method = trim(scf_env%mixing_store%iter_method)//"/Dav."
1743 END IF
1744
1745 my_check_moconv_only = .false.
1746 IF (PRESENT(check_moconv_only)) my_check_moconv_only = check_moconv_only
1747 do_prec = .false.
1748 IF (scf_env%block_davidson_env(1)%prec_type /= 0 .AND. &
1749 scf_env%iter_count >= scf_env%block_davidson_env(1)%first_prec) THEN
1750 do_prec = .true.
1751 END IF
1752
1753 nspins = SIZE(matrix_ks)
1754
1755 IF (do_prec .AND. (scf_env%iter_count == scf_env%block_davidson_env(1)%first_prec .OR. &
1756 modulo(scf_env%iter_count, scf_env%block_davidson_env(1)%niter_new_prec) == 0)) THEN
1757 CALL restart_preconditioner(qs_env, scf_env%ot_preconditioner, &
1758 prec_type=scf_env%block_davidson_env(1)%prec_type, nspins=nspins)
1759 CALL prepare_preconditioner(qs_env, mos, matrix_ks, matrix_s, scf_env%ot_preconditioner, &
1760 scf_env%block_davidson_env(1)%prec_type, &
1761 scf_env%block_davidson_env(1)%solver_type, &
1762 scf_env%block_davidson_env(1)%energy_gap, nspins, &
1763 convert_to_dbcsr=scf_env%block_davidson_env(1)%use_sparse_mos, &
1764 full_mo_set=.true.)
1765 END IF
1766
1767 DO ispin = 1, nspins
1768 IF (scf_env%block_davidson_env(ispin)%use_sparse_mos) THEN
1769 IF (.NOT. do_prec) THEN
1770 CALL generate_extended_space_sparse(scf_env%block_davidson_env(ispin), mos(ispin), &
1771 matrix_ks(ispin)%matrix, matrix_s(1)%matrix, output_unit)
1772 ELSE
1773 CALL generate_extended_space_sparse(scf_env%block_davidson_env(ispin), mos(ispin), &
1774 matrix_ks(ispin)%matrix, matrix_s(1)%matrix, output_unit, &
1775 scf_env%ot_preconditioner(ispin)%preconditioner)
1776 END IF
1777
1778 ELSE
1779 IF (.NOT. do_prec) THEN
1780 CALL generate_extended_space(scf_env%block_davidson_env(ispin), mos(ispin), &
1781 matrix_ks(ispin)%matrix, matrix_s(1)%matrix, output_unit)
1782 ELSE
1783 CALL generate_extended_space(scf_env%block_davidson_env(ispin), mos(ispin), &
1784 matrix_ks(ispin)%matrix, matrix_s(1)%matrix, output_unit, &
1785 scf_env%ot_preconditioner(ispin)%preconditioner)
1786 END IF
1787 END IF
1788 END DO !ispin
1789
1790 CALL set_mo_occupation(mo_array=mos, &
1791 smear=scf_control%smear)
1792
1793 DO ispin = 1, nspins
1794 ! does not yet handle k-points
1795 CALL calculate_density_matrix(mos(ispin), &
1796 scf_env%p_mix_new(ispin, 1)%matrix)
1797 END DO
1798
1799 IF (output_unit > 0) THEN
1800 WRITE (output_unit, "(T15,A/)") '<<<<<<<<< END DAVIDSON ITERATION <<<<<<<<<<'
1801 END IF
1802
1803 CALL cp_print_key_finished_output(output_unit, logger, scf_section, &
1804 "PRINT%DAVIDSON")
1805
1806 CALL timestop(handle)
1807
1808 END SUBROUTINE do_block_davidson_diag
1809
1810END MODULE qs_scf_diagonalization
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
various utilities that regard array of different kinds: output, allocation,... maybe it is not a good...
Basic linear algebra operations for complex full matrices.
subroutine, public cp_cfm_scale_and_add(alpha, matrix_a, beta, matrix_b)
Scale and add two BLACS matrices (a = alpha*a + beta*b).
subroutine, public cp_cfm_scale_and_add_fm(alpha, matrix_a, beta, matrix_b)
Scale and add two BLACS matrices (a = alpha*a + beta*b). where b is a real matrix (adapted from cp_cf...
used for collecting diagonalization schemes available for cp_cfm_type
Definition cp_cfm_diag.F:14
subroutine, public cp_cfm_geeig(amatrix, bmatrix, eigenvectors, eigenvalues, work)
General Eigenvalue Problem AX = BXE Single option version: Cholesky decomposition of B.
subroutine, public cp_cfm_geeig_canon(amatrix, bmatrix, eigenvectors, eigenvalues, work, epseig)
General Eigenvalue Problem AX = BXE Use canonical orthogonalization.
Represents a complex full matrix distributed on many processors.
subroutine, public cp_cfm_create(matrix, matrix_struct, name)
Creates a new full matrix with the given structure.
subroutine, public cp_cfm_release(matrix)
Releases a full matrix.
subroutine, public cp_cfm_to_fm(msource, mtargetr, mtargeti)
Copy real and imaginary parts of a complex full matrix into separate real-value full matrices.
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_deallocate_matrix(matrix)
...
subroutine, public dbcsr_desymmetrize(matrix_a, matrix_b)
...
subroutine, public dbcsr_copy(matrix_b, matrix_a, name, keep_sparsity, keep_imaginary)
...
subroutine, public dbcsr_set(matrix, alpha)
...
DBCSR operations in CP2K.
subroutine, public cp_dbcsr_sm_fm_multiply(matrix, fm_in, fm_out, ncol, alpha, beta)
multiply a dbcsr with a fm matrix
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
basic linear algebra operations for full matrices
subroutine, public cp_fm_cholesky_restore(fm_matrix, neig, fm_matrixb, fm_matrixout, op, pos, transa)
...
subroutine, public cp_fm_uplo_to_full(matrix, work, uplo)
given a triangular matrix according to uplo, computes the corresponding full matrix
subroutine, public cp_fm_symm(side, uplo, m, n, alpha, matrix_a, matrix_b, beta, matrix_c)
computes matrix_c = beta * matrix_c + alpha * matrix_a * matrix_b computes matrix_c = beta * matrix_c...
various cholesky decomposition related routines
subroutine, public cp_fm_cholesky_reduce(matrix, matrixb, itype)
reduce a matrix pencil A,B to normal form B has to be cholesky decomposed with cp_fm_cholesky_decompo...
used for collecting some of the diagonalization schemes available for cp_fm_type. cp_fm_power also mo...
Definition cp_fm_diag.F:17
subroutine, public cp_fm_geeig(amatrix, bmatrix, eigenvectors, eigenvalues, work)
General Eigenvalue Problem AX = BXE Single option version: Cholesky decomposition of B.
subroutine, public choose_eigv_solver(matrix, eigenvectors, eigenvalues, info)
Choose the Eigensolver depending on which library is available ELPA seems to be unstable for small sy...
Definition cp_fm_diag.F:229
subroutine, public cp_fm_geeig_canon(amatrix, bmatrix, eigenvectors, eigenvalues, work, epseig)
General Eigenvalue Problem AX = BXE Use canonical diagonalization : U*s**(-1/2)
pool for for elements that are retained and released
subroutine, public fm_pool_create_fm(pool, element, name)
returns an element, allocating it if none is in the pool
subroutine, public fm_pool_give_back_fm(pool, element)
returns the element to the pool
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_start_copy_general(source, destination, para_env, info)
Initiates the copy operation: get distribution data, post MPI isend and irecvs.
subroutine, public cp_fm_cleanup_copy_general(info)
Completes the copy operation: wait for comms clean up MPI state.
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
subroutine, public cp_fm_add_to_element(matrix, irow_global, icol_global, alpha)
...
subroutine, public cp_fm_finish_copy_general(destination, info)
Completes the copy operation: wait for comms, unpack, clean up MPI state.
subroutine, public cp_fm_create(matrix, matrix_struct, name, use_sp)
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 ...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
integer function, public cp_print_key_unit_nr(logger, basis_section, print_key_path, extension, middle_name, local, log_filename, ignore_should_output, file_form, file_position, file_action, file_status, do_backup, on_file, is_new_file, mpi_io, fout)
...
subroutine, public cp_print_key_finished_output(unit_nr, logger, basis_section, print_key_path, local, ignore_should_output, on_file, mpi_io)
should be called after you finish working with a unit obtained with cp_print_key_unit_nr,...
collects all constants needed in input so that they can be used without circular dependencies
integer, parameter, public core_guess
integer, parameter, public cholesky_restore
integer, parameter, public cholesky_dbcsr
integer, parameter, public cholesky_off
integer, parameter, public cholesky_reduce
integer, parameter, public high_spin_roks
integer, parameter, public cholesky_inverse
integer, parameter, public general_roks
integer, parameter, public restart_guess
objects that represent the structure of input sections and the data contained in an input section
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Defines the basic variable types.
Definition kinds.F:23
integer, parameter, public dp
Definition kinds.F:34
Routines needed for kpoint calculation.
subroutine, public rskp_transform(rmatrix, cmatrix, rsmat, ispin, xkp, cell_to_index, sab_nl, is_complex, rs_sign)
Transformation of real space matrices to a kpoint.
subroutine, public kpoint_density_transform(kpoint, denmat, wtype, tempmat, sab_nl, fmwork, for_aux_fit, pmat_ext)
generate real space density matrices in DBCSR format
subroutine, public kpoint_density_matrices(kpoint, energy_weighted, for_aux_fit)
Calculate kpoint density matrices (rho(k), owned by kpoint groups)
subroutine, public kpoint_set_mo_occupation(kpoint, smear)
Given the eigenvalues of all kpoints, calculates the occupation numbers.
Types and basic routines needed for a kpoint calculation.
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym)
Retrieve information from a kpoint environment.
Machine interface based on Fortran 2003 and POSIX.
Definition machine.F:17
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Definition machine.F:130
real(kind=dp) function, public m_walltime()
returns time from a real-time clock, protected against rolling early/easily
Definition machine.F:147
Interface to the message passing library MPI.
basic linear algebra operations for full matrixes
computes preconditioners, and implements methods to apply them currently used in qs_ot
subroutine, public restart_preconditioner(qs_env, preconditioner, prec_type, nspins)
Allows for a restart of the preconditioner depending on the method it purges all arrays or keeps them...
subroutine, public prepare_preconditioner(qs_env, mos, matrix_ks, matrix_s, ot_preconditioner, prec_type, solver_type, energy_gap, nspins, has_unit_metric, convert_to_dbcsr, chol_type, full_mo_set)
...
collects routines that calculate density matrices
module that contains the definitions of the scf types
integer, parameter, public direct_mixing_nr
integer, parameter, public gspace_mixing_nr
Apply the direct inversion in the iterative subspace (DIIS) of Pulay in the framework of an SCF itera...
Definition qs_diis.F:21
subroutine, public qs_diis_b_info_kp(diis_buffer, ib, nb)
Update info about the current buffer step ib and the current number of buffers nb.
Definition qs_diis.F:989
subroutine, public qs_diis_b_calc_err_kp(diis_buffer, ib, mos, kc, sc, ispin, ikp, nkp_local, scf_section)
Calculate and store the error for a given k-point.
Definition qs_diis.F:1013
subroutine, public qs_diis_b_step_kp(diis_buffer, coeffs, ib, nb, delta, error_max, diis_step, eps_diis, nspin, nkp, nkp_local, nmixing, scf_section, para_env)
Update the SCF DIIS buffer, and if appropriate does a diis step, for k-points.
Definition qs_diis.F:1100
subroutine, public qs_diis_b_step(diis_buffer, mo_array, kc, sc, delta, error_max, diis_step, eps_diis, nmixing, s_matrix, scf_section, roks)
Update the SCF DIIS buffer, and if appropriate does a diis step.
Definition qs_diis.F:230
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_pp, 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, 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)
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 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...
wrapper for the pools of matrixes
subroutine, public mpools_get(mpools, ao_mo_fm_pools, ao_ao_fm_pools, mo_mo_fm_pools, ao_mosub_fm_pools, mosub_mosub_fm_pools, maxao_maxmo_fm_pool, maxao_maxao_fm_pool, maxmo_maxmo_fm_pool)
returns various attributes of the mpools (notably the pools contained in it)
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
subroutine, public self_consistency_check(rho_ao, p_delta, para_env, p_out, delta)
...
collects routines that perform operations directly related to MOs
Set occupation of molecular orbitals.
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Define the neighbor list data types and the corresponding functionality.
an eigen-space solver for the generalised symmetric eigenvalue problem for sparse matrices,...
subroutine, public ot_eigensolver(matrix_h, matrix_s, matrix_orthogonal_space_fm, matrix_c_fm, preconditioner, eps_gradient, iter_max, size_ortho_space, silent, ot_settings)
...
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 algorithms to perform an itrative diagonalization by the block-Davidson appr...
subroutine, public generate_extended_space_sparse(bdav_env, mo_set, matrix_h, matrix_s, output_unit, preconditioner)
...
subroutine, public generate_extended_space(bdav_env, mo_set, matrix_h, matrix_s, output_unit, preconditioner)
...
Different diagonalization schemes that can be used for the iterative solution of the eigenvalue probl...
subroutine, public general_eigenproblem(scf_env, mos, matrix_ks, matrix_s, scf_control, scf_section, diis_step)
the inner loop of scf, specific to diagonalization with S matrix basically, in goes the ks matrix out...
subroutine, public diag_subspace_allocate(subspace_env, qs_env, mos)
...
subroutine, public do_ot_diag(scf_env, mos, matrix_ks, matrix_s, scf_control, scf_section, diis_step)
the inner loop of scf, specific to iterative diagonalization using OT with S matrix; basically,...
subroutine, public do_block_davidson_diag(qs_env, scf_env, mos, matrix_ks, matrix_s, scf_control, scf_section, check_moconv_only)
iterative diagonalization using the block davidson space approach
subroutine, public do_roks_diag(scf_env, mos, matrix_ks, matrix_s, scf_control, scf_section, diis_step, orthogonal_basis)
Solve a set restricted open Kohn-Sham (ROKS) equations based on the alpha and beta Kohn-Sham matrices...
subroutine, public do_general_diag_kp(matrix_ks, matrix_s, kpoints, scf_env, scf_control, update_p, diis_step, diis_error, qs_env)
Kpoint diagonalization routine Transforms matrices to kpoint, distributes kpoint groups,...
subroutine, public do_scf_diag_subspace(qs_env, scf_env, subspace_env, mos, rho, ks_env, scf_section, scf_control)
inner loop within MOS subspace, to refine occupation and density, before next diagonalization of the ...
subroutine, public do_block_krylov_diag(scf_env, mos, matrix_ks, scf_control, scf_section, check_moconv_only)
iterative diagonalization using the block Krylov-space approach
subroutine, public do_special_diag(scf_env, mos, matrix_ks, scf_control, scf_section, diis_step)
the inner loop of scf, specific to diagonalization without S matrix basically, in goes the ks matrix ...
subroutine, public do_general_diag(scf_env, mos, matrix_ks, matrix_s, scf_control, scf_section, diis_step)
...
module that contains the algorithms to perform an itrative diagonalization by the block-Lanczos appro...
subroutine, public lanczos_refinement(krylov_space, ks, c0, c1, eval, nao, eps_iter, ispin, check_moconv_only)
lanczos refinement by blocks of not-converged MOs
subroutine, public lanczos_refinement_2v(krylov_space, ks, c0, c1, eval, nao, eps_iter, ispin, check_moconv_only)
...
groups fairly general SCF methods, so that modules other than qs_scf can use them too split off from ...
subroutine, public eigensolver_simple(matrix_ks, mo_set, work, do_level_shift, level_shift, use_jacobi, jacobi_threshold)
...
subroutine, public eigensolver_dbcsr(matrix_ks, matrix_ks_fm, mo_set, ortho_dbcsr, ksbuf1, ksbuf2)
...
subroutine, public scf_env_density_mixing(p_mix_new, mixing_store, rho_ao, para_env, iter_delta, iter_count, diis, invert)
perform (if requested) a density mixing
subroutine, public eigensolver(matrix_ks_fm, mo_set, ortho, work, cholesky_method, do_level_shift, level_shift, matrix_u_fm, use_jacobi)
Diagonalise the Kohn-Sham matrix to get a new set of MO eigen- vectors and MO eigenvalues....
subroutine, public eigensolver_symm(matrix_ks_fm, mo_set, ortho, work, do_level_shift, level_shift, matrix_u_fm, use_jacobi, jacobi_threshold, ortho_red, work_red, matrix_ks_fm_red, matrix_u_fm_red)
...
module that contains the definitions of the scf types
parameters that control an scf iteration
represent a pointer to a 1d array
Represent a complex full matrix.
keeps the information about the structure of a full matrix
Stores the state of a copy between cp_fm_start_copy_general and cp_fm_finish_copy_general.
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...
Keeps information about a specific k-point.
Contains information about kpoints.
stores all the informations relevant to an mpi environment
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
container for the pools of matrixes used by qs
keeps the density in various representations, keeping track of which ones are valid.