79#include "./base/base_uses.f90"
85 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_mo_io'
108 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(IN) :: mo_array
111 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
113 CHARACTER(LEN=*),
PARAMETER :: routinen =
'write_mo_set_to_restart'
114 CHARACTER(LEN=30),
DIMENSION(2),
PARAMETER :: &
115 keys = (/
"SCF%PRINT%RESTART_HISTORY",
"SCF%PRINT%RESTART "/)
117 INTEGER :: handle, ikey, ires, ispin
120 CALL timeset(routinen, handle)
128 IF (mo_array(1)%use_mo_coeff_b)
THEN
131 DO ispin = 1,
SIZE(mo_array)
132 IF (.NOT.
ASSOCIATED(mo_array(ispin)%mo_coeff_b))
THEN
136 mo_array(ispin)%mo_coeff)
140 DO ikey = 1,
SIZE(keys)
143 dft_section, keys(ikey)),
cp_p_file))
THEN
146 extension=
".wfn", file_status=
"REPLACE", file_action=
"WRITE", &
147 do_backup=.true., file_form=
"UNFORMATTED")
150 qs_kind_set=qs_kind_set, ires=ires)
157 CALL timestop(handle)
171 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(IN) :: mo_array
173 TYPE(
dbcsr_p_type),
DIMENSION(:),
POINTER :: tmpl_matrix
175 CHARACTER(LEN=*),
PARAMETER :: routinen =
'write_dm_binary_restart'
177 CHARACTER(LEN=default_path_length) :: file_name, project_name
178 INTEGER :: handle, ispin, unit_nr
179 LOGICAL :: do_dm_restart
180 REAL(kind=
dp) :: cs_pos
184 CALL timeset(routinen, handle)
186 IF (logger%para_env%is_source())
THEN
192 project_name = logger%iter_info%project_name
194 NULLIFY (matrix_p_tmp)
196 IF (do_dm_restart)
THEN
197 ALLOCATE (matrix_p_tmp)
198 DO ispin = 1,
SIZE(mo_array)
199 CALL dbcsr_create(matrix_p_tmp, template=tmpl_matrix(ispin)%matrix, name=
"DM RESTART")
201 IF (.NOT.
ASSOCIATED(mo_array(ispin)%mo_coeff_b)) cpabort(
"mo_coeff_b NOT ASSOCIATED")
205 use_dbcsr=.true., retain_sparsity=.false.)
207 WRITE (file_name,
'(A,I0,A)') trim(project_name)//
"_SCF_DM_SPIN_", ispin,
"_RESTART.dm"
209 IF (unit_nr > 0)
THEN
210 WRITE (unit_nr,
'(T2,A,E20.8)')
"Writing restart DM "//trim(file_name)//
" with checksum: ", cs_pos
216 DEALLOCATE (matrix_p_tmp)
219 CALL timestop(handle)
234 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(IN) :: mo_array
235 TYPE(
cp_fm_type),
DIMENSION(:),
INTENT(IN) :: rt_mos
238 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
240 CHARACTER(LEN=*),
PARAMETER :: routinen =
'write_rt_mos_to_restart'
241 CHARACTER(LEN=43),
DIMENSION(2),
PARAMETER :: keys = (/ &
242 "REAL_TIME_PROPAGATION%PRINT%RESTART_HISTORY", &
243 "REAL_TIME_PROPAGATION%PRINT%RESTART "/)
245 INTEGER :: handle, ikey, ires
248 CALL timeset(routinen, handle)
256 DO ikey = 1,
SIZE(keys)
259 dft_section, keys(ikey)),
cp_p_file))
THEN
262 extension=
".rtpwfn", file_status=
"REPLACE", file_action=
"WRITE", &
263 do_backup=.true., file_form=
"UNFORMATTED")
266 particle_set=particle_set, ires=ires)
272 CALL timestop(handle)
286 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(IN) :: mo_array
287 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
293 CHARACTER(LEN=*),
PARAMETER :: routinen =
'write_mo_set_low'
295 INTEGER :: handle, iatom, ikind, imat, iset, &
296 ishell, ispin, lmax, lshell, &
297 max_block, nao, natom, nmo, nset, &
298 nset_max, nshell_max, nspin
299 INTEGER,
DIMENSION(:),
POINTER :: nset_info, nshell
300 INTEGER,
DIMENSION(:, :),
POINTER :: l, nshell_info
301 INTEGER,
DIMENSION(:, :, :),
POINTER :: nso_info
305 CALL timeset(routinen, handle)
306 nspin =
SIZE(mo_array)
307 nao = mo_array(1)%nao
311 natom =
SIZE(particle_set, 1)
316 NULLIFY (orb_basis_set, dftb_parameter)
317 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
319 basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
320 IF (
ASSOCIATED(orb_basis_set))
THEN
325 nset_max = max(nset_max, nset)
327 nshell_max = max(nshell_max, nshell(iset))
329 ELSEIF (
ASSOCIATED(dftb_parameter))
THEN
331 nset_max = max(nset_max, 1)
332 nshell_max = max(nshell_max, lmax + 1)
339 ALLOCATE (nso_info(nshell_max, nset_max, natom))
340 nso_info(:, :, :) = 0
342 ALLOCATE (nshell_info(nset_max, natom))
343 nshell_info(:, :) = 0
345 ALLOCATE (nset_info(natom))
349 NULLIFY (orb_basis_set, dftb_parameter)
350 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
352 basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
353 IF (
ASSOCIATED(orb_basis_set))
THEN
358 nset_info(iatom) = nset
360 nshell_info(iset, iatom) = nshell(iset)
361 DO ishell = 1, nshell(iset)
362 lshell = l(ishell, iset)
363 nso_info(ishell, iset, iatom) =
nso(lshell)
366 ELSEIF (
ASSOCIATED(dftb_parameter))
THEN
369 nshell_info(1, iatom) = lmax + 1
370 DO ishell = 1, lmax + 1
372 nso_info(ishell, 1, iatom) =
nso(lshell)
380 WRITE (ires) natom, nspin, nao, nset_max, nshell_max
381 WRITE (ires) nset_info
382 WRITE (ires) nshell_info
383 WRITE (ires) nso_info
385 DEALLOCATE (nset_info)
387 DEALLOCATE (nshell_info)
389 DEALLOCATE (nso_info)
395 nmo = mo_array(ispin)%nmo
396 IF ((ires > 0) .AND. (nmo > 0))
THEN
398 mo_array(ispin)%homo, &
399 mo_array(ispin)%lfomo, &
400 mo_array(ispin)%nelectron
401 WRITE (ires) mo_array(ispin)%eigenvalues(1:nmo), &
402 mo_array(ispin)%occupation_numbers(1:nmo)
404 IF (
PRESENT(rt_mos))
THEN
405 DO imat = 2*ispin - 1, 2*ispin
413 CALL timestop(handle)
428 CHARACTER(LEN=default_path_length),
INTENT(OUT) :: filename
429 LOGICAL,
INTENT(OUT) :: exist
432 LOGICAL,
INTENT(IN),
OPTIONAL :: kp, xas, rtp
435 LOGICAL :: my_kp, my_rtp, my_xas
441 IF (
PRESENT(kp)) my_kp = kp
442 IF (
PRESENT(xas)) my_xas = xas
443 IF (
PRESENT(rtp)) my_rtp = rtp
447 IF (n_rep_val > 0)
THEN
454 extension=
"", my_local=.false.)
455 ELSE IF (my_rtp)
THEN
459 extension=
".rtpwfn", my_local=.false.)
464 extension=
".kp", my_local=.false.)
469 extension=
".wfn", my_local=.false.)
472 IF (.NOT. my_xas)
THEN
473 INQUIRE (file=filename, exist=exist)
493 para_env, id_nr, multiplicity, dft_section, natom_mismatch, &
496 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(INOUT) :: mo_array
498 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
501 INTEGER,
INTENT(IN) :: id_nr, multiplicity
503 LOGICAL,
INTENT(OUT),
OPTIONAL :: natom_mismatch
504 LOGICAL,
INTENT(IN),
OPTIONAL :: cdft
505 INTEGER,
INTENT(IN),
OPTIONAL :: out_unit
507 CHARACTER(LEN=*),
PARAMETER :: routinen =
'read_mo_set_from_restart'
509 CHARACTER(LEN=default_path_length) :: file_name
510 INTEGER :: handle, ispin, my_out_unit, natom, &
512 LOGICAL :: exist, my_cdft
515 CALL timeset(routinen, handle)
518 IF (
PRESENT(cdft)) my_cdft = cdft
520 IF (
PRESENT(out_unit)) my_out_unit = out_unit
522 nspin =
SIZE(mo_array)
525 IF (para_env%is_source())
THEN
527 natom =
SIZE(particle_set, 1)
531 file_name = trim(file_name)//
".bak-"//adjustl(
cp_to_string(id_nr))
535 file_action=
"READ", &
536 file_form=
"UNFORMATTED", &
538 unit_number=restart_unit)
543 particle_set=particle_set, natom=natom, &
544 rst_unit=restart_unit, multiplicity=multiplicity, natom_mismatch=natom_mismatch)
546 IF (
PRESENT(natom_mismatch))
THEN
548 CALL para_env%bcast(natom_mismatch)
549 IF (natom_mismatch)
THEN
550 IF (para_env%is_source())
CALL close_file(unit_number=restart_unit)
551 CALL timestop(handle)
557 IF (para_env%is_source())
THEN
558 IF (my_out_unit > 0)
THEN
559 WRITE (unit=my_out_unit, fmt=
"(T2,A)") &
560 "WFN_RESTART| Restart file "//trim(file_name)//
" read"
566 IF (.NOT. my_cdft)
THEN
569 particle_set, dft_section, 4, 0, final_mos=.false.)
573 CALL timestop(handle)
590 particle_set, para_env, id_nr, multiplicity, dft_section)
592 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(INOUT) :: mo_array
593 TYPE(
cp_fm_type),
DIMENSION(:),
POINTER :: rt_mos
595 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
598 INTEGER,
INTENT(IN) :: id_nr, multiplicity
601 CHARACTER(LEN=*),
PARAMETER :: routinen =
'read_rt_mos_from_restart'
603 CHARACTER(LEN=default_path_length) :: file_name
604 INTEGER :: handle, ispin, natom, nspin, &
605 restart_unit, unit_nr
609 CALL timeset(routinen, handle)
612 nspin =
SIZE(mo_array)
615 IF (para_env%is_source())
THEN
617 natom =
SIZE(particle_set, 1)
621 file_name = trim(file_name)//
".bak-"//adjustl(
cp_to_string(id_nr))
625 IF (unit_nr > 0)
THEN
626 WRITE (unit_nr,
'(T2,A)')
"Read RTP restart from the file: "//trim(file_name)
630 file_action=
"READ", &
631 file_form=
"UNFORMATTED", &
633 unit_number=restart_unit)
638 particle_set=particle_set, qs_kind_set=qs_kind_set, natom=natom, &
639 rst_unit=restart_unit, multiplicity=multiplicity)
642 IF (para_env%is_source())
CALL close_file(unit_number=restart_unit)
646 particle_set, dft_section, 4, 0, final_mos=.false.)
649 CALL timestop(handle)
669 multiplicity, rt_mos, natom_mismatch)
671 TYPE(
mo_set_type),
DIMENSION(:),
INTENT(INOUT) :: mos
673 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
675 INTEGER,
INTENT(IN) :: natom, rst_unit
676 INTEGER,
INTENT(in),
OPTIONAL :: multiplicity
677 TYPE(
cp_fm_type),
DIMENSION(:),
OPTIONAL,
POINTER :: rt_mos
678 LOGICAL,
INTENT(OUT),
OPTIONAL :: natom_mismatch
680 INTEGER :: homo, homo_read, i, iatom, ikind, imat, irow, iset, iset_read, ishell, &
681 ishell_read, iso, ispin, lfomo_read, lmax, lshell, my_mult, nao, nao_read, natom_read, &
682 nelectron, nelectron_read, nmo, nmo_read, nnshell, nset, nset_max, nshell_max, nspin, &
683 nspin_read, offset_read
684 INTEGER,
DIMENSION(:),
POINTER :: nset_info, nshell
685 INTEGER,
DIMENSION(:, :),
POINTER :: l, nshell_info
686 INTEGER,
DIMENSION(:, :, :),
POINTER :: nso_info, offset_info
687 LOGICAL :: minbas, natom_match, use_this
688 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: eig_read, occ_read
689 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: vecbuffer, vecbuffer_read
699 IF (
PRESENT(multiplicity)) my_mult = multiplicity
701 IF (para_env%is_source())
THEN
702 READ (rst_unit) natom_read, nspin_read, nao_read, nset_max, nshell_max
703 IF (
PRESENT(rt_mos))
THEN
704 IF (nspin_read /= nspin)
THEN
705 cpabort(
"To change nspin is not possible. ")
709 IF (nspin_read /= nspin)
THEN
711 "READ RESTART : WARNING : nspin is not equal "
714 IF (nspin_read > nspin)
THEN
715 cpabort(
"Reducing nspin is not possible. ")
719 natom_match = (natom_read == natom)
721 IF (natom_match)
THEN
724 ALLOCATE (nso_info(nshell_max, nset_max, natom_read))
725 ALLOCATE (nshell_info(nset_max, natom_read))
726 ALLOCATE (nset_info(natom_read))
727 ALLOCATE (offset_info(nshell_max, nset_max, natom_read))
729 IF (nao_read /= nao)
THEN
731 " READ RESTART : WARNING : DIFFERENT # AOs ", nao, nao_read
732 IF (
PRESENT(rt_mos)) &
733 cpabort(
"To change basis is not possible. ")
736 READ (rst_unit) nset_info
737 READ (rst_unit) nshell_info
738 READ (rst_unit) nso_info
742 DO iset = 1, nset_info(iatom)
743 DO ishell = 1, nshell_info(iset, iatom)
744 offset_info(ishell, iset, iatom) = i
745 i = i + nso_info(ishell, iset, iatom)
750 ALLOCATE (vecbuffer_read(1, nao_read))
756 CALL para_env%bcast(natom_match)
757 IF (
PRESENT(natom_mismatch)) natom_mismatch = .NOT. natom_match
759 IF (.NOT. natom_match)
THEN
760 IF (
PRESENT(natom_mismatch))
THEN
762 " READ RESTART : WARNING : DIFFERENT natom, returning ", natom, natom_read
765 cpabort(
"Incorrect number of atoms in restart file. ")
769 CALL para_env%bcast(nspin_read)
771 ALLOCATE (vecbuffer(1, nao))
776 homo = mos(ispin)%homo
777 mos(ispin)%eigenvalues(:) = 0.0_dp
778 mos(ispin)%occupation_numbers(:) = 0.0_dp
781 IF (para_env%is_source() .AND. (nmo > 0))
THEN
782 READ (rst_unit) nmo_read, homo_read, lfomo_read, nelectron_read
783 ALLOCATE (eig_read(nmo_read), occ_read(nmo_read))
787 nmo = min(nmo, nmo_read)
788 IF (nmo_read < nmo) &
789 CALL cp_warn(__location__, &
790 "The number of MOs on the restart unit is smaller than the number of "// &
791 "the allocated MOs. The MO set will be padded with zeros!")
792 IF (nmo_read > nmo) &
793 CALL cp_warn(__location__, &
794 "The number of MOs on the restart unit is greater than the number of "// &
795 "the allocated MOs. The read MO set will be truncated!")
797 READ (rst_unit) eig_read(1:nmo_read), occ_read(1:nmo_read)
798 mos(ispin)%eigenvalues(1:nmo) = eig_read(1:nmo)
799 mos(ispin)%occupation_numbers(1:nmo) = occ_read(1:nmo)
800 DEALLOCATE (eig_read, occ_read)
802 mos(ispin)%homo = homo_read
803 mos(ispin)%lfomo = lfomo_read
804 IF (homo_read > nmo)
THEN
805 IF (nelectron_read == mos(ispin)%nelectron)
THEN
806 CALL cp_warn(__location__, &
807 "The number of occupied MOs on the restart unit is larger than "// &
808 "the allocated MOs. The read MO set will be truncated and the occupation numbers recalculated!")
813 cpabort(
"Number of occupied MOs on restart unit larger than allocated MOs. ")
818 CALL para_env%bcast(nmo)
819 CALL para_env%bcast(mos(ispin)%homo)
820 CALL para_env%bcast(mos(ispin)%lfomo)
821 CALL para_env%bcast(mos(ispin)%nelectron)
822 CALL para_env%bcast(mos(ispin)%eigenvalues)
823 CALL para_env%bcast(mos(ispin)%occupation_numbers)
825 IF (
PRESENT(rt_mos))
THEN
826 DO imat = 2*ispin - 1, 2*ispin
828 IF (para_env%is_source())
THEN
829 READ (rst_unit) vecbuffer
831 vecbuffer(1, :) = 0.0_dp
833 CALL para_env%bcast(vecbuffer)
835 vecbuffer, 1, i, nao, 1, transpose=.true.)
840 IF (para_env%is_source())
THEN
841 READ (rst_unit) vecbuffer_read
846 NULLIFY (orb_basis_set, dftb_parameter, l, nshell)
847 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
849 basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
850 IF (
ASSOCIATED(orb_basis_set))
THEN
856 ELSEIF (
ASSOCIATED(dftb_parameter))
THEN
873 nnshell = nshell(iset)
875 DO ishell = 1, nnshell
879 lshell = l(ishell, iset)
881 IF (iset_read > nset_info(iatom)) use_this = .false.
883 IF (
nso(lshell) == nso_info(ishell_read, iset_read, iatom))
THEN
884 offset_read = offset_info(ishell_read, iset_read, iatom)
885 ishell_read = ishell_read + 1
886 IF (ishell_read > nshell_info(iset, iatom))
THEN
888 iset_read = iset_read + 1
894 DO iso = 1,
nso(lshell)
896 IF (offset_read - 1 + iso .LT. 1 .OR. offset_read - 1 + iso .GT. nao_read)
THEN
897 vecbuffer(1, irow) = 0.0_dp
899 vecbuffer(1, irow) = vecbuffer_read(1, offset_read - 1 + iso)
902 vecbuffer(1, irow) = 0.0_dp
913 vecbuffer(1, :) = 0.0_dp
917 CALL para_env%bcast(vecbuffer)
919 vecbuffer, 1, i, nao, 1, transpose=.true.)
923 IF (para_env%is_source())
THEN
926 DO i = nmo + 1, nmo_read
927 READ (rst_unit) vecbuffer_read
932 IF (.NOT.
PRESENT(rt_mos))
THEN
933 IF (ispin == 1 .AND. nspin_read < nspin)
THEN
935 mos(ispin + 1)%homo = mos(ispin)%homo
936 mos(ispin + 1)%lfomo = mos(ispin)%lfomo
937 nelectron = mos(ispin)%nelectron
938 IF (my_mult .NE. 1)
THEN
939 CALL cp_abort(__location__, &
940 "Restarting an LSD calculation from an LDA wfn only works for multiplicity=1 (singlets).")
942 IF (mos(ispin + 1)%nelectron < 0)
THEN
943 cpabort(
"LSD: too few electrons for this multiplisity. ")
945 mos(ispin + 1)%eigenvalues = mos(ispin)%eigenvalues
946 mos(ispin)%occupation_numbers = mos(ispin)%occupation_numbers/2.0_dp
947 mos(ispin + 1)%occupation_numbers = mos(ispin)%occupation_numbers
948 CALL cp_fm_to_fm(mos(ispin)%mo_coeff, mos(ispin + 1)%mo_coeff)
954 DEALLOCATE (vecbuffer)
956 IF (para_env%is_source())
THEN
957 DEALLOCATE (vecbuffer_read)
958 DEALLOCATE (offset_info)
959 DEALLOCATE (nso_info)
960 DEALLOCATE (nshell_info)
961 DEALLOCATE (nset_info)
993 dft_section, before, kpoint, final_mos, spin, &
994 solver_method, rtp, cpart, sim_step, umo_set)
998 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
1001 INTEGER,
INTENT(IN) :: before, kpoint
1002 LOGICAL,
INTENT(IN),
OPTIONAL :: final_mos
1003 CHARACTER(LEN=*),
INTENT(IN),
OPTIONAL :: spin
1004 CHARACTER(LEN=2),
INTENT(IN),
OPTIONAL :: solver_method
1005 LOGICAL,
INTENT(IN),
OPTIONAL :: rtp
1006 INTEGER,
INTENT(IN),
OPTIONAL :: cpart, sim_step
1007 TYPE(
mo_set_type),
INTENT(IN),
OPTIONAL :: umo_set
1009 CHARACTER(LEN=12) :: symbol
1010 CHARACTER(LEN=12),
DIMENSION(:),
POINTER :: bcgf_symbol
1011 CHARACTER(LEN=14) :: fmtstr5
1012 CHARACTER(LEN=15) :: energy_str, orbital_str, step_string
1013 CHARACTER(LEN=2) :: element_symbol, my_solver_method
1014 CHARACTER(LEN=2*default_string_length) :: name
1015 CHARACTER(LEN=21) :: vector_str
1016 CHARACTER(LEN=22) :: fmtstr4
1017 CHARACTER(LEN=24) :: fmtstr2
1018 CHARACTER(LEN=25) :: fmtstr1
1019 CHARACTER(LEN=29) :: fmtstr6
1020 CHARACTER(LEN=4) :: reim
1021 CHARACTER(LEN=40) :: fmtstr3
1022 CHARACTER(LEN=6),
DIMENSION(:),
POINTER :: bsgf_symbol
1023 INTEGER :: after, first_mo, from, homo, iatom, icgf, ico, icol, ikind, imo, irow, iset, &
1024 isgf, ishell, iso, iw, jcol, last_mo, left, lmax, lshell, nao, natom, ncgf, ncol, nmo, &
1025 nset, nsgf, numo, right, scf_step, to, width
1026 INTEGER,
DIMENSION(:),
POINTER :: mo_index_range, nshell
1027 INTEGER,
DIMENSION(:, :),
POINTER :: l
1028 LOGICAL :: ionode, my_final, my_rtp, &
1029 print_cartesian, print_eigvals, &
1030 print_eigvecs, print_occup, &
1032 REAL(kind=
dp) :: gap, maxocc
1033 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:) :: mo_eigenvalues, mo_occupation_numbers
1034 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: cmatrix, smatrix
1035 REAL(kind=
dp),
DIMENSION(:),
POINTER :: eigenvalues, occupation_numbers
1036 TYPE(
cp_fm_type),
POINTER :: mo_coeff, umo_coeff
1041 NULLIFY (bcgf_symbol)
1042 NULLIFY (bsgf_symbol)
1044 NULLIFY (mo_index_range)
1049 ionode = logger%para_env%is_source()
1056 after = min(max(after, 1), 16)
1059 IF (
PRESENT(final_mos))
THEN
1060 my_final = final_mos
1067 IF (
PRESENT(rtp))
THEN
1072 .OR. (sim_step == 1)
1078 IF ((.NOT. should_output) .OR. (.NOT. (print_eigvals .OR. print_eigvecs .OR. print_occup)))
RETURN
1081 cpassert(
PRESENT(sim_step))
1082 cpassert(
PRESENT(cpart))
1084 IF (cpart == 0)
THEN
1089 print_eigvals = .false.
1091 scf_step = max(0, logger%iter_info%iteration(logger%iter_info%n_rlevel) - 1)
1094 IF (.NOT. my_final)
THEN
1095 IF (.NOT. my_rtp)
THEN
1096 step_string =
" AFTER SCF STEP"
1098 step_string =
" AFTER RTP STEP"
1102 IF (
PRESENT(solver_method))
THEN
1103 my_solver_method = solver_method
1106 my_solver_method =
"TD"
1111 mo_coeff=mo_coeff, &
1112 eigenvalues=eigenvalues, &
1113 occupation_numbers=occupation_numbers, &
1118 IF (
PRESENT(umo_set))
THEN
1120 mo_coeff=umo_coeff, &
1126 ALLOCATE (mo_eigenvalues(nmo))
1127 mo_eigenvalues(:) = 0.0_dp
1128 mo_eigenvalues(1:nmo - numo) = eigenvalues(1:nmo - numo)
1129 ALLOCATE (mo_occupation_numbers(nmo))
1130 mo_occupation_numbers(:) = 0.0_dp
1131 mo_occupation_numbers(1:nmo - numo) = occupation_numbers(1:nmo - numo)
1134 eigenvalues=eigenvalues)
1135 mo_eigenvalues(nmo - numo + 1:nmo) = eigenvalues(1:numo)
1138 IF (print_eigvecs)
THEN
1139 ALLOCATE (smatrix(nao, nmo))
1144 IF (.NOT. ionode)
THEN
1145 DEALLOCATE (smatrix)
1150 ignore_should_output=should_output, &
1160 fmtstr1 =
"(T2,A,21X, ( X,I5, X))"
1161 fmtstr2 =
"(T2,A,21X, (1X,F . ))"
1162 fmtstr3 =
"(T2,A,I5,1X,I5,1X,A,1X,A6, (1X,F . ))"
1164 width = before + after + 3
1165 ncol = int(56/width)
1167 right = max((after - 2), 1)
1168 left = width - right - 5
1170 WRITE (unit=fmtstr1(11:12), fmt=
"(I2)") ncol
1171 WRITE (unit=fmtstr1(14:15), fmt=
"(I2)") left
1172 WRITE (unit=fmtstr1(21:22), fmt=
"(I2)") right
1174 WRITE (unit=fmtstr2(11:12), fmt=
"(I2)") ncol
1175 WRITE (unit=fmtstr2(18:19), fmt=
"(I2)") width - 1
1176 WRITE (unit=fmtstr2(21:22), fmt=
"(I2)") after
1178 WRITE (unit=fmtstr3(27:28), fmt=
"(I2)") ncol
1179 WRITE (unit=fmtstr3(34:35), fmt=
"(I2)") width - 1
1180 WRITE (unit=fmtstr3(37:38), fmt=
"(I2)") after
1182 IF (my_final .OR. (my_solver_method ==
"TD"))
THEN
1183 energy_str =
"EIGENVALUES"
1184 vector_str =
"EIGENVECTORS"
1186 energy_str =
"ENERGIES"
1187 vector_str =
"COEFFICIENTS"
1191 energy_str =
"ZEROS"
1192 vector_str = trim(reim)//
" RTP COEFFICIENTS"
1195 IF (print_eigvecs)
THEN
1197 IF (print_cartesian)
THEN
1199 orbital_str =
"CARTESIAN"
1201 ALLOCATE (cmatrix(ncgf, ncgf))
1208 NULLIFY (orb_basis_set, dftb_parameter)
1209 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
1211 basis_set=orb_basis_set, &
1212 dftb_parameter=dftb_parameter)
1213 IF (
ASSOCIATED(orb_basis_set))
THEN
1219 DO ishell = 1, nshell(iset)
1220 lshell = l(ishell, iset)
1221 CALL dgemm(
"T",
"N",
nco(lshell), nmo,
nso(lshell), 1.0_dp, &
1223 smatrix(isgf, 1), nsgf, 0.0_dp, &
1224 cmatrix(icgf, 1), ncgf)
1225 icgf = icgf +
nco(lshell)
1226 isgf = isgf +
nso(lshell)
1229 ELSE IF (
ASSOCIATED(dftb_parameter))
THEN
1231 DO ishell = 1, lmax + 1
1233 CALL dgemm(
"T",
"N",
nco(lshell), nsgf,
nso(lshell), 1.0_dp, &
1235 smatrix(isgf, 1), nsgf, 0.0_dp, &
1236 cmatrix(icgf, 1), ncgf)
1237 icgf = icgf +
nco(lshell)
1238 isgf = isgf +
nso(lshell)
1248 orbital_str =
"SPHERICAL"
1252 name = trim(energy_str)//
", OCCUPATION NUMBERS, AND "// &
1253 trim(orbital_str)//
" "//trim(vector_str)
1255 IF (.NOT. my_final) &
1256 WRITE (unit=name, fmt=
"(A,1X,I0)") trim(name)//step_string, scf_step
1258 ELSE IF (print_occup .OR. print_eigvals)
THEN
1259 name = trim(energy_str)//
" AND OCCUPATION NUMBERS"
1261 IF (.NOT. my_final) &
1262 WRITE (unit=name, fmt=
"(A,1X,I0)") trim(name)//step_string, scf_step
1266 IF (
PRESENT(spin) .AND. (kpoint > 0))
THEN
1267 WRITE (unit=iw, fmt=
"(/,T2,A,I0)") &
1268 "MO| "//trim(spin)//
" "//trim(name)//
" FOR K POINT ", kpoint
1269 ELSE IF (
PRESENT(spin))
THEN
1270 WRITE (unit=iw, fmt=
"(/,T2,A)") &
1271 "MO| "//trim(spin)//
" "//trim(name)
1272 ELSE IF (kpoint > 0)
THEN
1273 WRITE (unit=iw, fmt=
"(/,T2,A,I0)") &
1274 "MO| "//trim(name)//
" FOR K POINT ", kpoint
1276 WRITE (unit=iw, fmt=
"(/,T2,A)") &
1281 IF ((mo_index_range(1) > 0) .AND. &
1282 (mo_index_range(2) > 0) .AND. &
1283 (mo_index_range(2) >= mo_index_range(1)))
THEN
1284 first_mo = max(1, mo_index_range(1))
1285 last_mo = min(nmo, mo_index_range(2))
1291 IF (print_eigvecs)
THEN
1295 DO icol = first_mo, last_mo, ncol
1298 to = min((from + ncol - 1), last_mo)
1300 WRITE (unit=iw, fmt=
"(T2,A)")
"MO|"
1301 WRITE (unit=iw, fmt=fmtstr1) &
1302 "MO|", (jcol, jcol=from, to)
1303 WRITE (unit=iw, fmt=fmtstr2) &
1304 "MO|", (mo_eigenvalues(jcol), jcol=from, to)
1305 WRITE (unit=iw, fmt=
"(T2,A)")
"MO|"
1306 WRITE (unit=iw, fmt=fmtstr2) &
1307 "MO|", (mo_occupation_numbers(jcol), jcol=from, to)
1308 WRITE (unit=iw, fmt=
"(T2,A)")
"MO|"
1314 IF (iatom /= 1)
WRITE (unit=iw, fmt=
"(T2,A)")
"MO|"
1316 NULLIFY (orb_basis_set, dftb_parameter)
1318 element_symbol=element_symbol, kind_number=ikind)
1320 basis_set=orb_basis_set, &
1321 dftb_parameter=dftb_parameter)
1323 IF (print_cartesian)
THEN
1325 IF (
ASSOCIATED(orb_basis_set))
THEN
1334 DO ishell = 1, nshell(iset)
1335 lshell = l(ishell, iset)
1336 DO ico = 1,
nco(lshell)
1337 WRITE (unit=iw, fmt=fmtstr3) &
1338 "MO|", irow, iatom, adjustr(element_symbol), bcgf_symbol(icgf), &
1339 (cmatrix(irow, jcol), jcol=from, to)
1345 ELSE IF (
ASSOCIATED(dftb_parameter))
THEN
1348 DO ishell = 1, lmax + 1
1350 DO ico = 1,
nco(lshell)
1353 WRITE (unit=iw, fmt=fmtstr3) &
1354 "MO|", irow, iatom, adjustr(element_symbol), symbol, &
1355 (cmatrix(irow, jcol), jcol=from, to)
1367 IF (
ASSOCIATED(orb_basis_set))
THEN
1375 DO ishell = 1, nshell(iset)
1376 lshell = l(ishell, iset)
1377 DO iso = 1,
nso(lshell)
1378 WRITE (unit=iw, fmt=fmtstr3) &
1379 "MO|", irow, iatom, adjustr(element_symbol), bsgf_symbol(isgf), &
1380 (smatrix(irow, jcol), jcol=from, to)
1386 ELSE IF (
ASSOCIATED(dftb_parameter))
THEN
1389 DO ishell = 1, lmax + 1
1391 DO iso = 1,
nso(lshell)
1392 symbol =
sgf_symbol(1, lshell, -lshell + iso - 1)
1394 WRITE (unit=iw, fmt=fmtstr3) &
1395 "MO|", irow, iatom, adjustr(element_symbol), symbol, &
1396 (smatrix(irow, jcol), jcol=from, to)
1412 WRITE (unit=iw, fmt=
"(T2,A)")
"MO|"
1416 DEALLOCATE (smatrix)
1417 IF (print_cartesian)
THEN
1418 DEALLOCATE (cmatrix)
1421 ELSE IF (print_occup .OR. print_eigvals)
THEN
1423 WRITE (unit=iw, fmt=
"(T2,A)")
"MO|"
1424 fmtstr4 =
"(T2,A,I7,3(1X,F22. ))"
1425 WRITE (unit=fmtstr4(19:20), fmt=
"(I2)") after
1426 IF (my_final .OR. (my_solver_method ==
"TD"))
THEN
1427 WRITE (unit=iw, fmt=
"(A)") &
1428 " MO| Index Eigenvalue [a.u.] Eigenvalue [eV] Occupation"
1430 WRITE (unit=iw, fmt=
"(A)") &
1431 " MO| Index Energy [a.u.] Energy [eV] Occupation"
1433 DO imo = first_mo, last_mo
1434 WRITE (unit=iw, fmt=fmtstr4) &
1435 "MO|", imo, mo_eigenvalues(imo), &
1436 mo_eigenvalues(imo)*
evolt, &
1437 mo_occupation_numbers(imo)
1439 fmtstr5 =
"(A,T59,F22. )"
1440 WRITE (unit=fmtstr5(12:13), fmt=
"(I2)") after
1441 WRITE (unit=iw, fmt=fmtstr5) &
1446 IF (.NOT. my_rtp)
THEN
1447 fmtstr6 =
"(A,T18,F17. ,A,T41,F17. ,A)"
1448 WRITE (unit=fmtstr6(12:13), fmt=
"(I2)") after
1449 WRITE (unit=fmtstr6(25:26), fmt=
"(I2)") after
1450 WRITE (unit=iw, fmt=fmtstr6) &
1451 " MO| E(Fermi):", mo_set%mu,
" a.u.", mo_set%mu*
evolt,
" eV"
1453 IF ((homo > 0) .AND. .NOT. my_rtp)
THEN
1454 IF ((mo_occupation_numbers(homo) == maxocc) .AND. (last_mo > homo))
THEN
1455 gap = mo_eigenvalues(homo + 1) - &
1456 mo_eigenvalues(homo)
1457 WRITE (unit=iw, fmt=fmtstr6) &
1458 " MO| Band gap:", gap,
" a.u.", gap*
evolt,
" eV"
1461 WRITE (unit=iw, fmt=
"(A)")
""
1465 IF (
ALLOCATED(mo_eigenvalues))
DEALLOCATE (mo_eigenvalues)
1466 IF (
ALLOCATED(mo_occupation_numbers))
DEALLOCATE (mo_occupation_numbers)
1469 ignore_should_output=should_output)
static void dgemm(const char transa, const char transb, const int m, const int n, const int k, const double alpha, const double *a, const int lda, const double *b, const int ldb, const double beta, double *c, const int ldc)
Convenient wrapper to hide Fortran nature of dgemm_, swapping a and b.
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
subroutine, public get_gto_basis_set(gto_basis_set, name, aliases, norm_type, kind_radius, ncgf, nset, nsgf, cgf_symbol, sgf_symbol, norm_cgf, set_radius, lmax, lmin, lx, ly, lz, m, ncgf_set, npgf, nsgf_set, nshell, cphi, pgf_radius, sphi, scon, zet, first_cgf, first_sgf, l, last_cgf, last_sgf, n, gcc, maxco, maxl, maxpgf, maxsgf_set, maxshell, maxso, nco_sum, npgf_sum, nshell_sum, maxder, short_kind_radius, npgf_seg_sum)
...
subroutine, public dbcsr_binary_write(matrix, filepath)
...
subroutine, public dbcsr_release(matrix)
...
real(kind=dp) function, public dbcsr_checksum(matrix, pos)
Calculates the checksum of a DBCSR matrix.
DBCSR operations in CP2K.
subroutine, public copy_dbcsr_to_fm(matrix, fm)
Copy a DBCSR matrix to a BLACS matrix.
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
Utility routines to open and close files. Tracking of preconnections.
subroutine, public open_file(file_name, file_status, file_form, file_action, file_position, file_pad, unit_number, debug, skip_get_unit_number, file_access)
Opens the requested file using a free unit number.
subroutine, public close_file(unit_number, file_status, keep_preconnection)
Close an open file given by its logical unit number. Optionally, keep the file and unit preconnected.
represent a full matrix distributed on many processors
subroutine, public cp_fm_get_info(matrix, name, nrow_global, ncol_global, nrow_block, ncol_block, nrow_local, ncol_local, row_indices, col_indices, local_data, context, nrow_locals, ncol_locals, matrix_struct, para_env)
returns all kind of information about the full matrix
subroutine, public cp_fm_write_unformatted(fm, unit)
...
subroutine, public cp_fm_set_submatrix(fm, new_values, start_row, start_col, n_rows, n_cols, alpha, beta, transpose)
sets a submatrix of a full matrix fm(start_row:start_row+n_rows,start_col:start_col+n_cols) = alpha*o...
subroutine, public cp_fm_set_all(matrix, alpha, beta)
set all elements of a matrix to the same value, and optionally the diagonal to a different one
subroutine, public cp_fm_get_submatrix(fm, target_m, start_row, start_col, n_rows, n_cols, transpose)
gets a submatrix of a full matrix op(target_m)(1:n_rows,1:n_cols) =fm(start_row:start_row+n_rows,...
various routines to log and control the output. The idea is that decisions about where to log should ...
recursive integer function, public cp_logger_get_default_unit_nr(logger, local, skip_not_ionode)
asks the default unit number of the given logger. try to use cp_logger_get_unit_nr
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)
...
character(len=default_path_length) function, public cp_print_key_generate_filename(logger, print_key, middle_name, extension, my_local)
Utility function that returns a unit number to write the print key. Might open a file with a unique f...
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...
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
integer, parameter, public default_path_length
Interface to the message passing library MPI.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public nco
integer, dimension(:, :), allocatable, public indco
integer, dimension(:), allocatable, public nso
character(len=12) function, public cgf_symbol(n, lxyz)
Build a Cartesian orbital symbol (orbital labels for printing).
character(len=6) function, public sgf_symbol(n, l, m)
Build a spherical orbital symbol (orbital labels for printing).
Define the data structure for the particle information.
Definition of physical constants:
real(kind=dp), parameter, public evolt
collects routines that calculate density matrices
Definition of the DFTB parameter types.
Working with the DFTB parameter types.
subroutine, public get_dftb_atom_param(dftb_parameter, name, typ, defined, z, zeff, natorb, lmax, skself, occupation, eta, energy, cutoff, xi, di, rcdisp, dudq)
...
Define the quickstep kind type and their sub types.
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
subroutine, public get_qs_kind_set(qs_kind_set, all_potential_present, tnadd_potential_present, gth_potential_present, sgp_potential_present, paw_atom_present, dft_plus_u_atom_present, maxcgf, maxsgf, maxco, maxco_proj, maxgtops, maxlgto, maxlprj, maxnset, maxsgf_set, ncgf, npgf, nset, nsgf, nshell, maxpol, maxlppl, maxlppnl, maxppnl, nelectron, maxder, max_ngrid_rad, max_sph_harm, maxg_iso_not0, lmax_rho0, basis_rcut, basis_type, total_zeff_corr, npgf_seg)
Get attributes of an atomic kind set.
Definition and initialisation of the mo data type.
subroutine, public wfn_restart_file_name(filename, exist, section, logger, kp, xas, rtp)
...
subroutine, public write_mo_set_low(mo_array, qs_kind_set, particle_set, ires, rt_mos)
...
subroutine, public read_rt_mos_from_restart(mo_array, rt_mos, atomic_kind_set, qs_kind_set, particle_set, para_env, id_nr, multiplicity, dft_section)
...
subroutine, public write_dm_binary_restart(mo_array, dft_section, tmpl_matrix)
calculates density matrix from mo set and writes the density matrix into a binary restart file
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)
subroutine, public write_mo_set_to_restart(mo_array, particle_set, dft_section, qs_kind_set)
...
subroutine, public read_mo_set_from_restart(mo_array, atomic_kind_set, qs_kind_set, particle_set, para_env, id_nr, multiplicity, dft_section, natom_mismatch, cdft, out_unit)
...
subroutine, public read_mos_restart_low(mos, para_env, qs_kind_set, particle_set, natom, rst_unit, multiplicity, rt_mos, natom_mismatch)
Reading the mos from apreviously defined restart file.
subroutine, public write_rt_mos_to_restart(mo_array, rt_mos, particle_set, dft_section, qs_kind_set)
...
Set occupation of molecular orbitals.
Definition and initialisation of the mo data type.
subroutine, public get_mo_set(mo_set, maxocc, homo, lfomo, nao, nelectron, n_el_f, nmo, eigenvalues, occupation_numbers, mo_coeff, mo_coeff_b, uniform_occupation, kts, mu, flexible_electron_count)
Get the components of a MO set data structure.
Provides all information about an atomic kind.
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment
Provides all information about a quickstep kind.