83#include "./base/base_uses.f90"
89 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_scf_output'
109 INTEGER,
INTENT(IN) :: output_unit
112 INTEGER :: nelectron_total
113 LOGICAL :: gapw, gapw_xc, qmmm
120 NULLIFY (rho, energy, dft_control, scf_env, qs_charges)
121 CALL get_qs_env(qs_env=qs_env, rho=rho, energy=energy, dft_control=dft_control, &
122 scf_env=scf_env, qs_charges=qs_charges)
124 gapw = dft_control%qs_control%gapw
125 gapw_xc = dft_control%qs_control%gapw_xc
127 nelectron_total = scf_env%nelectron
129 CALL qs_scf_print_scf_summary(output_unit, rho, qs_charges, energy, nelectron_total, &
130 dft_control, qmmm, qs_env, gapw, gapw_xc)
142 INTEGER :: output_unit
145 INTEGER,
INTENT(IN) :: ndep
147 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_initial_info'
149 INTEGER :: handle, homo, ispin, nao, &
152 CALL timeset(routinen, handle)
154 IF (output_unit > 0)
THEN
155 DO ispin = 1, dft_control%nspins
158 nelectron=nelectron_spin, &
161 IF (dft_control%nspins > 1)
THEN
162 WRITE (unit=output_unit, fmt=
"(/,T2,A,I2)")
"Spin", ispin
164 WRITE (unit=output_unit, fmt=
"(/,(T2,A,T71,I10))") &
165 "Number of electrons:", nelectron_spin, &
166 "Number of occupied orbitals:", homo, &
167 "Number of molecular orbitals:", nmo
169 WRITE (unit=output_unit, fmt=
"(/,(T2,A,T71,I10))") &
170 "Number of orbital functions:", nao, &
171 "Number of independent orbital functions:", nao - ndep
174 CALL timestop(handle)
189 LOGICAL,
INTENT(IN) :: final_mos
191 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_write_mos'
193 CHARACTER(LEN=2) :: solver_method
194 CHARACTER(LEN=3*default_string_length) :: message
195 CHARACTER(LEN=5) :: spin
196 CHARACTER(LEN=default_string_length), &
197 DIMENSION(:),
POINTER :: tmpstringlist
198 INTEGER :: handle, homo, ikp, ispin, iw, kpoint, &
199 nao, nelectron, nkp, nmo, nspin, numo
200 INTEGER,
DIMENSION(2) :: nmos_occ
201 INTEGER,
DIMENSION(:),
POINTER :: mo_index_range
202 LOGICAL :: do_kpoints, do_printout, print_eigvals, &
203 print_eigvecs, print_mo_info, &
204 print_occup, print_occup_stats
205 REAL(kind=
dp) :: flexible_electron_count, maxocc, n_el_f, &
206 occup_stats_occ_threshold
207 REAL(kind=
dp),
DIMENSION(:),
POINTER :: mo_eigenvalues, umo_eigenvalues
211 TYPE(
cp_fm_type),
POINTER :: mo_coeff, umo_coeff
214 TYPE(
dbcsr_type),
POINTER :: matrix_ks, matrix_s
222 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
226 CALL timeset(routinen, handle)
228 cpassert(
ASSOCIATED(qs_env))
232 blacs_env=blacs_env, &
233 dft_control=dft_control, &
234 do_kpoints=do_kpoints, &
236 qs_kind_set=qs_kind_set, &
238 particle_set=particle_set, &
239 scf_control=scf_control)
246 CALL section_vals_val_get(dft_section,
"PRINT%MO%OCCUPATION_NUMBERS_STATS", c_vals=tmpstringlist)
248 print_occup_stats = .false.
249 occup_stats_occ_threshold = 1e-6_dp
250 IF (
SIZE(tmpstringlist) > 0)
THEN
251 print_occup_stats = .true.
252 IF (len_trim(tmpstringlist(1)) > 0) &
253 READ (tmpstringlist(1), *) print_occup_stats
255 IF (
SIZE(tmpstringlist) > 1) &
256 READ (tmpstringlist(2), *) occup_stats_occ_threshold
261 IF ((.NOT. print_mo_info) .OR. (.NOT. (print_eigvals .OR. print_eigvecs .OR. print_occup .OR. print_occup_stats)))
THEN
262 CALL timestop(handle)
266 NULLIFY (fm_struct_tmp)
268 NULLIFY (mo_eigenvalues)
271 NULLIFY (umo_eigenvalues)
275 nspin = dft_control%nspins
281 nkp =
SIZE(kpoints%kp_env)
283 CALL get_qs_env(qs_env, matrix_ks=ks, matrix_s=s)
284 cpassert(
ASSOCIATED(ks))
285 cpassert(
ASSOCIATED(s))
289 kp_loop:
DO ikp = 1, nkp
292 mos => kpoints%kp_env(ikp)%kpoint_env%mos(1, :)
295 CALL get_qs_env(qs_env, matrix_ks=ks, mos=mos)
298 cpassert(
ASSOCIATED(mos))
309 cpabort(
"The OT method is not implemented for k points")
314 matrix_ks => ks(ispin)%matrix
315 matrix_s => s(1)%matrix
318 IF (dft_control%do_admm)
THEN
326 eigenvalues=mo_eigenvalues, &
329 nelectron=nelectron, &
333 flexible_electron_count=flexible_electron_count)
337 cpassert(
ASSOCIATED(mo_index_range))
338 IF (mo_index_range(2) == -1)
THEN
341 numo = min(mo_index_range(2) - homo, nao - homo)
360 flexible_electron_count=flexible_electron_count)
362 fm_struct=fm_struct_tmp, &
363 name=
"Temporary MO set (unoccupied MOs only) for printout")
366 mo_coeff=umo_coeff, &
367 eigenvalues=umo_eigenvalues)
373 NULLIFY (local_preconditioner)
374 IF (
ASSOCIATED(scf_env%ot_preconditioner))
THEN
375 local_preconditioner => scf_env%ot_preconditioner(1)%preconditioner
377 NULLIFY (local_preconditioner)
384 matrix_c_fm=umo_coeff, &
385 matrix_orthogonal_space_fm=mo_coeff, &
386 eps_gradient=scf_control%eps_lumos, &
388 iter_max=scf_control%max_iter_lumos, &
389 size_ortho_space=nmo)
392 ks_matrix=matrix_ks, &
393 evals_arg=umo_eigenvalues, &
398 IF (dft_control%do_admm)
THEN
406 message =
"The MO information is only calculated after SCF convergence "// &
407 "is achieved when the orbital transformation (OT) method is used"
408 cpwarn(trim(message))
409 do_printout = .false.
430 cpabort(
"Invalid spin")
432 IF (
ASSOCIATED(umo_set))
THEN
434 final_mos=final_mos, spin=trim(spin), solver_method=solver_method, &
438 final_mos=final_mos, spin=trim(spin), solver_method=solver_method)
441 IF (
ASSOCIATED(umo_set))
THEN
443 final_mos=final_mos, solver_method=solver_method, &
447 final_mos=final_mos, solver_method=solver_method)
451 nmos_occ(ispin) = max(nmos_occ(ispin), count(mo_set%occupation_numbers > occup_stats_occ_threshold))
455 IF (
ASSOCIATED(umo_set))
THEN
468 IF (do_printout .AND. print_mo_info .AND. print_occup_stats)
THEN
470 ignore_should_output=print_mo_info, &
473 IF (
SIZE(mos) > 1)
THEN
474 WRITE (unit=iw, fmt=
"(A,I4)")
" MO| Total occupied (ALPHA):", nmos_occ(1)
475 WRITE (unit=iw, fmt=
"(A,I4)")
" MO| Total occupied (BETA): ", nmos_occ(2)
477 WRITE (unit=iw, fmt=
"(A,I4)")
" MO| Total occupied: ", nmos_occ(1)
479 WRITE (unit=iw, fmt=
"(A)")
""
482 ignore_should_output=print_mo_info)
485 CALL timestop(handle)
500 energy, total_steps, should_stop, outer_loop_converged)
501 INTEGER :: output_unit
505 INTEGER :: total_steps
506 LOGICAL,
INTENT(IN) :: should_stop, outer_loop_converged
508 REAL(kind=
dp) :: outer_loop_eps
510 outer_loop_eps = sqrt(maxval(scf_env%outer_scf%gradient(:, scf_env%outer_scf%iter_count)**2))
511 IF (output_unit > 0)
WRITE (output_unit,
'(/,T3,A,I4,A,E10.2,A,F22.10)') &
512 "outer SCF iter = ", scf_env%outer_scf%iter_count, &
513 " RMS gradient = ", outer_loop_eps,
" energy =", energy%total
515 IF (outer_loop_converged)
THEN
516 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
517 "outer SCF loop converged in", scf_env%outer_scf%iter_count, &
518 " iterations or ", total_steps,
" steps"
519 ELSE IF (scf_env%outer_scf%iter_count > scf_control%outer_scf%max_scf &
520 .OR. should_stop)
THEN
521 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
522 "outer SCF loop FAILED to converge after ", &
523 scf_env%outer_scf%iter_count,
" iterations or ", total_steps,
" steps"
540 INTEGER :: output_unit
541 LOGICAL :: just_energy
542 REAL(kind=
dp) :: t1, t2
545 IF ((output_unit > 0) .AND. scf_env%print_iter_line)
THEN
546 IF (just_energy)
THEN
547 WRITE (unit=output_unit, &
548 fmt=
"(T2,A,1X,A,T20,E8.2,1X,F6.1,16X,F20.10)") &
549 " -", trim(scf_env%iter_method), scf_env%iter_param, t2 - t1, energy%total
551 IF ((abs(scf_env%iter_delta) < 1.0e-8_dp) .OR. &
552 (abs(scf_env%iter_delta) >= 1.0e5_dp))
THEN
553 WRITE (unit=output_unit, &
554 fmt=
"(T2,I5,1X,A,T20,E8.2,1X,F6.1,1X,ES14.4,1X,F20.10,1X,ES9.2)") &
555 scf_env%iter_count, trim(scf_env%iter_method), scf_env%iter_param, &
556 t2 - t1, scf_env%iter_delta, energy%total, energy%total - energy%tot_old
558 WRITE (unit=output_unit, &
559 fmt=
"(T2,I5,1X,A,T20,E8.2,1X,F6.1,1X,F14.8,1X,F20.10,1X,ES9.2)") &
560 scf_env%iter_count, trim(scf_env%iter_method), scf_env%iter_param, &
561 t2 - t1, scf_env%iter_delta, energy%total, energy%total - energy%tot_old
586 SUBROUTINE qs_scf_print_scf_summary(output_unit, rho, qs_charges, energy, nelectron_total, &
587 dft_control, qmmm, qs_env, gapw, gapw_xc)
588 INTEGER,
INTENT(IN) :: output_unit
592 INTEGER,
INTENT(IN) :: nelectron_total
594 LOGICAL,
INTENT(IN) :: qmmm
596 LOGICAL,
INTENT(IN) :: gapw, gapw_xc
598 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_print_scf_summary'
600 INTEGER :: bc, handle, ispin, psolver
601 REAL(kind=
dp) :: exc1_energy, exc_energy, &
602 implicit_ps_ehartree, tot1_h, tot1_s
603 REAL(kind=
dp),
DIMENSION(:),
POINTER :: tot_rho_r
606 NULLIFY (tot_rho_r, pw_env)
607 CALL timeset(routinen, handle)
610 psolver = pw_env%poisson_env%parameters%solver
612 IF (output_unit > 0)
THEN
614 IF (.NOT. (dft_control%qs_control%semi_empirical .OR. &
615 dft_control%qs_control%xtb .OR. &
616 dft_control%qs_control%dftb))
THEN
617 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T41,2F20.10))") &
618 "Electronic density on regular grids: ", &
621 "Core density on regular grids:", &
622 qs_charges%total_rho_core_rspace, &
623 qs_charges%total_rho_core_rspace + &
624 qs_charges%total_rho1_hard_nuc - &
625 REAL(nelectron_total + dft_control%charge,
dp)
627 IF (dft_control%correct_surf_dip)
THEN
628 WRITE (unit=output_unit, fmt=
"((T3,A,/,T3,A,T41,F20.10))") &
629 "Total dipole moment perpendicular to ", &
630 "the slab [electrons-Angstroem]: ", &
631 qs_env%surface_dipole_moment
635 tot1_h = qs_charges%total_rho1_hard(1)
636 tot1_s = qs_charges%total_rho1_soft(1)
637 DO ispin = 2, dft_control%nspins
638 tot1_h = tot1_h + qs_charges%total_rho1_hard(ispin)
639 tot1_s = tot1_s + qs_charges%total_rho1_soft(ispin)
641 WRITE (unit=output_unit, fmt=
"((T3,A,T41,2F20.10))") &
642 "Hard and soft densities (Lebedev):", &
644 WRITE (unit=output_unit, fmt=
"(T3,A,T41,F20.10)") &
645 "Total Rho_soft + Rho1_hard - Rho1_soft (r-space): ", &
647 "Total charge density (r-space): ", &
649 + qs_charges%total_rho_core_rspace &
650 + qs_charges%total_rho1_hard_nuc
651 IF (qs_charges%total_rho1_hard_nuc /= 0.0_dp)
THEN
652 WRITE (unit=output_unit, fmt=
"(T3,A,T41,F20.10)") &
653 "Total CNEO nuc. char. den. (Lebedev): ", &
654 qs_charges%total_rho1_hard_nuc, &
655 "Total CNEO soft char. den. (Lebedev): ", &
656 qs_charges%total_rho1_soft_nuc_lebedev, &
657 "Total CNEO soft char. den. (r-space): ", &
658 qs_charges%total_rho1_soft_nuc_rspace, &
659 "Total soft Rho_e+n+0 (g-space):", &
660 qs_charges%total_rho_gspace
662 WRITE (unit=output_unit, fmt=
"(T3,A,T41,F20.10)") &
663 "Total Rho_soft + Rho0_soft (g-space):", &
664 qs_charges%total_rho_gspace
668 WRITE (unit=output_unit, fmt=
"(T3,A,T41,F20.10)") &
669 "Total charge density on r-space grids: ", &
671 qs_charges%total_rho_core_rspace, &
672 "Total charge density g-space grids: ", &
673 qs_charges%total_rho_gspace
676 IF (dft_control%qs_control%semi_empirical)
THEN
677 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
678 "Core-core repulsion energy [eV]: ", energy%core_overlap*
evolt, &
679 "Core Hamiltonian energy [eV]: ", energy%core*
evolt, &
680 "Two-electron integral energy [eV]: ", energy%hartree*
evolt, &
681 "Electronic energy [eV]: ", &
682 (energy%core + 0.5_dp*energy%hartree)*
evolt
683 IF (energy%dispersion /= 0.0_dp) &
684 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
685 "Dispersion energy [eV]: ", energy%dispersion*
evolt
686 ELSEIF (dft_control%qs_control%dftb)
THEN
687 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
688 "Core Hamiltonian energy: ", energy%core, &
689 "Repulsive potential energy: ", energy%repulsive, &
690 "Electronic energy: ", energy%hartree, &
691 "Dispersion energy: ", energy%dispersion
692 IF (energy%dftb3 /= 0.0_dp) &
693 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
694 "DFTB3 3rd order energy: ", energy%dftb3
695 IF (energy%efield /= 0.0_dp) &
696 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
697 "Electric field interaction energy: ", energy%efield
698 ELSEIF (dft_control%qs_control%xtb)
THEN
699 IF (dft_control%qs_control%xtb_control%do_tblite)
THEN
700 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
701 "Core Hamiltonian energy: ", energy%core, &
702 "Repulsive potential energy: ", energy%repulsive, &
703 "Electrostatic energy: ", energy%el_stat, &
704 "Self-consistent dispersion energy: ", energy%dispersion_sc, &
705 "Non-self consistent dispersion energy: ", energy%dispersion, &
706 "Correction for halogen bonding: ", energy%xtb_xb_inter
708 IF (dft_control%qs_control%xtb_control%gfn_type == 0)
THEN
709 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
710 "Core Hamiltonian energy: ", energy%core, &
711 "Repulsive potential energy: ", energy%repulsive, &
712 "SRB Correction energy: ", energy%srb, &
713 "Charge equilibration energy: ", energy%eeq, &
714 "Dispersion energy: ", energy%dispersion
715 ELSEIF (dft_control%qs_control%xtb_control%gfn_type == 1)
THEN
716 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
717 "Core Hamiltonian energy: ", energy%core, &
718 "Repulsive potential energy: ", energy%repulsive, &
719 "Electronic energy: ", energy%hartree, &
720 "DFTB3 3rd order energy: ", energy%dftb3, &
721 "Dispersion energy: ", energy%dispersion
722 IF (dft_control%qs_control%xtb_control%xb_interaction) &
723 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
724 "Correction for halogen bonding: ", energy%xtb_xb_inter
725 ELSEIF (dft_control%qs_control%xtb_control%gfn_type == 2)
THEN
726 cpabort(
"gfn_typ 2 NYA")
728 cpabort(
"invalid gfn_typ")
731 IF (dft_control%qs_control%xtb_control%do_nonbonded) &
732 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
733 "Correction for nonbonded interactions: ", energy%xtb_nonbonded
734 IF (energy%efield /= 0.0_dp) &
735 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
736 "Electric field interaction energy: ", energy%efield
738 IF (dft_control%do_admm)
THEN
739 exc_energy = energy%exc + energy%exc_aux_fit
740 IF (gapw .OR. gapw_xc) exc1_energy = energy%exc1 + energy%exc1_aux_fit
742 exc_energy = energy%exc
743 IF (gapw .OR. gapw_xc) exc1_energy = energy%exc1
747 implicit_ps_ehartree = pw_env%poisson_env%implicit_env%ehartree
748 bc = pw_env%poisson_env%parameters%ps_implicit_params%boundary_condition
751 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
752 "Overlap energy of the core charge distribution:", energy%core_overlap, &
753 "Self energy of the core charge distribution: ", energy%core_self, &
754 "Core Hamiltonian energy: ", energy%core, &
755 "Hartree energy: ", implicit_ps_ehartree, &
756 "Electric enthalpy: ", energy%hartree, &
757 "Exchange-correlation energy: ", exc_energy
759 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
760 "Overlap energy of the core charge distribution:", energy%core_overlap, &
761 "Self energy of the core charge distribution: ", energy%core_self, &
762 "Core Hamiltonian energy: ", energy%core, &
763 "Hartree energy: ", energy%hartree, &
764 "Exchange-correlation energy: ", exc_energy
767 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
768 "Overlap energy of the core charge distribution:", energy%core_overlap, &
769 "Self energy of the core charge distribution: ", energy%core_self, &
770 "Core Hamiltonian energy: ", energy%core, &
771 "Hartree energy: ", energy%hartree, &
772 "Exchange-correlation energy: ", exc_energy
774 IF (energy%e_hartree /= 0.0_dp) &
775 WRITE (unit=output_unit, fmt=
"(T3,A,/,T3,A,T56,F25.14)") &
776 "Coulomb Electron-Electron Interaction Energy ", &
777 "- Already included in the total Hartree term ", energy%e_hartree
778 IF (energy%ex /= 0.0_dp) &
779 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
780 "Hartree-Fock Exchange energy: ", energy%ex
781 IF (energy%dispersion /= 0.0_dp) &
782 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
783 "Dispersion energy: ", energy%dispersion
784 IF (energy%gcp /= 0.0_dp) &
785 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
786 "gCP energy: ", energy%gcp
787 IF (energy%efield /= 0.0_dp) &
788 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
789 "Electric field interaction energy: ", energy%efield
791 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
792 "GAPW| Exc from hard and soft atomic rho1: ", exc1_energy, &
793 "GAPW| local Eh = 1 center integrals: ", energy%hartree_1c
796 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
797 "GAPW_XC| Exc from hard and soft atomic rho1: ", exc1_energy
799 IF (energy%core_cneo /= 0.0_dp)
THEN
800 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
801 "CNEO| quantum nuclear core energy: ", energy%core_cneo
804 IF (dft_control%hairy_probes .EQV. .true.)
THEN
805 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
806 "Electronic entropic energy:", energy%kTS
807 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
808 "Fermi energy:", energy%efermi
810 IF (dft_control%smear)
THEN
811 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
812 "Electronic entropic energy:", energy%kTS
813 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
814 "Fermi energy:", energy%efermi
816 IF (dft_control%dft_plus_u)
THEN
817 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
818 "DFT+U energy:", energy%dft_plus_u
820 IF (dft_control%do_sccs)
THEN
821 WRITE (unit=output_unit, fmt=
"(A)")
""
825 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
826 "QM/MM Electrostatic energy: ", energy%qmmm_el
827 IF (qs_env%qmmm_env_qm%image_charge)
THEN
828 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
829 "QM/MM image charge energy: ", energy%image_charge
832 IF (dft_control%qs_control%mulliken_restraint)
THEN
833 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
834 "Mulliken restraint energy: ", energy%mulliken
836 IF (dft_control%qs_control%semi_empirical)
THEN
837 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
838 "Total energy [eV]: ", energy%total*
evolt
839 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
840 "Atomic reference energy [eV]: ", energy%core_self*
evolt, &
841 "Heat of formation [kcal/mol]: ", &
842 (energy%total + energy%core_self)*
kcalmol
844 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
845 "Total energy: ", energy%total
848 IF (qs_env%qmmm_env_qm%image_charge)
THEN
855 CALL timestop(handle)
857 END SUBROUTINE qs_scf_print_scf_summary
872 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_loop_print'
874 INTEGER :: after, handle, ic, ispin, iw
875 LOGICAL :: do_kpoints, omit_headers
876 REAL(kind=
dp) :: mo_mag_max, mo_mag_min, orthonormality
878 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_ks, matrix_p, matrix_s
885 CALL timeset(routinen, handle)
887 CALL get_qs_env(qs_env=qs_env, input=input, dft_control=dft_control, &
888 do_kpoints=do_kpoints)
894 DO ispin = 1, dft_control%nspins
897 dft_section,
"PRINT%AO_MATRICES/DENSITY"),
cp_p_file))
THEN
903 after = min(max(after, 1), 16)
904 DO ic = 1,
SIZE(matrix_p, 2)
906 output_unit=iw, omit_headers=omit_headers)
909 "PRINT%AO_MATRICES/DENSITY")
913 dft_section,
"PRINT%AO_MATRICES/KOHN_SHAM_MATRIX"),
cp_p_file))
THEN
917 after = min(max(after, 1), 16)
918 CALL get_qs_env(qs_env=qs_env, matrix_ks_kp=matrix_ks)
919 DO ic = 1,
SIZE(matrix_ks, 2)
920 IF (dft_control%qs_control%semi_empirical)
THEN
922 scale=
evolt, output_unit=iw, omit_headers=omit_headers)
925 output_unit=iw, omit_headers=omit_headers)
929 "PRINT%AO_MATRICES/KOHN_SHAM_MATRIX")
935 scf_section,
"PRINT%MO_ORTHONORMALITY"),
cp_p_file))
THEN
940 WRITE (iw,
'(T8,A)') &
941 " K-points: Maximum deviation from MO S-orthonormality not determined"
944 "PRINT%MO_ORTHONORMALITY")
950 CALL get_qs_env(qs_env=qs_env, matrix_s_kp=matrix_s)
956 WRITE (iw,
'(T8,A,T61,E20.4)') &
957 " Maximum deviation from MO S-orthonormality", orthonormality
960 "PRINT%MO_ORTHONORMALITY")
964 scf_section,
"PRINT%MO_MAGNITUDE"),
cp_p_file))
THEN
969 WRITE (iw,
'(T8,A)') &
970 " K-points: Minimum/Maximum MO magnitude not determined"
973 "PRINT%MO_MAGNITUDE")
980 WRITE (iw,
'(T8,A,T41,2E20.4)') &
981 " Minimum/Maximum MO magnitude ", mo_mag_min, mo_mag_max
984 "PRINT%MO_MAGNITUDE")
988 CALL timestop(handle)
1007 energy, total_steps, should_stop, outer_loop_converged, &
1009 INTEGER :: output_unit
1014 INTEGER :: total_steps
1015 LOGICAL,
INTENT(IN) :: should_stop, outer_loop_converged, &
1018 REAL(kind=
dp) :: outer_loop_eps
1021 outer_loop_eps = sqrt(maxval(scf_env%outer_scf%gradient(:, scf_env%outer_scf%iter_count)**2))
1022 IF (output_unit > 0)
WRITE (output_unit,
'(/,T3,A,I4,A,E10.2,A,F22.10)') &
1023 "CDFT SCF iter = ", scf_env%outer_scf%iter_count, &
1024 " RMS gradient = ", outer_loop_eps,
" energy =", energy%total
1025 IF (outer_loop_converged)
THEN
1026 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
1027 "CDFT SCF loop converged in", scf_env%outer_scf%iter_count, &
1028 " iterations or ", total_steps,
" steps"
1030 IF ((scf_env%outer_scf%iter_count > scf_control%outer_scf%max_scf .OR. should_stop) &
1031 .AND. .NOT. outer_loop_converged)
THEN
1032 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
1033 "CDFT SCF loop FAILED to converge after ", &
1034 scf_env%outer_scf%iter_count,
" iterations or ", total_steps,
" steps"
1049 INTEGER :: output_unit
1052 IF (output_unit > 0)
THEN
1053 WRITE (output_unit,
'(/,A)') &
1054 " ---------------------------------- CDFT --------------------------------------"
1055 WRITE (output_unit,
'(A)') &
1056 " Optimizing a density constraint in an external SCF loop "
1057 WRITE (output_unit,
'(A)')
" "
1058 SELECT CASE (cdft_control%type)
1060 WRITE (output_unit,
'(A)')
" Type of constraint: Hirshfeld"
1062 WRITE (output_unit,
'(A)')
" Type of constraint: Becke"
1064 WRITE (output_unit,
'(A,I8)')
" Number of constraints: ",
SIZE(cdft_control%group)
1065 WRITE (output_unit,
'(A,L8)')
" Using fragment densities:", cdft_control%fragment_density
1066 WRITE (output_unit,
'(A)')
" "
1067 IF (cdft_control%atomic_charges)
WRITE (output_unit,
'(A,/)')
" Calculating atomic CDFT charges"
1068 SELECT CASE (cdft_control%constraint_control%optimizer)
1070 WRITE (output_unit,
'(A)') &
1071 " Minimizer : SD : steepest descent"
1073 WRITE (output_unit,
'(A)') &
1074 " Minimizer : DIIS : direct inversion"
1075 WRITE (output_unit,
'(A)') &
1076 " in the iterative subspace"
1077 WRITE (output_unit,
'(A,I3,A)') &
1079 cdft_control%constraint_control%diis_buffer_length,
" DIIS vectors"
1081 WRITE (output_unit,
'(A)') &
1082 " Minimizer : BISECT : gradient bisection"
1083 WRITE (output_unit,
'(A,I3)') &
1084 " using a trust count of", &
1085 cdft_control%constraint_control%bisect_trust_count
1089 cdft_control%constraint_control%optimizer, output_unit)
1091 WRITE (output_unit,
'(A)')
" Minimizer : Secant"
1095 WRITE (output_unit,
'(/,A,L7)') &
1096 " Reusing OT preconditioner: ", cdft_control%reuse_precond
1097 IF (cdft_control%reuse_precond)
THEN
1098 WRITE (output_unit,
'(A,I3,A,I3,A)') &
1099 " using old preconditioner for up to ", &
1100 cdft_control%max_reuse,
" subsequent CDFT SCF"
1101 WRITE (output_unit,
'(A,I3,A,I3,A)') &
1102 " iterations if the relevant loop converged in less than ", &
1103 cdft_control%precond_freq,
" steps"
1105 SELECT CASE (cdft_control%type)
1107 WRITE (output_unit,
'(/,A)')
" Hirshfeld constraint settings"
1108 WRITE (output_unit,
'(A)')
" "
1109 SELECT CASE (cdft_control%hirshfeld_control%shape_function)
1111 WRITE (output_unit,
'(A, A8)') &
1112 " Shape function type: ",
"Gaussian"
1113 WRITE (output_unit,
'(A)', advance=
'NO') &
1114 " Type of Gaussian: "
1115 SELECT CASE (cdft_control%hirshfeld_control%gaussian_shape)
1117 WRITE (output_unit,
'(A13)')
"Default"
1119 WRITE (output_unit,
'(A13)')
"Covalent"
1121 WRITE (output_unit,
'(A13)')
"Fixed radius"
1123 WRITE (output_unit,
'(A13)')
"Van der Waals"
1125 WRITE (output_unit,
'(A13)')
"User-defined"
1129 WRITE (output_unit,
'(A, A8)') &
1130 " Shape function type: ",
"Density"
1133 WRITE (output_unit,
'(/, A)')
" Becke constraint settings"
1134 WRITE (output_unit,
'(A)')
" "
1135 SELECT CASE (cdft_control%becke_control%cutoff_type)
1137 WRITE (output_unit,
'(A,F8.3,A)') &
1138 " Cutoff for partitioning :",
cp_unit_from_cp2k(cdft_control%becke_control%rglobal, &
1139 "angstrom"),
" angstrom"
1141 WRITE (output_unit,
'(A)') &
1142 " Using element specific cutoffs for partitioning"
1144 WRITE (output_unit,
'(A,L7)') &
1145 " Skipping distant gpoints: ", cdft_control%becke_control%should_skip
1146 WRITE (output_unit,
'(A,L7)') &
1147 " Precompute gradients : ", cdft_control%becke_control%in_memory
1148 WRITE (output_unit,
'(A)')
" "
1149 IF (cdft_control%becke_control%adjust) &
1150 WRITE (output_unit,
'(A)') &
1151 " Using atomic radii to generate a heteronuclear charge partitioning"
1152 WRITE (output_unit,
'(A)')
" "
1153 IF (.NOT. cdft_control%becke_control%cavity_confine)
THEN
1154 WRITE (output_unit,
'(A)') &
1155 " No confinement is active"
1157 WRITE (output_unit,
'(A)')
" Confinement using a Gaussian shaped cavity is active"
1158 SELECT CASE (cdft_control%becke_control%cavity_shape)
1160 WRITE (output_unit,
'(A,F8.4, A)') &
1161 " Type of Gaussian : Fixed radius: ", &
1164 WRITE (output_unit,
'(A)') &
1165 " Type of Gaussian : Covalent radius "
1167 WRITE (output_unit,
'(A)') &
1168 " Type of Gaussian : vdW radius "
1170 WRITE (output_unit,
'(A)') &
1171 " Type of Gaussian : User radius "
1173 WRITE (output_unit,
'(A,ES12.4)') &
1174 " Cavity threshold : ", cdft_control%becke_control%eps_cavity
1177 WRITE (output_unit,
'(/,A)') &
1178 " ---------------------------------- CDFT --------------------------------------"
1191 INTEGER :: output_unit
1196 IF (output_unit > 0)
THEN
1197 SELECT CASE (cdft_control%type)
1199 WRITE (output_unit,
'(/,T3,A,T60)') &
1200 '------------------- Hirshfeld constraint information -------------------'
1202 WRITE (output_unit,
'(/,T3,A,T60)') &
1203 '--------------------- Becke constraint information ---------------------'
1205 cpabort(
"Unknown CDFT constraint.")
1207 DO igroup = 1,
SIZE(cdft_control%target)
1208 IF (igroup > 1)
WRITE (output_unit,
'(T3,A)')
' '
1209 WRITE (output_unit,
'(T3,A,T54,(3X,I18))') &
1210 'Atomic group :', igroup
1211 SELECT CASE (cdft_control%group(igroup)%constraint_type)
1213 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1214 WRITE (output_unit,
'(T3,A,T42,A)') &
1215 'Type of constraint :', adjustr(
'Charge density constraint (frag.)')
1217 WRITE (output_unit,
'(T3,A,T50,A)') &
1218 'Type of constraint :', adjustr(
'Charge density constraint')
1221 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1222 WRITE (output_unit,
'(T3,A,T35,A)') &
1223 'Type of constraint :', adjustr(
'Magnetization density constraint (frag.)')
1225 WRITE (output_unit,
'(T3,A,T43,A)') &
1226 'Type of constraint :', adjustr(
'Magnetization density constraint')
1229 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1230 WRITE (output_unit,
'(T3,A,T38,A)') &
1231 'Type of constraint :', adjustr(
'Alpha spin density constraint (frag.)')
1233 WRITE (output_unit,
'(T3,A,T46,A)') &
1234 'Type of constraint :', adjustr(
'Alpha spin density constraint')
1237 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1238 WRITE (output_unit,
'(T3,A,T39,A)') &
1239 'Type of constraint :', adjustr(
'Beta spin density constraint (frag.)')
1241 WRITE (output_unit,
'(T3,A,T47,A)') &
1242 'Type of constraint :', adjustr(
'Beta spin density constraint')
1245 cpabort(
"Unknown constraint type.")
1247 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1248 'Target value of constraint :', cdft_control%target(igroup)
1249 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1250 'Current value of constraint :', cdft_control%value(igroup)
1251 WRITE (output_unit,
'(T3,A,T59,(3X,ES13.3))') &
1252 'Deviation from target :', cdft_control%value(igroup) - cdft_control%target(igroup)
1253 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1254 'Strength of constraint :', cdft_control%strength(igroup)
1256 WRITE (output_unit,
'(T3,A)') &
1257 '------------------------------------------------------------------------'
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)
...
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_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, rhoz_cneo_set, ecoul_1c, rho0_s_rs, rho0_s_gs, rhoz_cneo_s_rs, rhoz_cneo_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs, do_rixs, tb_tblite)
Get the QUICKSTEP environment.
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, 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_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_initial_info(output_unit, mos, dft_control, ndep)
writes basic information at the beginning of an scf run
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
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.