28 USE dbcsr_api,
ONLY: dbcsr_p_type,&
84#include "./base/base_uses.f90"
90 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_scf_output'
110 INTEGER,
INTENT(IN) :: output_unit
113 INTEGER :: nelectron_total
114 LOGICAL :: gapw, gapw_xc, qmmm
121 NULLIFY (rho, energy, dft_control, scf_env, qs_charges)
122 CALL get_qs_env(qs_env=qs_env, rho=rho, energy=energy, dft_control=dft_control, &
123 scf_env=scf_env, qs_charges=qs_charges)
125 gapw = dft_control%qs_control%gapw
126 gapw_xc = dft_control%qs_control%gapw_xc
128 nelectron_total = scf_env%nelectron
130 CALL qs_scf_print_scf_summary(output_unit, rho, qs_charges, energy, nelectron_total, &
131 dft_control, qmmm, qs_env, gapw, gapw_xc)
142 INTEGER :: output_unit
146 INTEGER :: homo, ispin, nao, nelectron_spin, nmo
148 IF (output_unit > 0)
THEN
149 DO ispin = 1, dft_control%nspins
152 nelectron=nelectron_spin, &
155 IF (dft_control%nspins > 1)
THEN
156 WRITE (unit=output_unit, fmt=
"(/,T2,A,I2)")
"Spin", ispin
158 WRITE (unit=output_unit, fmt=
"(/,(T2,A,T71,I10))") &
159 "Number of electrons:", nelectron_spin, &
160 "Number of occupied orbitals:", homo, &
161 "Number of molecular orbitals:", nmo
163 WRITE (unit=output_unit, fmt=
"(/,T2,A,T71,I10)") &
164 "Number of orbital functions:", nao
180 LOGICAL,
INTENT(IN) :: final_mos
182 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_write_mos'
184 CHARACTER(LEN=2) :: solver_method
185 CHARACTER(LEN=3*default_string_length) :: message
186 CHARACTER(LEN=5) :: spin
187 CHARACTER(LEN=default_string_length), &
188 DIMENSION(:),
POINTER :: tmpstringlist
189 INTEGER :: handle, homo, ikp, ispin, iw, kpoint, &
190 nao, nelectron, nkp, nmo, nspin, numo
191 INTEGER,
DIMENSION(2) :: nmos_occ
192 INTEGER,
DIMENSION(:),
POINTER :: mo_index_range
193 LOGICAL :: do_kpoints, print_eigvals, &
194 print_eigvecs, print_mo_info, &
195 print_occup, print_occup_stats
196 REAL(kind=
dp) :: flexible_electron_count, maxocc, n_el_f, &
197 occup_stats_occ_threshold
198 REAL(kind=
dp),
DIMENSION(:),
POINTER :: mo_eigenvalues, umo_eigenvalues
203 TYPE(
cp_fm_type),
POINTER :: mo_coeff, umo_coeff
205 TYPE(dbcsr_p_type),
DIMENSION(:),
POINTER :: ks, s
206 TYPE(dbcsr_type),
POINTER :: matrix_ks, matrix_s
214 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
218 CALL timeset(routinen, handle)
220 cpassert(
ASSOCIATED(qs_env))
224 atomic_kind_set=atomic_kind_set, &
225 blacs_env=blacs_env, &
226 dft_control=dft_control, &
227 do_kpoints=do_kpoints, &
229 qs_kind_set=qs_kind_set, &
231 particle_set=particle_set, &
232 scf_control=scf_control)
239 CALL section_vals_val_get(dft_section,
"PRINT%MO%OCCUPATION_NUMBERS_STATS", c_vals=tmpstringlist)
241 print_occup_stats = .false.
242 occup_stats_occ_threshold = 1e-6_dp
243 IF (
SIZE(tmpstringlist) > 0)
THEN
244 print_occup_stats = .true.
245 IF (len_trim(tmpstringlist(1)) > 0) &
246 READ (tmpstringlist(1), *) print_occup_stats
248 IF (
SIZE(tmpstringlist) > 1) &
249 READ (tmpstringlist(2), *) occup_stats_occ_threshold
254 IF ((.NOT. print_mo_info) .OR. (.NOT. (print_eigvals .OR. print_eigvecs .OR. print_occup .OR. print_occup_stats)))
THEN
255 CALL timestop(handle)
259 NULLIFY (fm_struct_tmp)
261 NULLIFY (mo_eigenvalues)
264 NULLIFY (umo_eigenvalues)
267 nspin = dft_control%nspins
273 nkp =
SIZE(kpoints%kp_env)
275 CALL get_qs_env(qs_env, matrix_ks=ks, matrix_s=s)
276 cpassert(
ASSOCIATED(ks))
277 cpassert(
ASSOCIATED(s))
284 mos => kpoints%kp_env(ikp)%kpoint_env%mos(1, :)
287 CALL get_qs_env(qs_env, matrix_ks=ks, mos=mos)
290 cpassert(
ASSOCIATED(mos))
301 cpabort(
"The OT method is not implemented for k points")
304 matrix_ks => ks(ispin)%matrix
305 matrix_s => s(1)%matrix
308 IF (dft_control%do_admm)
THEN
316 eigenvalues=mo_eigenvalues, &
319 nelectron=nelectron, &
323 flexible_electron_count=flexible_electron_count)
327 cpassert(
ASSOCIATED(mo_index_range))
328 numo = min(mo_index_range(2) - homo, nao - homo)
330 IF (.NOT. final_mos)
THEN
332 message =
"The MO information for unoccupied MOs is only calculated after "// &
333 "SCF convergence is achieved when the orbital transformation (OT) "// &
335 cpwarn(trim(message))
354 flexible_electron_count=flexible_electron_count)
356 fm_struct=fm_struct_tmp, &
357 name=
"Temporary MO set (unoccupied MOs only)for printout")
360 mo_coeff=umo_coeff, &
361 eigenvalues=umo_eigenvalues)
368 NULLIFY (local_preconditioner)
369 IF (
ASSOCIATED(scf_env%ot_preconditioner))
THEN
370 local_preconditioner => scf_env%ot_preconditioner(1)%preconditioner
372 NULLIFY (local_preconditioner)
377 matrix_c_fm=umo_coeff, &
378 matrix_orthogonal_space_fm=mo_coeff, &
379 eps_gradient=scf_control%eps_lumos, &
381 iter_max=scf_control%max_iter_lumos, &
382 size_ortho_space=nmo)
386 ks_matrix=matrix_ks, &
387 evals_arg=umo_eigenvalues, &
392 IF (dft_control%do_admm)
THEN
418 cpabort(
"Invalid spin")
420 IF (
ASSOCIATED(umo_set))
THEN
422 dft_section, 4, kpoint, final_mos=final_mos, spin=trim(spin), &
423 solver_method=solver_method, umo_set=umo_set)
426 dft_section, 4, kpoint, final_mos=final_mos, spin=trim(spin), &
427 solver_method=solver_method)
430 IF (
ASSOCIATED(umo_set))
THEN
432 dft_section, 4, kpoint, final_mos=final_mos, &
433 solver_method=solver_method, umo_set=umo_set)
436 dft_section, 4, kpoint, final_mos=final_mos, &
437 solver_method=solver_method)
441 nmos_occ(ispin) = max(nmos_occ(ispin), count(mo_set%occupation_numbers > occup_stats_occ_threshold))
445 IF (
ASSOCIATED(umo_set))
THEN
458 IF (print_mo_info .AND. print_occup_stats)
THEN
460 ignore_should_output=print_mo_info, &
463 IF (
SIZE(mos) > 1)
THEN
464 WRITE (unit=iw, fmt=
"(A,I4)")
" MO| Total occupied (ALPHA):", nmos_occ(1)
465 WRITE (unit=iw, fmt=
"(A,I4)")
" MO| Total occupied (BETA): ", nmos_occ(2)
467 WRITE (unit=iw, fmt=
"(A,I4)")
" MO| Total occupied: ", nmos_occ(1)
469 WRITE (unit=iw, fmt=
"(A)")
""
472 ignore_should_output=print_mo_info)
475 CALL timestop(handle)
490 energy, total_steps, should_stop, outer_loop_converged)
491 INTEGER :: output_unit
495 INTEGER :: total_steps
496 LOGICAL,
INTENT(IN) :: should_stop, outer_loop_converged
498 REAL(kind=
dp) :: outer_loop_eps
500 outer_loop_eps = sqrt(maxval(scf_env%outer_scf%gradient(:, scf_env%outer_scf%iter_count)**2))
501 IF (output_unit > 0)
WRITE (output_unit,
'(/,T3,A,I4,A,E10.2,A,F22.10)') &
502 "outer SCF iter = ", scf_env%outer_scf%iter_count, &
503 " RMS gradient = ", outer_loop_eps,
" energy =", energy%total
505 IF (outer_loop_converged)
THEN
506 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
507 "outer SCF loop converged in", scf_env%outer_scf%iter_count, &
508 " iterations or ", total_steps,
" steps"
509 ELSE IF (scf_env%outer_scf%iter_count > scf_control%outer_scf%max_scf &
510 .OR. should_stop)
THEN
511 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
512 "outer SCF loop FAILED to converge after ", &
513 scf_env%outer_scf%iter_count,
" iterations or ", total_steps,
" steps"
530 INTEGER :: output_unit
531 LOGICAL :: just_energy
532 REAL(kind=
dp) :: t1, t2
535 IF ((output_unit > 0) .AND. scf_env%print_iter_line)
THEN
536 IF (just_energy)
THEN
537 WRITE (unit=output_unit, &
538 fmt=
"(T2,I5,1X,A,T20,E8.2,1X,F6.1,16X,F20.10)") &
539 scf_env%iter_count, trim(scf_env%iter_method), scf_env%iter_param, &
540 t2 - t1, energy%total
542 IF ((abs(scf_env%iter_delta) < 1.0e-8_dp) .OR. &
543 (abs(scf_env%iter_delta) >= 1.0e5_dp))
THEN
544 WRITE (unit=output_unit, &
545 fmt=
"(T2,I5,1X,A,T20,E8.2,1X,F6.1,1X,ES14.4,1X,F20.10,1X,ES9.2)") &
546 scf_env%iter_count, trim(scf_env%iter_method), scf_env%iter_param, &
547 t2 - t1, scf_env%iter_delta, energy%total, energy%total - energy%tot_old
549 WRITE (unit=output_unit, &
550 fmt=
"(T2,I5,1X,A,T20,E8.2,1X,F6.1,1X,F14.8,1X,F20.10,1X,ES9.2)") &
551 scf_env%iter_count, trim(scf_env%iter_method), scf_env%iter_param, &
552 t2 - t1, scf_env%iter_delta, energy%total, energy%total - energy%tot_old
577 SUBROUTINE qs_scf_print_scf_summary(output_unit, rho, qs_charges, energy, nelectron_total, &
578 dft_control, qmmm, qs_env, gapw, gapw_xc)
579 INTEGER,
INTENT(IN) :: output_unit
583 INTEGER,
INTENT(IN) :: nelectron_total
585 LOGICAL,
INTENT(IN) :: qmmm
587 LOGICAL,
INTENT(IN) :: gapw, gapw_xc
589 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_print_scf_summary'
591 INTEGER :: bc, handle, ispin, psolver
592 REAL(kind=
dp) :: exc1_energy, exc_energy, &
593 implicit_ps_ehartree, tot1_h, tot1_s
594 REAL(kind=
dp),
DIMENSION(:),
POINTER :: tot_rho_r
597 NULLIFY (tot_rho_r, pw_env)
598 CALL timeset(routinen, handle)
601 psolver = pw_env%poisson_env%parameters%solver
603 IF (output_unit > 0)
THEN
605 IF (.NOT. (dft_control%qs_control%semi_empirical .OR. &
606 dft_control%qs_control%xtb .OR. &
607 dft_control%qs_control%dftb))
THEN
608 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T41,2F20.10))") &
609 "Electronic density on regular grids: ", &
612 "Core density on regular grids:", &
613 qs_charges%total_rho_core_rspace, &
614 qs_charges%total_rho_core_rspace - real(nelectron_total + dft_control%charge,
dp)
616 IF (dft_control%correct_surf_dip)
THEN
617 WRITE (unit=output_unit, fmt=
"((T3,A,/,T3,A,T41,F20.10))") &
618 "Total dipole moment perpendicular to ", &
619 "the slab [electrons-Angstroem]: ", &
620 qs_env%surface_dipole_moment
624 tot1_h = qs_charges%total_rho1_hard(1)
625 tot1_s = qs_charges%total_rho1_soft(1)
626 DO ispin = 2, dft_control%nspins
627 tot1_h = tot1_h + qs_charges%total_rho1_hard(ispin)
628 tot1_s = tot1_s + qs_charges%total_rho1_soft(ispin)
630 WRITE (unit=output_unit, fmt=
"((T3,A,T41,2F20.10))") &
631 "Hard and soft densities (Lebedev):", &
633 WRITE (unit=output_unit, fmt=
"(T3,A,T41,F20.10)") &
634 "Total Rho_soft + Rho1_hard - Rho1_soft (r-space): ", &
636 "Total charge density (r-space): ", &
638 + qs_charges%total_rho_core_rspace, &
639 "Total Rho_soft + Rho0_soft (g-space):", &
640 qs_charges%total_rho_gspace
642 WRITE (unit=output_unit, fmt=
"(T3,A,T41,F20.10)") &
643 "Total charge density on r-space grids: ", &
645 qs_charges%total_rho_core_rspace, &
646 "Total charge density g-space grids: ", &
647 qs_charges%total_rho_gspace
650 IF (dft_control%qs_control%semi_empirical)
THEN
651 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
652 "Core-core repulsion energy [eV]: ", energy%core_overlap*
evolt, &
653 "Core Hamiltonian energy [eV]: ", energy%core*
evolt, &
654 "Two-electron integral energy [eV]: ", energy%hartree*
evolt, &
655 "Electronic energy [eV]: ", &
656 (energy%core + 0.5_dp*energy%hartree)*
evolt
657 IF (energy%dispersion /= 0.0_dp) &
658 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
659 "Dispersion energy [eV]: ", energy%dispersion*
evolt
660 ELSEIF (dft_control%qs_control%dftb)
THEN
661 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
662 "Core Hamiltonian energy: ", energy%core, &
663 "Repulsive potential energy: ", energy%repulsive, &
664 "Electronic energy: ", energy%hartree, &
665 "Dispersion energy: ", energy%dispersion
666 IF (energy%dftb3 /= 0.0_dp) &
667 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
668 "DFTB3 3rd order energy: ", energy%dftb3
669 IF (energy%efield /= 0.0_dp) &
670 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
671 "Electric field interaction energy: ", energy%efield
672 ELSEIF (dft_control%qs_control%xtb)
THEN
673 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
674 "Core Hamiltonian energy: ", energy%core, &
675 "Repulsive potential energy: ", energy%repulsive, &
676 "Electronic energy: ", energy%hartree, &
677 "DFTB3 3rd order energy: ", energy%dftb3, &
678 "Dispersion energy: ", energy%dispersion
679 IF (dft_control%qs_control%xtb_control%xb_interaction) &
680 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
681 "Correction for halogen bonding: ", energy%xtb_xb_inter
682 IF (dft_control%qs_control%xtb_control%do_nonbonded) &
683 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
684 "Correction for nonbonded interactions: ", energy%xtb_nonbonded
685 IF (energy%efield /= 0.0_dp) &
686 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
687 "Electric field interaction energy: ", energy%efield
689 IF (dft_control%do_admm)
THEN
690 exc_energy = energy%exc + energy%exc_aux_fit
691 IF (gapw .OR. gapw_xc) exc1_energy = energy%exc1 + energy%exc1_aux_fit
693 exc_energy = energy%exc
694 IF (gapw .OR. gapw_xc) exc1_energy = energy%exc1
698 implicit_ps_ehartree = pw_env%poisson_env%implicit_env%ehartree
699 bc = pw_env%poisson_env%parameters%ps_implicit_params%boundary_condition
702 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
703 "Overlap energy of the core charge distribution:", energy%core_overlap, &
704 "Self energy of the core charge distribution: ", energy%core_self, &
705 "Core Hamiltonian energy: ", energy%core, &
706 "Hartree energy: ", implicit_ps_ehartree, &
707 "Electric enthalpy: ", energy%hartree, &
708 "Exchange-correlation energy: ", exc_energy
710 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
711 "Overlap energy of the core charge distribution:", energy%core_overlap, &
712 "Self energy of the core charge distribution: ", energy%core_self, &
713 "Core Hamiltonian energy: ", energy%core, &
714 "Hartree energy: ", energy%hartree, &
715 "Exchange-correlation energy: ", exc_energy
718 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
719 "Overlap energy of the core charge distribution:", energy%core_overlap, &
720 "Self energy of the core charge distribution: ", energy%core_self, &
721 "Core Hamiltonian energy: ", energy%core, &
722 "Hartree energy: ", energy%hartree, &
723 "Exchange-correlation energy: ", exc_energy
725 IF (energy%e_hartree /= 0.0_dp) &
726 WRITE (unit=output_unit, fmt=
"(T3,A,/,T3,A,T56,F25.14)") &
727 "Coulomb Electron-Electron Interaction Energy ", &
728 "- Already included in the total Hartree term ", energy%e_hartree
729 IF (energy%ex /= 0.0_dp) &
730 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
731 "Hartree-Fock Exchange energy: ", energy%ex
732 IF (energy%dispersion /= 0.0_dp) &
733 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
734 "Dispersion energy: ", energy%dispersion
735 IF (energy%gcp /= 0.0_dp) &
736 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
737 "gCP energy: ", energy%gcp
738 IF (energy%efield /= 0.0_dp) &
739 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
740 "Electric field interaction energy: ", energy%efield
742 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
743 "GAPW| Exc from hard and soft atomic rho1: ", exc1_energy, &
744 "GAPW| local Eh = 1 center integrals: ", energy%hartree_1c
747 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
748 "GAPW_XC| Exc from hard and soft atomic rho1: ", exc1_energy
751 IF (dft_control%smear)
THEN
752 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
753 "Electronic entropic energy:", energy%kTS
754 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
755 "Fermi energy:", energy%efermi
757 IF (dft_control%dft_plus_u)
THEN
758 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
759 "DFT+U energy:", energy%dft_plus_u
761 IF (dft_control%do_sccs)
THEN
762 WRITE (unit=output_unit, fmt=
"(A)")
""
766 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
767 "QM/MM Electrostatic energy: ", energy%qmmm_el
768 IF (qs_env%qmmm_env_qm%image_charge)
THEN
769 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
770 "QM/MM image charge energy: ", energy%image_charge
773 IF (dft_control%qs_control%mulliken_restraint)
THEN
774 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
775 "Mulliken restraint energy: ", energy%mulliken
777 IF (dft_control%qs_control%semi_empirical)
THEN
778 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
779 "Total energy [eV]: ", energy%total*
evolt
780 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
781 "Atomic reference energy [eV]: ", energy%core_self*
evolt, &
782 "Heat of formation [kcal/mol]: ", &
783 (energy%total + energy%core_self)*
kcalmol
785 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
786 "Total energy: ", energy%total
789 IF (qs_env%qmmm_env_qm%image_charge)
THEN
796 CALL timestop(handle)
798 END SUBROUTINE qs_scf_print_scf_summary
813 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_loop_print'
815 INTEGER :: after, handle, ic, ispin, iw
816 LOGICAL :: do_kpoints, omit_headers
817 REAL(kind=
dp) :: mo_mag_max, mo_mag_min, orthonormality
819 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_ks, matrix_p, matrix_s
826 CALL timeset(routinen, handle)
828 CALL get_qs_env(qs_env=qs_env, input=input, dft_control=dft_control, &
829 do_kpoints=do_kpoints)
835 DO ispin = 1, dft_control%nspins
838 dft_section,
"PRINT%AO_MATRICES/DENSITY"),
cp_p_file))
THEN
844 after = min(max(after, 1), 16)
845 DO ic = 1,
SIZE(matrix_p, 2)
847 output_unit=iw, omit_headers=omit_headers)
850 "PRINT%AO_MATRICES/DENSITY")
854 dft_section,
"PRINT%AO_MATRICES/KOHN_SHAM_MATRIX"),
cp_p_file))
THEN
858 after = min(max(after, 1), 16)
859 CALL get_qs_env(qs_env=qs_env, matrix_ks_kp=matrix_ks)
860 DO ic = 1,
SIZE(matrix_ks, 2)
861 IF (dft_control%qs_control%semi_empirical)
THEN
863 scale=
evolt, output_unit=iw, omit_headers=omit_headers)
866 output_unit=iw, omit_headers=omit_headers)
870 "PRINT%AO_MATRICES/KOHN_SHAM_MATRIX")
876 scf_section,
"PRINT%MO_ORTHONORMALITY"),
cp_p_file))
THEN
881 WRITE (iw,
'(T8,A)') &
882 " K-points: Maximum deviation from MO S-orthonormality not determined"
885 "PRINT%MO_ORTHONORMALITY")
891 CALL get_qs_env(qs_env=qs_env, matrix_s_kp=matrix_s)
897 WRITE (iw,
'(T8,A,T61,E20.4)') &
898 " Maximum deviation from MO S-orthonormality", orthonormality
901 "PRINT%MO_ORTHONORMALITY")
905 scf_section,
"PRINT%MO_MAGNITUDE"),
cp_p_file))
THEN
910 WRITE (iw,
'(T8,A)') &
911 " K-points: Minimum/Maximum MO magnitude not determined"
914 "PRINT%MO_MAGNITUDE")
921 WRITE (iw,
'(T8,A,T41,2E20.4)') &
922 " Minimum/Maximum MO magnitude ", mo_mag_min, mo_mag_max
925 "PRINT%MO_MAGNITUDE")
929 CALL timestop(handle)
948 energy, total_steps, should_stop, outer_loop_converged, &
950 INTEGER :: output_unit
955 INTEGER :: total_steps
956 LOGICAL,
INTENT(IN) :: should_stop, outer_loop_converged, &
959 REAL(kind=
dp) :: outer_loop_eps
962 outer_loop_eps = sqrt(maxval(scf_env%outer_scf%gradient(:, scf_env%outer_scf%iter_count)**2))
963 IF (output_unit > 0)
WRITE (output_unit,
'(/,T3,A,I4,A,E10.2,A,F22.10)') &
964 "CDFT SCF iter = ", scf_env%outer_scf%iter_count, &
965 " RMS gradient = ", outer_loop_eps,
" energy =", energy%total
966 IF (outer_loop_converged)
THEN
967 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
968 "CDFT SCF loop converged in", scf_env%outer_scf%iter_count, &
969 " iterations or ", total_steps,
" steps"
971 IF ((scf_env%outer_scf%iter_count > scf_control%outer_scf%max_scf .OR. should_stop) &
972 .AND. .NOT. outer_loop_converged)
THEN
973 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
974 "CDFT SCF loop FAILED to converge after ", &
975 scf_env%outer_scf%iter_count,
" iterations or ", total_steps,
" steps"
990 INTEGER :: output_unit
993 IF (output_unit > 0)
THEN
994 WRITE (output_unit,
'(/,A)') &
995 " ---------------------------------- CDFT --------------------------------------"
996 WRITE (output_unit,
'(A)') &
997 " Optimizing a density constraint in an external SCF loop "
998 WRITE (output_unit,
'(A)')
" "
999 SELECT CASE (cdft_control%type)
1001 WRITE (output_unit,
'(A)')
" Type of constraint: Hirshfeld"
1003 WRITE (output_unit,
'(A)')
" Type of constraint: Becke"
1005 WRITE (output_unit,
'(A,I8)')
" Number of constraints: ",
SIZE(cdft_control%group)
1006 WRITE (output_unit,
'(A,L8)')
" Using fragment densities:", cdft_control%fragment_density
1007 WRITE (output_unit,
'(A)')
" "
1008 IF (cdft_control%atomic_charges)
WRITE (output_unit,
'(A,/)')
" Calculating atomic CDFT charges"
1009 SELECT CASE (cdft_control%constraint_control%optimizer)
1011 WRITE (output_unit,
'(A)') &
1012 " Minimizer : SD : steepest descent"
1014 WRITE (output_unit,
'(A)') &
1015 " Minimizer : DIIS : direct inversion"
1016 WRITE (output_unit,
'(A)') &
1017 " in the iterative subspace"
1018 WRITE (output_unit,
'(A,I3,A)') &
1020 cdft_control%constraint_control%diis_buffer_length,
" DIIS vectors"
1022 WRITE (output_unit,
'(A)') &
1023 " Minimizer : BISECT : gradient bisection"
1024 WRITE (output_unit,
'(A,I3)') &
1025 " using a trust count of", &
1026 cdft_control%constraint_control%bisect_trust_count
1030 cdft_control%constraint_control%optimizer, output_unit)
1032 WRITE (output_unit,
'(A)')
" Minimizer : Secant"
1036 WRITE (output_unit,
'(/,A,L7)') &
1037 " Reusing OT preconditioner: ", cdft_control%reuse_precond
1038 IF (cdft_control%reuse_precond)
THEN
1039 WRITE (output_unit,
'(A,I3,A,I3,A)') &
1040 " using old preconditioner for up to ", &
1041 cdft_control%max_reuse,
" subsequent CDFT SCF"
1042 WRITE (output_unit,
'(A,I3,A,I3,A)') &
1043 " iterations if the relevant loop converged in less than ", &
1044 cdft_control%precond_freq,
" steps"
1046 SELECT CASE (cdft_control%type)
1048 WRITE (output_unit,
'(/,A)')
" Hirshfeld constraint settings"
1049 WRITE (output_unit,
'(A)')
" "
1050 SELECT CASE (cdft_control%hirshfeld_control%shape_function)
1052 WRITE (output_unit,
'(A, A8)') &
1053 " Shape function type: ",
"Gaussian"
1054 WRITE (output_unit,
'(A)', advance=
'NO') &
1055 " Type of Gaussian: "
1056 SELECT CASE (cdft_control%hirshfeld_control%gaussian_shape)
1058 WRITE (output_unit,
'(A13)')
"Default"
1060 WRITE (output_unit,
'(A13)')
"Covalent"
1062 WRITE (output_unit,
'(A13)')
"Fixed radius"
1064 WRITE (output_unit,
'(A13)')
"Van der Waals"
1066 WRITE (output_unit,
'(A13)')
"User-defined"
1070 WRITE (output_unit,
'(A, A8)') &
1071 " Shape function type: ",
"Density"
1074 WRITE (output_unit,
'(/, A)')
" Becke constraint settings"
1075 WRITE (output_unit,
'(A)')
" "
1076 SELECT CASE (cdft_control%becke_control%cutoff_type)
1078 WRITE (output_unit,
'(A,F8.3,A)') &
1079 " Cutoff for partitioning :",
cp_unit_from_cp2k(cdft_control%becke_control%rglobal, &
1080 "angstrom"),
" angstrom"
1082 WRITE (output_unit,
'(A)') &
1083 " Using element specific cutoffs for partitioning"
1085 WRITE (output_unit,
'(A,L7)') &
1086 " Skipping distant gpoints: ", cdft_control%becke_control%should_skip
1087 WRITE (output_unit,
'(A,L7)') &
1088 " Precompute gradients : ", cdft_control%becke_control%in_memory
1089 WRITE (output_unit,
'(A)')
" "
1090 IF (cdft_control%becke_control%adjust) &
1091 WRITE (output_unit,
'(A)') &
1092 " Using atomic radii to generate a heteronuclear charge partitioning"
1093 WRITE (output_unit,
'(A)')
" "
1094 IF (.NOT. cdft_control%becke_control%cavity_confine)
THEN
1095 WRITE (output_unit,
'(A)') &
1096 " No confinement is active"
1098 WRITE (output_unit,
'(A)')
" Confinement using a Gaussian shaped cavity is active"
1099 SELECT CASE (cdft_control%becke_control%cavity_shape)
1101 WRITE (output_unit,
'(A,F8.4, A)') &
1102 " Type of Gaussian : Fixed radius: ", &
1105 WRITE (output_unit,
'(A)') &
1106 " Type of Gaussian : Covalent radius "
1108 WRITE (output_unit,
'(A)') &
1109 " Type of Gaussian : vdW radius "
1111 WRITE (output_unit,
'(A)') &
1112 " Type of Gaussian : User radius "
1114 WRITE (output_unit,
'(A,ES12.4)') &
1115 " Cavity threshold : ", cdft_control%becke_control%eps_cavity
1118 WRITE (output_unit,
'(/,A)') &
1119 " ---------------------------------- CDFT --------------------------------------"
1132 INTEGER :: output_unit
1137 IF (output_unit > 0)
THEN
1138 SELECT CASE (cdft_control%type)
1140 WRITE (output_unit,
'(/,T3,A,T60)') &
1141 '------------------- Hirshfeld constraint information -------------------'
1143 WRITE (output_unit,
'(/,T3,A,T60)') &
1144 '--------------------- Becke constraint information ---------------------'
1146 cpabort(
"Unknown CDFT constraint.")
1148 DO igroup = 1,
SIZE(cdft_control%target)
1149 IF (igroup > 1)
WRITE (output_unit,
'(T3,A)')
' '
1150 WRITE (output_unit,
'(T3,A,T54,(3X,I18))') &
1151 'Atomic group :', igroup
1152 SELECT CASE (cdft_control%group(igroup)%constraint_type)
1154 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1155 WRITE (output_unit,
'(T3,A,T42,A)') &
1156 'Type of constraint :', adjustr(
'Charge density constraint (frag.)')
1158 WRITE (output_unit,
'(T3,A,T50,A)') &
1159 'Type of constraint :', adjustr(
'Charge density constraint')
1162 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1163 WRITE (output_unit,
'(T3,A,T35,A)') &
1164 'Type of constraint :', adjustr(
'Magnetization density constraint (frag.)')
1166 WRITE (output_unit,
'(T3,A,T43,A)') &
1167 'Type of constraint :', adjustr(
'Magnetization density constraint')
1170 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1171 WRITE (output_unit,
'(T3,A,T38,A)') &
1172 'Type of constraint :', adjustr(
'Alpha spin density constraint (frag.)')
1174 WRITE (output_unit,
'(T3,A,T46,A)') &
1175 'Type of constraint :', adjustr(
'Alpha spin density constraint')
1178 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1179 WRITE (output_unit,
'(T3,A,T39,A)') &
1180 'Type of constraint :', adjustr(
'Beta spin density constraint (frag.)')
1182 WRITE (output_unit,
'(T3,A,T47,A)') &
1183 'Type of constraint :', adjustr(
'Beta spin density constraint')
1186 cpabort(
"Unknown constraint type.")
1188 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1189 'Target value of constraint :', cdft_control%target(igroup)
1190 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1191 'Current value of constraint :', cdft_control%value(igroup)
1192 WRITE (output_unit,
'(T3,A,T59,(3X,ES13.3))') &
1193 'Deviation from target :', cdft_control%value(igroup) - cdft_control%target(igroup)
1194 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1195 'Strength of constraint :', cdft_control%strength(igroup)
1197 WRITE (output_unit,
'(T3,A)') &
1198 '------------------------------------------------------------------------'
Types and set/get functions for auxiliary density matrix methods.
Contains methods used in the context of density fitting.
subroutine, public admm_uncorrect_for_eigenvalues(ispin, admm_env, ks_matrix)
...
subroutine, public admm_correct_for_eigenvalues(ispin, admm_env, ks_matrix)
...
Define the atomic kind types and their sub types.
methods related to the blacs parallel environment
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public cp_dbcsr_write_sparse_matrix(sparse_matrix, before, after, qs_env, para_env, first_row, last_row, first_col, last_col, scale, output_unit, omit_headers)
...
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_init_random(matrix, ncol, start_col)
fills a matrix with random numbers
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,...
integer, parameter, public cp_p_file
integer function, public cp_print_key_should_output(iteration_info, basis_section, print_key_path, used_print_key, first_time)
returns what should be done with the given property if btest(res,cp_p_store) then the property should...
real(kind=dp) function, public cp_unit_from_cp2k(value, unit_str, defaults, power)
converts from the internal cp2k units to the given unit
sums arrays of real/complex numbers with much reduced round-off as compared to a naive implementation...
Defines the basic variable types.
integer, parameter, public dp
integer, parameter, public default_string_length
Types and basic routines needed for a kpoint calculation.
Machine interface based on Fortran 2003 and POSIX.
subroutine, public m_flush(lunit)
flushes units if the &GLOBAL flag is set accordingly
Interface to the message passing library MPI.
Define the data structure for the particle information.
Definition of physical constants:
real(kind=dp), parameter, public kcalmol
real(kind=dp), parameter, public evolt
computes preconditioners, and implements methods to apply them currently used in qs_ot
Types containing essential information for running implicit (iterative) Poisson solver.
integer, parameter, public neumann_bc
integer, parameter, public mixed_bc
integer, parameter, public mixed_periodic_bc
integer, parameter, public periodic_bc
container for various plainwaves related things
functions related to the poisson solver on regular grids
integer, parameter, public pw_poisson_implicit
Routines for image charge calculation within QM/MM.
subroutine, public print_image_coefficients(image_coeff, qs_env)
Print image coefficients.
Control parameters for optimizers that work with CDFT constraints.
subroutine, public cdft_opt_type_write(cdft_opt_control, optimizer, output_unit)
writes information about the CDFT optimizer object
Defines CDFT control structures.
container for information about total charges on the grids
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.
Define the quickstep kind type and their sub types.
Definition and initialisation of the mo data type.
subroutine, public write_mo_set_to_output_unit(mo_set, atomic_kind_set, qs_kind_set, particle_set, dft_section, before, kpoint, final_mos, spin, solver_method, rtp, cpart, sim_step, umo_set)
Write MO information to output file (eigenvalues, occupation numbers, coefficients)
collects routines that perform operations directly related to MOs
subroutine, public calculate_magnitude(mo_array, mo_mag_min, mo_mag_max)
...
subroutine, public calculate_orthonormality(orthonormality, mo_array, matrix_s)
...
Set occupation of molecular orbitals.
Definition and initialisation of the mo data type.
subroutine, public allocate_mo_set(mo_set, nao, nmo, nelectron, n_el_f, maxocc, flexible_electron_count)
Allocates a mo set and partially initializes it (nao,nmo,nelectron, and flexible_electron_count are v...
subroutine, public deallocate_mo_set(mo_set)
Deallocate a wavefunction data structure.
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.
subroutine, public init_mo_set(mo_set, fm_pool, fm_ref, fm_struct, name)
initializes an allocated mo_set. eigenvalues, mo_coeff, occupation_numbers are valid only after this ...
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)
...
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...
Self-consistent continuum solvation (SCCS) model implementation.
subroutine, public print_sccs_results(energy, sccs_control, output_unit)
Print SCCS results.
subroutine, public qs_scf_print_summary(output_unit, qs_env)
writes a summary of information after scf
subroutine, public qs_scf_cdft_constraint_info(output_unit, cdft_control)
writes CDFT constraint information
subroutine, public qs_scf_loop_print(qs_env, scf_env, para_env)
collects the 'heavy duty' printing tasks out of the SCF loop
subroutine, public qs_scf_initial_info(output_unit, mos, dft_control)
writes basic information at the beginning of an scf run
subroutine, public qs_scf_outer_loop_info(output_unit, scf_control, scf_env, energy, total_steps, should_stop, outer_loop_converged)
writes basic information obtained in a scf outer loop step
subroutine, public qs_scf_write_mos(qs_env, scf_env, final_mos)
Write the MO eigenvector, eigenvalues, and occupation numbers to the output unit.
subroutine, public qs_scf_loop_info(scf_env, output_unit, just_energy, t1, t2, energy)
writes basic information obtained in a scf step
subroutine, public qs_scf_cdft_info(output_unit, scf_control, scf_env, cdft_control, energy, total_steps, should_stop, outer_loop_converged, cdft_loop)
writes CDFT constraint information and optionally CDFT scf loop info
subroutine, public qs_scf_cdft_initial_info(output_unit, cdft_control)
writes information about the CDFT env
module that contains the definitions of the scf types
integer, parameter, public ot_method_nr
integer, parameter, public special_diag_method_nr
parameters that control an scf iteration
stores some data used in wavefunction fitting
Provides all information about an atomic kind.
represent a blacs multidimensional parallel environment (for the mpi corrispective see cp_paratypes/m...
keeps the information about the structure of a full matrix
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Contains information about kpoints.
stores all the informations relevant to an mpi environment
contained for different pw related things
Container for information about total charges on the grids.
Provides all information about a quickstep kind.
keeps the density in various representations, keeping track of which ones are valid.