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
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 IF (dft_control%qs_control%xtb_control%gfn_type == 0)
THEN
674 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
675 "Core Hamiltonian energy: ", energy%core, &
676 "Repulsive potential energy: ", energy%repulsive, &
677 "SRB Correction energy: ", energy%srb, &
678 "Charge equilibration energy: ", energy%eeq, &
679 "Dispersion energy: ", energy%dispersion
680 ELSEIF (dft_control%qs_control%xtb_control%gfn_type == 1)
THEN
681 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
682 "Core Hamiltonian energy: ", energy%core, &
683 "Repulsive potential energy: ", energy%repulsive, &
684 "Electronic energy: ", energy%hartree, &
685 "DFTB3 3rd order energy: ", energy%dftb3, &
686 "Dispersion energy: ", energy%dispersion
687 IF (dft_control%qs_control%xtb_control%xb_interaction) &
688 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
689 "Correction for halogen bonding: ", energy%xtb_xb_inter
690 ELSEIF (dft_control%qs_control%xtb_control%gfn_type == 2)
THEN
691 cpabort(
"gfn_typ 2 NYA")
693 cpabort(
"invalid gfn_typ")
695 IF (dft_control%qs_control%xtb_control%do_nonbonded) &
696 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
697 "Correction for nonbonded interactions: ", energy%xtb_nonbonded
698 IF (energy%efield /= 0.0_dp) &
699 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
700 "Electric field interaction energy: ", energy%efield
702 IF (dft_control%do_admm)
THEN
703 exc_energy = energy%exc + energy%exc_aux_fit
704 IF (gapw .OR. gapw_xc) exc1_energy = energy%exc1 + energy%exc1_aux_fit
706 exc_energy = energy%exc
707 IF (gapw .OR. gapw_xc) exc1_energy = energy%exc1
711 implicit_ps_ehartree = pw_env%poisson_env%implicit_env%ehartree
712 bc = pw_env%poisson_env%parameters%ps_implicit_params%boundary_condition
715 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
716 "Overlap energy of the core charge distribution:", energy%core_overlap, &
717 "Self energy of the core charge distribution: ", energy%core_self, &
718 "Core Hamiltonian energy: ", energy%core, &
719 "Hartree energy: ", implicit_ps_ehartree, &
720 "Electric enthalpy: ", energy%hartree, &
721 "Exchange-correlation energy: ", exc_energy
723 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
724 "Overlap energy of the core charge distribution:", energy%core_overlap, &
725 "Self energy of the core charge distribution: ", energy%core_self, &
726 "Core Hamiltonian energy: ", energy%core, &
727 "Hartree energy: ", energy%hartree, &
728 "Exchange-correlation energy: ", exc_energy
731 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
732 "Overlap energy of the core charge distribution:", energy%core_overlap, &
733 "Self energy of the core charge distribution: ", energy%core_self, &
734 "Core Hamiltonian energy: ", energy%core, &
735 "Hartree energy: ", energy%hartree, &
736 "Exchange-correlation energy: ", exc_energy
738 IF (energy%e_hartree /= 0.0_dp) &
739 WRITE (unit=output_unit, fmt=
"(T3,A,/,T3,A,T56,F25.14)") &
740 "Coulomb Electron-Electron Interaction Energy ", &
741 "- Already included in the total Hartree term ", energy%e_hartree
742 IF (energy%ex /= 0.0_dp) &
743 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
744 "Hartree-Fock Exchange energy: ", energy%ex
745 IF (energy%dispersion /= 0.0_dp) &
746 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
747 "Dispersion energy: ", energy%dispersion
748 IF (energy%gcp /= 0.0_dp) &
749 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
750 "gCP energy: ", energy%gcp
751 IF (energy%efield /= 0.0_dp) &
752 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
753 "Electric field interaction energy: ", energy%efield
755 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
756 "GAPW| Exc from hard and soft atomic rho1: ", exc1_energy, &
757 "GAPW| local Eh = 1 center integrals: ", energy%hartree_1c
760 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
761 "GAPW_XC| Exc from hard and soft atomic rho1: ", exc1_energy
764 IF (dft_control%smear)
THEN
765 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
766 "Electronic entropic energy:", energy%kTS
767 WRITE (unit=output_unit, fmt=
"((T3,A,T56,F25.14))") &
768 "Fermi energy:", energy%efermi
770 IF (dft_control%dft_plus_u)
THEN
771 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
772 "DFT+U energy:", energy%dft_plus_u
774 IF (dft_control%do_sccs)
THEN
775 WRITE (unit=output_unit, fmt=
"(A)")
""
779 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
780 "QM/MM Electrostatic energy: ", energy%qmmm_el
781 IF (qs_env%qmmm_env_qm%image_charge)
THEN
782 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
783 "QM/MM image charge energy: ", energy%image_charge
786 IF (dft_control%qs_control%mulliken_restraint)
THEN
787 WRITE (unit=output_unit, fmt=
"(T3,A,T56,F25.14)") &
788 "Mulliken restraint energy: ", energy%mulliken
790 IF (dft_control%qs_control%semi_empirical)
THEN
791 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
792 "Total energy [eV]: ", energy%total*
evolt
793 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
794 "Atomic reference energy [eV]: ", energy%core_self*
evolt, &
795 "Heat of formation [kcal/mol]: ", &
796 (energy%total + energy%core_self)*
kcalmol
798 WRITE (unit=output_unit, fmt=
"(/,(T3,A,T56,F25.14))") &
799 "Total energy: ", energy%total
802 IF (qs_env%qmmm_env_qm%image_charge)
THEN
809 CALL timestop(handle)
811 END SUBROUTINE qs_scf_print_scf_summary
826 CHARACTER(LEN=*),
PARAMETER :: routinen =
'qs_scf_loop_print'
828 INTEGER :: after, handle, ic, ispin, iw
829 LOGICAL :: do_kpoints, omit_headers
830 REAL(kind=
dp) :: mo_mag_max, mo_mag_min, orthonormality
832 TYPE(
dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_ks, matrix_p, matrix_s
839 CALL timeset(routinen, handle)
841 CALL get_qs_env(qs_env=qs_env, input=input, dft_control=dft_control, &
842 do_kpoints=do_kpoints)
848 DO ispin = 1, dft_control%nspins
851 dft_section,
"PRINT%AO_MATRICES/DENSITY"),
cp_p_file))
THEN
857 after = min(max(after, 1), 16)
858 DO ic = 1,
SIZE(matrix_p, 2)
860 output_unit=iw, omit_headers=omit_headers)
863 "PRINT%AO_MATRICES/DENSITY")
867 dft_section,
"PRINT%AO_MATRICES/KOHN_SHAM_MATRIX"),
cp_p_file))
THEN
871 after = min(max(after, 1), 16)
872 CALL get_qs_env(qs_env=qs_env, matrix_ks_kp=matrix_ks)
873 DO ic = 1,
SIZE(matrix_ks, 2)
874 IF (dft_control%qs_control%semi_empirical)
THEN
876 scale=
evolt, output_unit=iw, omit_headers=omit_headers)
879 output_unit=iw, omit_headers=omit_headers)
883 "PRINT%AO_MATRICES/KOHN_SHAM_MATRIX")
889 scf_section,
"PRINT%MO_ORTHONORMALITY"),
cp_p_file))
THEN
894 WRITE (iw,
'(T8,A)') &
895 " K-points: Maximum deviation from MO S-orthonormality not determined"
898 "PRINT%MO_ORTHONORMALITY")
904 CALL get_qs_env(qs_env=qs_env, matrix_s_kp=matrix_s)
910 WRITE (iw,
'(T8,A,T61,E20.4)') &
911 " Maximum deviation from MO S-orthonormality", orthonormality
914 "PRINT%MO_ORTHONORMALITY")
918 scf_section,
"PRINT%MO_MAGNITUDE"),
cp_p_file))
THEN
923 WRITE (iw,
'(T8,A)') &
924 " K-points: Minimum/Maximum MO magnitude not determined"
927 "PRINT%MO_MAGNITUDE")
934 WRITE (iw,
'(T8,A,T41,2E20.4)') &
935 " Minimum/Maximum MO magnitude ", mo_mag_min, mo_mag_max
938 "PRINT%MO_MAGNITUDE")
942 CALL timestop(handle)
961 energy, total_steps, should_stop, outer_loop_converged, &
963 INTEGER :: output_unit
968 INTEGER :: total_steps
969 LOGICAL,
INTENT(IN) :: should_stop, outer_loop_converged, &
972 REAL(kind=
dp) :: outer_loop_eps
975 outer_loop_eps = sqrt(maxval(scf_env%outer_scf%gradient(:, scf_env%outer_scf%iter_count)**2))
976 IF (output_unit > 0)
WRITE (output_unit,
'(/,T3,A,I4,A,E10.2,A,F22.10)') &
977 "CDFT SCF iter = ", scf_env%outer_scf%iter_count, &
978 " RMS gradient = ", outer_loop_eps,
" energy =", energy%total
979 IF (outer_loop_converged)
THEN
980 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
981 "CDFT SCF loop converged in", scf_env%outer_scf%iter_count, &
982 " iterations or ", total_steps,
" steps"
984 IF ((scf_env%outer_scf%iter_count > scf_control%outer_scf%max_scf .OR. should_stop) &
985 .AND. .NOT. outer_loop_converged)
THEN
986 IF (output_unit > 0)
WRITE (output_unit,
'(T3,A,I4,A,I4,A,/)') &
987 "CDFT SCF loop FAILED to converge after ", &
988 scf_env%outer_scf%iter_count,
" iterations or ", total_steps,
" steps"
1003 INTEGER :: output_unit
1006 IF (output_unit > 0)
THEN
1007 WRITE (output_unit,
'(/,A)') &
1008 " ---------------------------------- CDFT --------------------------------------"
1009 WRITE (output_unit,
'(A)') &
1010 " Optimizing a density constraint in an external SCF loop "
1011 WRITE (output_unit,
'(A)')
" "
1012 SELECT CASE (cdft_control%type)
1014 WRITE (output_unit,
'(A)')
" Type of constraint: Hirshfeld"
1016 WRITE (output_unit,
'(A)')
" Type of constraint: Becke"
1018 WRITE (output_unit,
'(A,I8)')
" Number of constraints: ",
SIZE(cdft_control%group)
1019 WRITE (output_unit,
'(A,L8)')
" Using fragment densities:", cdft_control%fragment_density
1020 WRITE (output_unit,
'(A)')
" "
1021 IF (cdft_control%atomic_charges)
WRITE (output_unit,
'(A,/)')
" Calculating atomic CDFT charges"
1022 SELECT CASE (cdft_control%constraint_control%optimizer)
1024 WRITE (output_unit,
'(A)') &
1025 " Minimizer : SD : steepest descent"
1027 WRITE (output_unit,
'(A)') &
1028 " Minimizer : DIIS : direct inversion"
1029 WRITE (output_unit,
'(A)') &
1030 " in the iterative subspace"
1031 WRITE (output_unit,
'(A,I3,A)') &
1033 cdft_control%constraint_control%diis_buffer_length,
" DIIS vectors"
1035 WRITE (output_unit,
'(A)') &
1036 " Minimizer : BISECT : gradient bisection"
1037 WRITE (output_unit,
'(A,I3)') &
1038 " using a trust count of", &
1039 cdft_control%constraint_control%bisect_trust_count
1043 cdft_control%constraint_control%optimizer, output_unit)
1045 WRITE (output_unit,
'(A)')
" Minimizer : Secant"
1049 WRITE (output_unit,
'(/,A,L7)') &
1050 " Reusing OT preconditioner: ", cdft_control%reuse_precond
1051 IF (cdft_control%reuse_precond)
THEN
1052 WRITE (output_unit,
'(A,I3,A,I3,A)') &
1053 " using old preconditioner for up to ", &
1054 cdft_control%max_reuse,
" subsequent CDFT SCF"
1055 WRITE (output_unit,
'(A,I3,A,I3,A)') &
1056 " iterations if the relevant loop converged in less than ", &
1057 cdft_control%precond_freq,
" steps"
1059 SELECT CASE (cdft_control%type)
1061 WRITE (output_unit,
'(/,A)')
" Hirshfeld constraint settings"
1062 WRITE (output_unit,
'(A)')
" "
1063 SELECT CASE (cdft_control%hirshfeld_control%shape_function)
1065 WRITE (output_unit,
'(A, A8)') &
1066 " Shape function type: ",
"Gaussian"
1067 WRITE (output_unit,
'(A)', advance=
'NO') &
1068 " Type of Gaussian: "
1069 SELECT CASE (cdft_control%hirshfeld_control%gaussian_shape)
1071 WRITE (output_unit,
'(A13)')
"Default"
1073 WRITE (output_unit,
'(A13)')
"Covalent"
1075 WRITE (output_unit,
'(A13)')
"Fixed radius"
1077 WRITE (output_unit,
'(A13)')
"Van der Waals"
1079 WRITE (output_unit,
'(A13)')
"User-defined"
1083 WRITE (output_unit,
'(A, A8)') &
1084 " Shape function type: ",
"Density"
1087 WRITE (output_unit,
'(/, A)')
" Becke constraint settings"
1088 WRITE (output_unit,
'(A)')
" "
1089 SELECT CASE (cdft_control%becke_control%cutoff_type)
1091 WRITE (output_unit,
'(A,F8.3,A)') &
1092 " Cutoff for partitioning :",
cp_unit_from_cp2k(cdft_control%becke_control%rglobal, &
1093 "angstrom"),
" angstrom"
1095 WRITE (output_unit,
'(A)') &
1096 " Using element specific cutoffs for partitioning"
1098 WRITE (output_unit,
'(A,L7)') &
1099 " Skipping distant gpoints: ", cdft_control%becke_control%should_skip
1100 WRITE (output_unit,
'(A,L7)') &
1101 " Precompute gradients : ", cdft_control%becke_control%in_memory
1102 WRITE (output_unit,
'(A)')
" "
1103 IF (cdft_control%becke_control%adjust) &
1104 WRITE (output_unit,
'(A)') &
1105 " Using atomic radii to generate a heteronuclear charge partitioning"
1106 WRITE (output_unit,
'(A)')
" "
1107 IF (.NOT. cdft_control%becke_control%cavity_confine)
THEN
1108 WRITE (output_unit,
'(A)') &
1109 " No confinement is active"
1111 WRITE (output_unit,
'(A)')
" Confinement using a Gaussian shaped cavity is active"
1112 SELECT CASE (cdft_control%becke_control%cavity_shape)
1114 WRITE (output_unit,
'(A,F8.4, A)') &
1115 " Type of Gaussian : Fixed radius: ", &
1118 WRITE (output_unit,
'(A)') &
1119 " Type of Gaussian : Covalent radius "
1121 WRITE (output_unit,
'(A)') &
1122 " Type of Gaussian : vdW radius "
1124 WRITE (output_unit,
'(A)') &
1125 " Type of Gaussian : User radius "
1127 WRITE (output_unit,
'(A,ES12.4)') &
1128 " Cavity threshold : ", cdft_control%becke_control%eps_cavity
1131 WRITE (output_unit,
'(/,A)') &
1132 " ---------------------------------- CDFT --------------------------------------"
1145 INTEGER :: output_unit
1150 IF (output_unit > 0)
THEN
1151 SELECT CASE (cdft_control%type)
1153 WRITE (output_unit,
'(/,T3,A,T60)') &
1154 '------------------- Hirshfeld constraint information -------------------'
1156 WRITE (output_unit,
'(/,T3,A,T60)') &
1157 '--------------------- Becke constraint information ---------------------'
1159 cpabort(
"Unknown CDFT constraint.")
1161 DO igroup = 1,
SIZE(cdft_control%target)
1162 IF (igroup > 1)
WRITE (output_unit,
'(T3,A)')
' '
1163 WRITE (output_unit,
'(T3,A,T54,(3X,I18))') &
1164 'Atomic group :', igroup
1165 SELECT CASE (cdft_control%group(igroup)%constraint_type)
1167 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1168 WRITE (output_unit,
'(T3,A,T42,A)') &
1169 'Type of constraint :', adjustr(
'Charge density constraint (frag.)')
1171 WRITE (output_unit,
'(T3,A,T50,A)') &
1172 'Type of constraint :', adjustr(
'Charge density constraint')
1175 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1176 WRITE (output_unit,
'(T3,A,T35,A)') &
1177 'Type of constraint :', adjustr(
'Magnetization density constraint (frag.)')
1179 WRITE (output_unit,
'(T3,A,T43,A)') &
1180 'Type of constraint :', adjustr(
'Magnetization density constraint')
1183 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1184 WRITE (output_unit,
'(T3,A,T38,A)') &
1185 'Type of constraint :', adjustr(
'Alpha spin density constraint (frag.)')
1187 WRITE (output_unit,
'(T3,A,T46,A)') &
1188 'Type of constraint :', adjustr(
'Alpha spin density constraint')
1191 IF (cdft_control%group(igroup)%is_fragment_constraint)
THEN
1192 WRITE (output_unit,
'(T3,A,T39,A)') &
1193 'Type of constraint :', adjustr(
'Beta spin density constraint (frag.)')
1195 WRITE (output_unit,
'(T3,A,T47,A)') &
1196 'Type of constraint :', adjustr(
'Beta spin density constraint')
1199 cpabort(
"Unknown constraint type.")
1201 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1202 'Target value of constraint :', cdft_control%target(igroup)
1203 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1204 'Current value of constraint :', cdft_control%value(igroup)
1205 WRITE (output_unit,
'(T3,A,T59,(3X,ES13.3))') &
1206 'Deviation from target :', cdft_control%value(igroup) - cdft_control%target(igroup)
1207 WRITE (output_unit,
'(T3,A,T54,(3X,F18.12))') &
1208 'Strength of constraint :', cdft_control%strength(igroup)
1210 WRITE (output_unit,
'(T3,A)') &
1211 '------------------------------------------------------------------------'
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_pp, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, harris_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, wanniercentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, eeq, rhs)
Get the QUICKSTEP environment.
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.