34 fm_pools_create_fm_vect,&
35 fm_pools_give_back_fm_vect
50 USE dbcsr_api,
ONLY: dbcsr_add,&
52 dbcsr_deallocate_matrix,&
91 #include "./base/base_uses.f90"
96 LOGICAL,
PRIVATE,
PARAMETER :: debug_this_module = .true.
97 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_wf_history_methods'
113 SUBROUTINE wfs_create(snapshot)
114 TYPE(qs_wf_snapshot_type),
INTENT(OUT) :: snapshot
116 NULLIFY (snapshot%wf, snapshot%rho_r, &
117 snapshot%rho_g, snapshot%rho_ao, snapshot%rho_ao_kp, &
118 snapshot%overlap, snapshot%rho_frozen)
120 END SUBROUTINE wfs_create
133 SUBROUTINE wfs_update(snapshot, wf_history, qs_env, dt)
134 TYPE(qs_wf_snapshot_type),
POINTER :: snapshot
135 TYPE(qs_wf_history_type),
POINTER :: wf_history
136 TYPE(qs_environment_type),
POINTER :: qs_env
137 REAL(KIND=
dp),
INTENT(in),
OPTIONAL :: dt
139 CHARACTER(len=*),
PARAMETER :: routineN =
'wfs_update'
141 INTEGER :: handle, img, ispin, nimg, nspins
142 TYPE(cp_fm_pool_p_type),
DIMENSION(:),
POINTER :: ao_mo_pools
143 TYPE(cp_fm_type),
POINTER :: mo_coeff
144 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s, rho_ao
145 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: rho_ao_kp
146 TYPE(dft_control_type),
POINTER :: dft_control
147 TYPE(mo_set_type),
DIMENSION(:),
POINTER :: mos
148 TYPE(pw_c1d_gs_type),
DIMENSION(:),
POINTER :: rho_g
149 TYPE(pw_env_type),
POINTER :: pw_env
150 TYPE(pw_pool_type),
POINTER :: auxbas_pw_pool
151 TYPE(pw_r3d_rs_type),
DIMENSION(:),
POINTER :: rho_r
152 TYPE(qs_rho_type),
POINTER :: rho
154 CALL timeset(routinen, handle)
156 NULLIFY (pw_env, auxbas_pw_pool, ao_mo_pools, dft_control, mos, mo_coeff, &
157 rho, rho_r, rho_g, rho_ao, matrix_s)
159 dft_control=dft_control, rho=rho)
160 CALL mpools_get(qs_env%mpools, ao_mo_fm_pools=ao_mo_pools)
161 CALL pw_env_get(pw_env, auxbas_pw_pool=auxbas_pw_pool)
163 cpassert(
ASSOCIATED(wf_history))
164 cpassert(
ASSOCIATED(dft_control))
165 IF (.NOT.
ASSOCIATED(snapshot))
THEN
167 CALL wfs_create(snapshot)
169 cpassert(wf_history%ref_count > 0)
171 nspins = dft_control%nspins
173 IF (
PRESENT(dt)) snapshot%dt = dt
174 IF (wf_history%store_wf)
THEN
176 IF (.NOT.
ASSOCIATED(snapshot%wf))
THEN
177 CALL fm_pools_create_fm_vect(ao_mo_pools, snapshot%wf, &
179 cpassert(nspins ==
SIZE(snapshot%wf))
182 CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff)
183 CALL cp_fm_to_fm(mo_coeff, snapshot%wf(ispin))
186 CALL fm_pools_give_back_fm_vect(ao_mo_pools, snapshot%wf)
189 IF (wf_history%store_rho_r)
THEN
191 cpassert(
ASSOCIATED(rho_r))
192 IF (.NOT.
ASSOCIATED(snapshot%rho_r))
THEN
193 ALLOCATE (snapshot%rho_r(nspins))
195 CALL auxbas_pw_pool%create_pw(snapshot%rho_r(ispin))
199 CALL pw_copy(rho_r(ispin), snapshot%rho_r(ispin))
201 ELSE IF (
ASSOCIATED(snapshot%rho_r))
THEN
202 DO ispin = 1,
SIZE(snapshot%rho_r)
203 CALL auxbas_pw_pool%give_back_pw(snapshot%rho_r(ispin))
205 DEALLOCATE (snapshot%rho_r)
208 IF (wf_history%store_rho_g)
THEN
210 cpassert(
ASSOCIATED(rho_g))
211 IF (.NOT.
ASSOCIATED(snapshot%rho_g))
THEN
212 ALLOCATE (snapshot%rho_g(nspins))
214 CALL auxbas_pw_pool%create_pw(snapshot%rho_g(ispin))
218 CALL pw_copy(rho_g(ispin), snapshot%rho_g(ispin))
220 ELSE IF (
ASSOCIATED(snapshot%rho_g))
THEN
221 DO ispin = 1,
SIZE(snapshot%rho_g)
222 CALL auxbas_pw_pool%give_back_pw(snapshot%rho_g(ispin))
224 DEALLOCATE (snapshot%rho_g)
227 IF (
ASSOCIATED(snapshot%rho_ao))
THEN
231 IF (wf_history%store_rho_ao)
THEN
233 cpassert(
ASSOCIATED(rho_ao))
237 ALLOCATE (snapshot%rho_ao(ispin)%matrix)
238 CALL dbcsr_copy(snapshot%rho_ao(ispin)%matrix, rho_ao(ispin)%matrix)
242 IF (
ASSOCIATED(snapshot%rho_ao_kp))
THEN
246 IF (wf_history%store_rho_ao_kp)
THEN
248 cpassert(
ASSOCIATED(rho_ao_kp))
250 nimg = dft_control%nimages
254 ALLOCATE (snapshot%rho_ao_kp(ispin, img)%matrix)
255 CALL dbcsr_copy(snapshot%rho_ao_kp(ispin, img)%matrix, &
256 rho_ao_kp(ispin, img)%matrix)
261 IF (
ASSOCIATED(snapshot%overlap))
THEN
263 CALL dbcsr_deallocate_matrix(snapshot%overlap)
265 IF (wf_history%store_overlap)
THEN
267 cpassert(
ASSOCIATED(matrix_s))
268 cpassert(
ASSOCIATED(matrix_s(1)%matrix))
269 ALLOCATE (snapshot%overlap)
270 CALL dbcsr_copy(snapshot%overlap, matrix_s(1)%matrix)
273 IF (wf_history%store_frozen_density)
THEN
278 CALL timestop(handle)
280 END SUBROUTINE wfs_update
294 SUBROUTINE wfi_create(wf_history, interpolation_method_nr, extrapolation_order, &
296 TYPE(qs_wf_history_type),
POINTER :: wf_history
297 INTEGER,
INTENT(in) :: interpolation_method_nr, &
299 LOGICAL,
INTENT(IN) :: has_unit_metric
301 CHARACTER(len=*),
PARAMETER :: routinen =
'wfi_create'
305 CALL timeset(routinen, handle)
307 ALLOCATE (wf_history)
308 wf_history%ref_count = 1
309 wf_history%memory_depth = 0
310 wf_history%snapshot_count = 0
311 wf_history%last_state_index = 1
312 wf_history%store_wf = .false.
313 wf_history%store_rho_r = .false.
314 wf_history%store_rho_g = .false.
315 wf_history%store_rho_ao = .false.
316 wf_history%store_rho_ao_kp = .false.
317 wf_history%store_overlap = .false.
318 wf_history%store_frozen_density = .false.
319 NULLIFY (wf_history%past_states)
321 wf_history%interpolation_method_nr = interpolation_method_nr
323 SELECT CASE (wf_history%interpolation_method_nr)
325 wf_history%memory_depth = 0
327 wf_history%memory_depth = 0
329 wf_history%memory_depth = 1
330 wf_history%store_rho_ao = .true.
332 wf_history%memory_depth = 1
333 wf_history%store_rho_ao = .true.
335 wf_history%memory_depth = 2
336 wf_history%store_wf = .true.
338 wf_history%memory_depth = 2
339 wf_history%store_rho_ao = .true.
341 wf_history%memory_depth = 2
342 wf_history%store_wf = .true.
343 IF (.NOT. has_unit_metric) wf_history%store_overlap = .true.
346 wf_history%memory_depth = extrapolation_order + 1
347 wf_history%store_wf = .true.
348 IF (.NOT. has_unit_metric) wf_history%store_overlap = .true.
350 wf_history%memory_depth = 1
351 wf_history%store_frozen_density = .true.
353 wf_history%memory_depth = extrapolation_order + 2
354 wf_history%store_wf = .true.
355 IF (.NOT. has_unit_metric) wf_history%store_overlap = .true.
357 CALL cp_abort(__location__, &
358 "Unknown interpolation method: "// &
359 trim(adjustl(cp_to_string(interpolation_method_nr))))
362 ALLOCATE (wf_history%past_states(wf_history%memory_depth))
364 DO i = 1,
SIZE(wf_history%past_states)
365 NULLIFY (wf_history%past_states(i)%snapshot)
368 CALL timestop(handle)
379 TYPE(qs_wf_history_type),
POINTER :: wf_history
381 cpassert(
ASSOCIATED(wf_history))
382 IF (wf_history%store_rho_ao)
THEN
383 wf_history%store_rho_ao_kp = .true.
384 wf_history%store_rho_ao = .false.
387 IF (wf_history%store_wf)
THEN
388 cpabort(
"WFN based interpolation method not possible for kpoints.")
390 IF (wf_history%store_frozen_density)
THEN
391 cpabort(
"Frozen density initialization method not possible for kpoints.")
393 IF (wf_history%store_overlap)
THEN
394 cpabort(
"Inconsistent interpolation method for kpoints.")
408 INTEGER,
INTENT(in) :: method_nr
409 CHARACTER(len=30) :: res
412 SELECT CASE (method_nr)
418 res =
"previous_rho_r"
420 res =
"initial_guess"
430 res =
"frozen density approximation"
434 CALL cp_abort(__location__, &
435 "Unknown interpolation method: "// &
436 trim(adjustl(cp_to_string(method_nr))))
456 TYPE(qs_wf_history_type),
POINTER :: wf_history
457 TYPE(qs_environment_type),
POINTER :: qs_env
458 REAL(kind=
dp),
INTENT(IN) :: dt
459 INTEGER,
INTENT(OUT),
OPTIONAL :: extrapolation_method_nr
460 LOGICAL,
INTENT(OUT),
OPTIONAL :: orthogonal_wf
462 CHARACTER(len=*),
PARAMETER :: routinen =
'wfi_extrapolate'
464 INTEGER :: actual_extrapolation_method_nr, handle, &
465 i, img, io_unit, ispin, k, n, nmo, &
467 LOGICAL :: do_kpoints, my_orthogonal_wf, use_overlap
468 REAL(kind=
dp) :: alpha, t0, t1, t2
469 TYPE(cp_fm_pool_p_type),
DIMENSION(:),
POINTER :: ao_mo_fm_pools
470 TYPE(cp_fm_struct_type),
POINTER :: matrix_struct, matrix_struct_new
471 TYPE(cp_fm_type) :: csc, fm_tmp
472 TYPE(cp_fm_type),
POINTER :: mo_coeff
473 TYPE(cp_logger_type),
POINTER :: logger
474 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: rho_ao, rho_frozen_ao
475 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: rho_ao_kp
476 TYPE(mo_set_type),
DIMENSION(:),
POINTER :: mos
477 TYPE(qs_rho_type),
POINTER :: rho
478 TYPE(qs_wf_snapshot_type),
POINTER :: t0_state, t1_state
480 NULLIFY (mos, ao_mo_fm_pools, t0_state, t1_state, mo_coeff, &
481 rho, rho_ao, rho_frozen_ao)
483 use_overlap = wf_history%store_overlap
485 CALL timeset(routinen, handle)
487 print_level = logger%iter_info%print_level
491 cpassert(
ASSOCIATED(wf_history))
492 cpassert(wf_history%ref_count > 0)
493 cpassert(
ASSOCIATED(qs_env))
494 CALL get_qs_env(qs_env, mos=mos, rho=rho, do_kpoints=do_kpoints)
495 CALL mpools_get(qs_env%mpools, ao_mo_fm_pools=ao_mo_fm_pools)
497 IF (wf_history%snapshot_count < 1)
THEN
500 actual_extrapolation_method_nr = wf_history%interpolation_method_nr
503 SELECT CASE (actual_extrapolation_method_nr)
505 IF (wf_history%snapshot_count < 2)
THEN
509 IF (wf_history%snapshot_count < 2)
THEN
517 IF (wf_history%snapshot_count < 2)
THEN
522 IF (
PRESENT(extrapolation_method_nr)) &
523 extrapolation_method_nr = actual_extrapolation_method_nr
524 my_orthogonal_wf = .false.
526 SELECT CASE (actual_extrapolation_method_nr)
528 cpassert(.NOT. do_kpoints)
530 cpassert(
ASSOCIATED(t0_state%rho_frozen))
532 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
533 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
535 CALL qs_rho_get(t0_state%rho_frozen, rho_ao=rho_frozen_ao)
537 DO ispin = 1,
SIZE(rho_frozen_ao)
538 CALL dbcsr_copy(rho_ao(ispin)%matrix, &
539 rho_frozen_ao(ispin)%matrix, &
540 keep_sparsity=.true.)
547 my_orthogonal_wf = .false.
550 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
551 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
553 cpassert(
ASSOCIATED(t0_state%rho_ao_kp))
555 DO ispin = 1,
SIZE(t0_state%rho_ao_kp, 1)
556 DO img = 1,
SIZE(t0_state%rho_ao_kp, 2)
557 IF (img >
SIZE(rho_ao_kp, 2))
THEN
558 cpwarn(
"Change in cell neighborlist: might affect quality of initial guess")
560 CALL dbcsr_copy(rho_ao_kp(ispin, img)%matrix, &
561 t0_state%rho_ao_kp(ispin, img)%matrix, &
562 keep_sparsity=.true.)
567 cpassert(
ASSOCIATED(t0_state%rho_ao))
569 DO ispin = 1,
SIZE(t0_state%rho_ao)
570 CALL dbcsr_copy(rho_ao(ispin)%matrix, &
571 t0_state%rho_ao(ispin)%matrix, &
572 keep_sparsity=.true.)
580 cpassert(.NOT. do_kpoints)
581 my_orthogonal_wf = .true.
582 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
583 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
585 DO ispin = 1,
SIZE(mos)
586 CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff, &
591 CALL calculate_density_matrix(mo_set=mos(ispin), &
592 density_matrix=rho_ao(ispin)%matrix)
599 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
600 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
602 cpassert(
ASSOCIATED(t0_state%rho_ao_kp))
604 DO ispin = 1,
SIZE(t0_state%rho_ao_kp, 1)
605 DO img = 1,
SIZE(t0_state%rho_ao_kp, 2)
606 IF (img >
SIZE(rho_ao_kp, 2))
THEN
607 cpwarn(
"Change in cell neighborlist: might affect quality of initial guess")
609 CALL dbcsr_copy(rho_ao_kp(ispin, img)%matrix, &
610 t0_state%rho_ao_kp(ispin, img)%matrix, &
611 keep_sparsity=.true.)
616 cpassert(
ASSOCIATED(t0_state%rho_ao))
618 DO ispin = 1,
SIZE(t0_state%rho_ao)
619 CALL dbcsr_copy(rho_ao(ispin)%matrix, &
620 t0_state%rho_ao(ispin)%matrix, &
621 keep_sparsity=.true.)
635 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
636 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
638 cpassert(.NOT. do_kpoints)
641 cpassert(
ASSOCIATED(t0_state))
642 cpassert(
ASSOCIATED(t1_state))
643 cpassert(
ASSOCIATED(t0_state%wf))
644 cpassert(
ASSOCIATED(t1_state%wf))
645 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
646 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
648 my_orthogonal_wf = .true.
653 DO ispin = 1,
SIZE(mos)
654 CALL get_mo_set(mos(ispin), mo_coeff=mo_coeff, &
658 matrix_b=t1_state%wf(ispin), &
659 beta=(t2 - t0)/(t1 - t0))
663 beta=(t1 - t2)/(t1 - t0), matrix_b=t0_state%wf(ispin))
667 CALL calculate_density_matrix(mo_set=mos(ispin), &
668 density_matrix=rho_ao(ispin)%matrix)
677 cpassert(
ASSOCIATED(t0_state))
678 cpassert(
ASSOCIATED(t1_state))
680 cpassert(
ASSOCIATED(t0_state%rho_ao_kp))
681 cpassert(
ASSOCIATED(t1_state%rho_ao_kp))
683 cpassert(
ASSOCIATED(t0_state%rho_ao))
684 cpassert(
ASSOCIATED(t1_state%rho_ao))
686 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
687 CALL wfi_set_history_variables(qs_env=qs_env, nvec=nvec)
694 DO ispin = 1,
SIZE(rho_ao_kp, 1)
695 DO img = 1,
SIZE(rho_ao_kp, 2)
696 IF (img >
SIZE(t0_state%rho_ao_kp, 2) .OR. &
697 img >
SIZE(t1_state%rho_ao_kp, 2))
THEN
698 cpwarn(
"Change in cell neighborlist: might affect quality of initial guess")
700 CALL dbcsr_add(rho_ao_kp(ispin, img)%matrix, t1_state%rho_ao_kp(ispin, img)%matrix, &
701 alpha_scalar=0.0_dp, beta_scalar=(t2 - t0)/(t1 - t0))
702 CALL dbcsr_add(rho_ao_kp(ispin, img)%matrix, t0_state%rho_ao_kp(ispin, img)%matrix, &
703 alpha_scalar=1.0_dp, beta_scalar=(t1 - t2)/(t1 - t0))
709 DO ispin = 1,
SIZE(rho_ao)
710 CALL dbcsr_add(rho_ao(ispin)%matrix, t1_state%rho_ao(ispin)%matrix, &
711 alpha_scalar=0.0_dp, beta_scalar=(t2 - t0)/(t1 - t0))
712 CALL dbcsr_add(rho_ao(ispin)%matrix, t0_state%rho_ao(ispin)%matrix, &
713 alpha_scalar=1.0_dp, beta_scalar=(t1 - t2)/(t1 - t0))
721 cpassert(.NOT. do_kpoints)
724 cpassert(
ASSOCIATED(t0_state))
725 cpassert(
ASSOCIATED(t1_state))
726 cpassert(
ASSOCIATED(t0_state%wf))
727 cpassert(
ASSOCIATED(t1_state%wf))
728 IF (wf_history%store_overlap)
THEN
729 cpassert(
ASSOCIATED(t0_state%overlap))
730 cpassert(
ASSOCIATED(t1_state%overlap))
732 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
733 IF (nvec >= wf_history%memory_depth)
THEN
734 IF ((qs_env%scf_control%max_scf_hist .NE. 0) .AND. (qs_env%scf_control%eps_scf_hist .NE. 0))
THEN
735 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
736 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
737 qs_env%scf_control%outer_scf%have_scf = .false.
738 ELSE IF (qs_env%scf_control%max_scf_hist .NE. 0)
THEN
739 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
740 qs_env%scf_control%outer_scf%have_scf = .false.
741 ELSE IF (qs_env%scf_control%eps_scf_hist .NE. 0)
THEN
742 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
746 my_orthogonal_wf = .true.
750 DO ispin = 1,
SIZE(mos)
751 NULLIFY (mo_coeff, matrix_struct, matrix_struct_new)
752 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
754 matrix_struct=matrix_struct)
756 nrow_global=k, ncol_global=k)
760 IF (use_overlap)
THEN
762 CALL parallel_gemm(
'T',
'N', k, k, n, 1.0_dp, t0_state%wf(ispin), mo_coeff, 0.0_dp, csc)
764 CALL parallel_gemm(
'T',
'N', k, k, n, 1.0_dp, t0_state%wf(ispin), &
765 t1_state%wf(ispin), 0.0_dp, csc)
767 CALL parallel_gemm(
'N',
'N', n, k, k, 1.0_dp, t0_state%wf(ispin), csc, 0.0_dp, mo_coeff)
768 CALL cp_fm_release(csc)
773 CALL calculate_density_matrix(mo_set=mos(ispin), &
774 density_matrix=rho_ao(ispin)%matrix)
780 cpassert(.NOT. do_kpoints)
782 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
783 cpassert(nvec .GT. 0)
784 IF (nvec >= wf_history%memory_depth)
THEN
785 IF ((qs_env%scf_control%max_scf_hist .NE. 0) .AND. (qs_env%scf_control%eps_scf_hist .NE. 0))
THEN
786 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
787 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
788 qs_env%scf_control%outer_scf%have_scf = .false.
789 ELSE IF (qs_env%scf_control%max_scf_hist .NE. 0)
THEN
790 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
791 qs_env%scf_control%outer_scf%have_scf = .false.
792 ELSE IF (qs_env%scf_control%eps_scf_hist .NE. 0)
THEN
793 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
797 my_orthogonal_wf = .true.
798 DO ispin = 1,
SIZE(mos)
799 NULLIFY (mo_coeff, matrix_struct, matrix_struct_new)
800 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
802 matrix_struct=matrix_struct)
805 nrow_global=k, ncol_global=k)
810 CALL cp_fm_to_fm(t1_state%wf(ispin), mo_coeff)
816 IF (use_overlap)
THEN
818 CALL parallel_gemm(
'T',
'N', k, k, n, 1.0_dp, t0_state%wf(ispin), fm_tmp, 0.0_dp, csc)
820 CALL parallel_gemm(
'T',
'N', k, k, n, 1.0_dp, t0_state%wf(ispin), &
821 t1_state%wf(ispin), 0.0_dp, csc)
823 CALL parallel_gemm(
'N',
'N', n, k, k, 1.0_dp, t0_state%wf(ispin), csc, 0.0_dp, fm_tmp)
824 alpha = -1.0_dp*alpha*real(nvec - i + 1,
dp)/real(i,
dp)
828 CALL cp_fm_release(csc)
829 CALL cp_fm_release(fm_tmp)
833 CALL calculate_density_matrix(mo_set=mos(ispin), &
834 density_matrix=rho_ao(ispin)%matrix)
840 cpassert(.NOT. do_kpoints)
844 nvec = min(wf_history%memory_depth, wf_history%snapshot_count)
845 cpassert(nvec .GT. 0)
846 IF (nvec >= wf_history%memory_depth)
THEN
847 IF ((qs_env%scf_control%max_scf_hist .NE. 0) .AND. &
848 (qs_env%scf_control%eps_scf_hist .NE. 0))
THEN
849 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
850 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
851 qs_env%scf_control%outer_scf%have_scf = .false.
852 ELSE IF (qs_env%scf_control%max_scf_hist .NE. 0)
THEN
853 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
854 qs_env%scf_control%outer_scf%have_scf = .false.
855 ELSE IF (qs_env%scf_control%eps_scf_hist .NE. 0)
THEN
856 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
860 my_orthogonal_wf = .true.
862 DO ispin = 1,
SIZE(mos)
863 NULLIFY (mo_coeff, matrix_struct, matrix_struct_new)
864 CALL get_mo_set(mo_set=mos(ispin), mo_coeff=mo_coeff)
868 matrix_struct=matrix_struct)
871 template_fmstruct=matrix_struct, &
879 CALL cp_fm_to_fm(t1_state%wf(ispin), mo_coeff)
880 alpha = real(4*nvec - 2, kind=
dp)/real(nvec + 1, kind=
dp)
882 WRITE (unit=io_unit, fmt=
"(/,T2,A,/,/,T3,A,I0,/,/,T3,A2,I0,A4,F10.6)") &
883 "Parameters for the always stable predictor-corrector (ASPC) method:", &
884 "ASPC order: ", max(nvec - 2, 0), &
885 "B(", 1,
") = ", alpha
891 IF (use_overlap)
THEN
893 CALL parallel_gemm(
'T',
'N', k, k, n, 1.0_dp, t0_state%wf(ispin), fm_tmp, 0.0_dp, csc)
895 CALL parallel_gemm(
'T',
'N', k, k, n, 1.0_dp, t0_state%wf(ispin), &
896 t1_state%wf(ispin), 0.0_dp, csc)
898 CALL parallel_gemm(
'N',
'N', n, k, k, 1.0_dp, t0_state%wf(ispin), csc, 0.0_dp, fm_tmp)
899 alpha = (-1.0_dp)**(i + 1)*real(i, kind=
dp)* &
902 WRITE (unit=io_unit, fmt=
"(T3,A2,I0,A4,F10.6)") &
903 "B(", i,
") = ", alpha
907 CALL cp_fm_release(csc)
908 CALL cp_fm_release(fm_tmp)
912 CALL calculate_density_matrix(mo_set=mos(ispin), &
913 density_matrix=rho_ao(ispin)%matrix)
919 CALL cp_abort(__location__, &
920 "Unknown interpolation method: "// &
921 trim(adjustl(cp_to_string(wf_history%interpolation_method_nr))))
923 IF (
PRESENT(orthogonal_wf)) orthogonal_wf = my_orthogonal_wf
925 "DFT%SCF%PRINT%PROGRAM_RUN_INFO")
926 CALL timestop(handle)
938 ELEMENTAL SUBROUTINE wfi_set_history_variables(qs_env, nvec)
939 TYPE(qs_environment_type),
INTENT(INOUT) :: qs_env
940 INTEGER,
INTENT(IN) :: nvec
942 IF (nvec >= qs_env%wf_history%memory_depth)
THEN
943 IF ((qs_env%scf_control%max_scf_hist .NE. 0) .AND. (qs_env%scf_control%eps_scf_hist .NE. 0))
THEN
944 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
945 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
946 qs_env%scf_control%outer_scf%have_scf = .false.
947 ELSE IF (qs_env%scf_control%max_scf_hist .NE. 0)
THEN
948 qs_env%scf_control%max_scf = qs_env%scf_control%max_scf_hist
949 qs_env%scf_control%outer_scf%have_scf = .false.
950 ELSE IF (qs_env%scf_control%eps_scf_hist .NE. 0)
THEN
951 qs_env%scf_control%eps_scf = qs_env%scf_control%eps_scf_hist
952 qs_env%scf_control%outer_scf%eps_scf = qs_env%scf_control%eps_scf_hist
956 END SUBROUTINE wfi_set_history_variables
968 TYPE(qs_wf_history_type),
POINTER :: wf_history
969 TYPE(qs_environment_type),
POINTER :: qs_env
970 REAL(kind=
dp),
INTENT(in) :: dt
972 cpassert(
ASSOCIATED(wf_history))
973 cpassert(wf_history%ref_count > 0)
974 cpassert(
ASSOCIATED(qs_env))
976 wf_history%snapshot_count = wf_history%snapshot_count + 1
977 IF (wf_history%memory_depth > 0)
THEN
978 wf_history%last_state_index =
modulo(wf_history%snapshot_count, &
979 wf_history%memory_depth) + 1
980 CALL wfs_update(snapshot=wf_history%past_states &
981 (wf_history%last_state_index)%snapshot, wf_history=wf_history, &
982 qs_env=qs_env, dt=dt)
996 TYPE(qs_environment_type),
POINTER :: qs_env
997 TYPE(cp_fm_type),
INTENT(IN) :: v_matrix
998 INTEGER,
INTENT(in),
OPTIONAL :: n_col
1000 CHARACTER(len=*),
PARAMETER :: routinen =
'reorthogonalize_vectors'
1002 INTEGER :: handle, my_n_col
1003 LOGICAL :: has_unit_metric, &
1004 ortho_contains_cholesky, &
1006 TYPE(cp_fm_pool_type),
POINTER :: maxao_maxmo_fm_pool
1007 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: matrix_s
1008 TYPE(dft_control_type),
POINTER :: dft_control
1009 TYPE(qs_matrix_pools_type),
POINTER :: mpools
1010 TYPE(qs_scf_env_type),
POINTER :: scf_env
1011 TYPE(scf_control_type),
POINTER :: scf_control
1013 NULLIFY (scf_env, scf_control, maxao_maxmo_fm_pool, matrix_s, mpools, dft_control)
1014 CALL timeset(routinen, handle)
1016 cpassert(
ASSOCIATED(qs_env))
1019 IF (
PRESENT(n_col)) my_n_col = n_col
1022 scf_control=scf_control, &
1023 matrix_s=matrix_s, &
1024 dft_control=dft_control)
1025 CALL mpools_get(mpools, maxao_maxmo_fm_pool=maxao_maxmo_fm_pool)
1026 IF (
ASSOCIATED(scf_env))
THEN
1027 ortho_contains_cholesky = (scf_env%method /=
ot_method_nr) .AND. &
1028 (scf_env%cholesky_method > 0) .AND. &
1029 ASSOCIATED(scf_env%ortho)
1031 ortho_contains_cholesky = .false.
1034 CALL get_qs_env(qs_env, has_unit_metric=has_unit_metric)
1035 smearing_is_used = .false.
1036 IF (dft_control%smear)
THEN
1037 smearing_is_used = .true.
1040 IF (has_unit_metric)
THEN
1042 ELSE IF (smearing_is_used)
THEN
1044 matrix_s=matrix_s(1)%matrix)
1045 ELSE IF (ortho_contains_cholesky)
THEN
1047 ortho=scf_env%ortho)
1051 CALL timestop(handle)
1063 TYPE(qs_environment_type),
POINTER :: qs_env
1065 CHARACTER(len=*),
PARAMETER :: routinen =
'wfi_purge_history'
1067 INTEGER :: handle, io_unit, print_level
1068 TYPE(cp_logger_type),
POINTER :: logger
1069 TYPE(dft_control_type),
POINTER :: dft_control
1070 TYPE(qs_wf_history_type),
POINTER :: wf_history
1072 NULLIFY (dft_control, wf_history)
1074 CALL timeset(routinen, handle)
1076 print_level = logger%iter_info%print_level
1078 extension=
".scfLog")
1080 cpassert(
ASSOCIATED(qs_env))
1081 cpassert(
ASSOCIATED(qs_env%wf_history))
1082 cpassert(qs_env%wf_history%ref_count > 0)
1083 CALL get_qs_env(qs_env, dft_control=dft_control)
1085 SELECT CASE (qs_env%wf_history%interpolation_method_nr)
1093 IF (qs_env%wf_history%snapshot_count .GE. 2)
THEN
1094 IF (debug_this_module .AND. io_unit > 0) &
1095 WRITE (io_unit, fmt=
"(T2,A)")
"QS| Purging WFN history"
1096 CALL wfi_create(wf_history, interpolation_method_nr= &
1097 dft_control%qs_control%wf_interpolation_method_nr, &
1098 extrapolation_order=dft_control%qs_control%wf_extrapolation_order, &
1099 has_unit_metric=qs_env%has_unit_metric)
1101 wf_history=wf_history)
1103 CALL wfi_update(qs_env%wf_history, qs_env=qs_env, dt=1.0_dp)
1106 cpabort(
"Unknown extrapolation method.")
1108 CALL timestop(handle)
static GRID_HOST_DEVICE int modulo(int a, int m)
Equivalent of Fortran's MODULO, which always return a positive number. https://gcc....
collects all references to literature in CP2K as new algorithms / method are included from literature...
integer, save, public vandevondele2005a
integer, save, public kuhne2007
integer, save, public kolafa2004
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
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
basic linear algebra operations for full matrices
subroutine, public cp_fm_scale_and_add(alpha, matrix_a, beta, matrix_b)
calc A <- alpha*A + beta*B optimized for alpha == 1.0 (just add beta*B) and beta == 0....
subroutine, public cp_fm_scale(alpha, matrix_a)
scales a matrix matrix_a = alpha * matrix_b
pool for for elements that are retained and released
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
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_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)
...
integer, parameter, public low_print_level
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,...
Defines the basic variable types.
integer, parameter, public dp
Collection of simple mathematical functions and subroutines.
elemental real(kind=dp) function, public binomial(n, k)
The binomial coefficient n over k for 0 <= k <= n is calculated, otherwise zero is returned.
basic linear algebra operations for full matrixes
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 ...
collects routines that calculate density matrices
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 set_qs_env(qs_env, super_cell, mos, qmmm, qmmm_periodic, ewald_env, ewald_pw, mpools, rho_external, external_vxc, mask, scf_control, rel_control, qs_charges, ks_env, ks_qmmm_env, wf_history, scf_env, active_space, input, oce, rho_atom_set, rho0_atom_set, rho0_mpole, run_rtp, rtp, rhoz_set, rhoz_tot, ecoul_1c, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, efield, linres_control, xas_env, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, ls_scf_env, do_transport, transport_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, mp2_env, bs_env, kg_env, force, kpoints, WannierCentres, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Set the QUICKSTEP environment.
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)
collects routines that perform operations directly related to MOs
subroutine, public make_basis_simple(vmatrix, ncol)
given a set of vectors, return an orthogonal (C^T C == 1) set spanning the same space (notice,...
subroutine, public make_basis_lowdin(vmatrix, ncol, matrix_s)
return a set of S orthonormal vectors (C^T S C == 1) where a Loedwin transformation is applied to kee...
subroutine, public make_basis_cholesky(vmatrix, ncol, ortho)
return a set of S orthonormal vectors (C^T S C == 1) where the cholesky decomposed form of S is passe...
subroutine, public make_basis_sm(vmatrix, ncol, matrix_s)
returns an S-orthonormal basis v (v^T S v ==1)
Definition and initialisation of the mo data type.
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.
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_set(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)
...
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
module that contains the definitions of the scf types
integer, parameter, public ot_method_nr
Storage of past states of the qs_env. Methods to interpolate (or actually normally extrapolate) the n...
subroutine, public wfi_create_for_kp(wf_history)
...
subroutine, public wfi_purge_history(qs_env)
purges wf_history retaining only the latest snapshot
character(len=30) function, public wfi_get_method_label(method_nr)
returns a string describing the interpolation method
subroutine, public reorthogonalize_vectors(qs_env, v_matrix, n_col)
reorthogonalizes the mos
subroutine, public wfi_update(wf_history, qs_env, dt)
updates the snapshot buffer, taking a new snapshot
subroutine, public wfi_extrapolate(wf_history, qs_env, dt, extrapolation_method_nr, orthogonal_wf)
calculates the new starting state for the scf for the next wf optimization
subroutine, public wfi_create(wf_history, interpolation_method_nr, extrapolation_order, has_unit_metric)
...
interpolate the wavefunctions to speed up the convergence when doing MD
type(qs_wf_snapshot_type) function, pointer, public wfi_get_snapshot(wf_history, wf_index)
returns a snapshot, the first being the latest snapshot
subroutine, public wfi_release(wf_history)
releases a wf_history of a wavefunction (see doc/ReferenceCounting.html)
parameters that control an scf iteration