(git:37370f3)
Loading...
Searching...
No Matches
qs_mo_io.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Definition and initialisation of the mo data type.
10!> \par History
11!> - adapted to the new QS environment data structure (02.04.2002,MK)
12!> - set_mo_occupation added (17.04.02,MK)
13!> - correct_mo_eigenvalues added (18.04.02,MK)
14!> - calculate_density_matrix moved from qs_scf to here (22.04.02,MK)
15!> - mo_set_p_type added (23.04.02,MK)
16!> - PRIVATE attribute set for TYPE mo_set_type (23.04.02,MK)
17!> - started conversion to LSD (1.2003, Joost VandeVondele)
18!> - Split of from qs_mo_types (07.2014, JGH)
19!> \author Matthias Krack (09.05.2001,MK)
20! **************************************************************************************************
21MODULE qs_mo_io
22
23 USE atomic_kind_types, ONLY: get_atomic_kind
24 USE basis_set_types, ONLY: get_gto_basis_set,&
25 gto_basis_set_p_type,&
26 gto_basis_set_type
27 USE cp_cfm_types, ONLY: cp_cfm_get_submatrix,&
28 cp_cfm_type
29 USE cp_dbcsr_api, ONLY: dbcsr_binary_write,&
30 dbcsr_create,&
31 dbcsr_p_type,&
32 dbcsr_release,&
33 dbcsr_type
34 USE cp_dbcsr_contrib, ONLY: dbcsr_checksum
35 USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
36 copy_fm_to_dbcsr,&
37 dbcsr_deallocate_matrix_set
38 USE cp_dbcsr_output, ONLY: cp_dbcsr_write_sparse_matrix
39 USE cp_files, ONLY: close_file,&
40 open_file
41 USE cp_fm_types, ONLY: cp_fm_get_info,&
42 cp_fm_get_submatrix,&
43 cp_fm_set_all,&
44 cp_fm_set_submatrix,&
45 cp_fm_to_fm,&
46 cp_fm_type,&
47 cp_fm_write_unformatted
48 USE cp_log_handling, ONLY: cp_get_default_logger,&
49 cp_logger_get_default_unit_nr,&
50 cp_logger_type,&
51 cp_to_string
52 USE cp_output_handling, ONLY: cp_p_file,&
53 cp_print_key_finished_output,&
54 cp_print_key_generate_filename,&
55 cp_print_key_should_output,&
56 cp_print_key_unit_nr
57 USE input_section_types, ONLY: section_vals_get_subs_vals,&
58 section_vals_type,&
59 section_vals_val_get
60 USE kahan_sum, ONLY: accurate_sum
61 USE kinds, ONLY: default_path_length,&
62 default_string_length,&
63 dp
64 USE message_passing, ONLY: mp_para_env_type
65 USE orbital_pointers, ONLY: indco,&
66 nco,&
67 nso
68 USE orbital_symbols, ONLY: cgf_symbol,&
69 sgf_symbol
70 USE orbital_transformation_matrices, ONLY: orbtramat
71 USE particle_types, ONLY: particle_type
72 USE physcon, ONLY: evolt
73 USE qs_density_matrices, ONLY: calculate_density_matrix
74 USE qs_dftb_types, ONLY: qs_dftb_atom_type
75 USE qs_dftb_utils, ONLY: get_dftb_atom_param
76 USE qs_environment_types, ONLY: get_qs_env,&
77 qs_environment_type
78 USE qs_kind_types, ONLY: get_qs_kind,&
79 get_qs_kind_set,&
80 qs_kind_type
81 USE qs_ks_types, ONLY: qs_ks_env_type
82 USE qs_mo_methods, ONLY: calculate_subspace_eigenvalues
83 USE qs_mo_occupation, ONLY: set_mo_occupation
84 USE qs_mo_types, ONLY: get_mo_set,&
85 mo_set_type
86 USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type,&
87 release_neighbor_list_sets
88 USE qs_neighbor_lists, ONLY: setup_neighbor_list
89 USE qs_overlap, ONLY: build_overlap_matrix_simple
90#include "./base/base_uses.f90"
91
92 IMPLICIT NONE
93
94 PRIVATE
95
96 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_mo_io'
97
98 PUBLIC :: wfn_restart_file_name, &
99 write_rt_mos_to_restart, &
100 read_rt_mos_from_restart, &
101 write_dm_binary_restart, &
102 write_mo_set_to_output_unit, &
103 write_mo_set_to_restart, &
104 read_mo_set_from_restart, &
105 read_mos_restart_low, &
106 write_mo_set_low
107
108CONTAINS
109
110! **************************************************************************************************
111!> \brief ...
112!> \param mo_array ...
113!> \param particle_set ...
114!> \param dft_section ...
115!> \param qs_kind_set ...
116!> \param matrix_ks ...
117! **************************************************************************************************
118 SUBROUTINE write_mo_set_to_restart(mo_array, particle_set, dft_section, qs_kind_set, matrix_ks)
119
120 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mo_array
121 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
122 TYPE(section_vals_type), POINTER :: dft_section
123 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
124 TYPE(dbcsr_p_type), DIMENSION(:), OPTIONAL, &
125 POINTER :: matrix_ks
126
127 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_mo_set_to_restart'
128 CHARACTER(LEN=30), DIMENSION(2), PARAMETER :: &
129 keys = ["SCF%PRINT%RESTART_HISTORY", "SCF%PRINT%RESTART "]
130
131 INTEGER :: handle, ikey, ires, ispin
132 TYPE(cp_logger_type), POINTER :: logger
133
134 CALL timeset(routinen, handle)
135
136 logger => cp_get_default_logger()
137
138 IF (btest(cp_print_key_should_output(logger%iter_info, &
139 dft_section, keys(1)), cp_p_file) .OR. &
140 btest(cp_print_key_should_output(logger%iter_info, &
141 dft_section, keys(2)), cp_p_file)) THEN
142
143 IF (mo_array(1)%use_mo_coeff_b) THEN
144 ! we are using the dbcsr mo_coeff
145 ! we copy it to the fm for anycase
146 DO ispin = 1, SIZE(mo_array)
147 cpassert(ASSOCIATED(mo_array(ispin)%mo_coeff_b))
148 CALL copy_dbcsr_to_fm(mo_array(ispin)%mo_coeff_b, &
149 mo_array(ispin)%mo_coeff) !fm->dbcsr
150 END DO
151 END IF
152
153 DO ikey = 1, SIZE(keys)
154 IF (btest(cp_print_key_should_output(logger%iter_info, &
155 dft_section, keys(ikey)), cp_p_file)) THEN
156 ires = cp_print_key_unit_nr(logger, dft_section, keys(ikey), &
157 extension=".wfn", file_status="REPLACE", file_action="WRITE", &
158 do_backup=.true., file_form="UNFORMATTED")
159 IF (PRESENT(matrix_ks)) THEN
160 CALL write_mo_set_low(mo_array, particle_set=particle_set, qs_kind_set=qs_kind_set, &
161 ires=ires, matrix_ks=matrix_ks)
162 ELSE
163 CALL write_mo_set_low(mo_array, particle_set=particle_set, qs_kind_set=qs_kind_set, &
164 ires=ires)
165 END IF
166 CALL cp_print_key_finished_output(ires, logger, dft_section, trim(keys(ikey)))
167 END IF
168 END DO
169 END IF
170
171 CALL timestop(handle)
172
173 END SUBROUTINE write_mo_set_to_restart
174
175! **************************************************************************************************
176!> \brief calculates density matrix from mo set and writes the density matrix
177!> into a binary restart file
178!> \param mo_array mos
179!> \param dft_section dft input section
180!> \param tmpl_matrix template dbcsr matrix
181!> \author Mohammad Hossein Bani-Hashemian
182! **************************************************************************************************
183 SUBROUTINE write_dm_binary_restart(mo_array, dft_section, tmpl_matrix)
184
185 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mo_array
186 TYPE(section_vals_type), POINTER :: dft_section
187 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: tmpl_matrix
188
189 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_dm_binary_restart'
190
191 CHARACTER(LEN=default_path_length) :: file_name, project_name
192 INTEGER :: handle, ispin, unit_nr
193 LOGICAL :: do_dm_restart
194 REAL(kind=dp) :: cs_pos
195 TYPE(cp_logger_type), POINTER :: logger
196 TYPE(dbcsr_type), POINTER :: matrix_p_tmp
197
198 CALL timeset(routinen, handle)
199 logger => cp_get_default_logger()
200 IF (logger%para_env%is_source()) THEN
201 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
202 ELSE
203 unit_nr = -1
204 END IF
205
206 project_name = logger%iter_info%project_name
207 CALL section_vals_val_get(dft_section, "SCF%PRINT%DM_RESTART_WRITE", l_val=do_dm_restart)
208 NULLIFY (matrix_p_tmp)
209
210 IF (do_dm_restart) THEN
211 ALLOCATE (matrix_p_tmp)
212 DO ispin = 1, SIZE(mo_array)
213 CALL dbcsr_create(matrix_p_tmp, template=tmpl_matrix(ispin)%matrix, name="DM RESTART")
214
215 IF (.NOT. ASSOCIATED(mo_array(ispin)%mo_coeff_b)) cpabort("mo_coeff_b NOT ASSOCIATED")
216
217 CALL copy_fm_to_dbcsr(mo_array(ispin)%mo_coeff, mo_array(ispin)%mo_coeff_b)
218 CALL calculate_density_matrix(mo_array(ispin), matrix_p_tmp, &
219 use_dbcsr=.true., retain_sparsity=.false.)
220
221 WRITE (file_name, '(A,I0,A)') trim(project_name)//"_SCF_DM_SPIN_", ispin, "_RESTART.dm"
222 cs_pos = dbcsr_checksum(matrix_p_tmp, pos=.true.)
223 IF (unit_nr > 0) THEN
224 WRITE (unit_nr, '(T2,A,E20.8)') "Writing restart DM "//trim(file_name)//" with checksum: ", cs_pos
225 END IF
226 CALL dbcsr_binary_write(matrix_p_tmp, file_name)
227
228 CALL dbcsr_release(matrix_p_tmp)
229 END DO
230 DEALLOCATE (matrix_p_tmp)
231 END IF
232
233 CALL timestop(handle)
234
235 END SUBROUTINE write_dm_binary_restart
236
237! **************************************************************************************************
238!> \brief ...
239!> \param mo_array ...
240!> \param rt_mos ...
241!> \param particle_set ...
242!> \param dft_section ...
243!> \param qs_kind_set ...
244! **************************************************************************************************
245 SUBROUTINE write_rt_mos_to_restart(mo_array, rt_mos, particle_set, dft_section, qs_kind_set)
246
247 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mo_array
248 TYPE(cp_fm_type), DIMENSION(:), INTENT(IN) :: rt_mos
249 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
250 TYPE(section_vals_type), POINTER :: dft_section
251 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
252
253 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_rt_mos_to_restart'
254 CHARACTER(LEN=43), DIMENSION(2), PARAMETER :: keys = [ &
255 "REAL_TIME_PROPAGATION%PRINT%RESTART_HISTORY", &
256 "REAL_TIME_PROPAGATION%PRINT%RESTART "]
257
258 INTEGER :: handle, ikey, ires
259 TYPE(cp_logger_type), POINTER :: logger
260
261 CALL timeset(routinen, handle)
262 logger => cp_get_default_logger()
263
264 IF (btest(cp_print_key_should_output(logger%iter_info, &
265 dft_section, keys(1)), cp_p_file) .OR. &
266 btest(cp_print_key_should_output(logger%iter_info, &
267 dft_section, keys(2)), cp_p_file)) THEN
268
269 DO ikey = 1, SIZE(keys)
270
271 IF (btest(cp_print_key_should_output(logger%iter_info, &
272 dft_section, keys(ikey)), cp_p_file)) THEN
273 ires = cp_print_key_unit_nr(logger, dft_section, keys(ikey), &
274 extension=".rtpwfn", file_status="REPLACE", file_action="WRITE", &
275 do_backup=.true., file_form="UNFORMATTED")
276 CALL write_mo_set_low(mo_array, qs_kind_set=qs_kind_set, particle_set=particle_set, &
277 ires=ires, rt_mos=rt_mos)
278 CALL cp_print_key_finished_output(ires, logger, dft_section, trim(keys(ikey)))
279 END IF
280 END DO
281 END IF
282
283 CALL timestop(handle)
284
285 END SUBROUTINE write_rt_mos_to_restart
286
287! **************************************************************************************************
288!> \brief ...
289!> \param mo_array ...
290!> \param qs_kind_set ...
291!> \param particle_set ...
292!> \param ires ...
293!> \param rt_mos ...
294!> \param matrix_ks ...
295! **************************************************************************************************
296 SUBROUTINE write_mo_set_low(mo_array, qs_kind_set, particle_set, ires, rt_mos, matrix_ks)
297
298 TYPE(mo_set_type), DIMENSION(:), INTENT(IN) :: mo_array
299 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
300 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
301 INTEGER :: ires
302 TYPE(cp_fm_type), DIMENSION(:), INTENT(IN), &
303 OPTIONAL :: rt_mos
304 TYPE(dbcsr_p_type), DIMENSION(:), OPTIONAL, &
305 POINTER :: matrix_ks
306
307 CHARACTER(LEN=*), PARAMETER :: routinen = 'write_mo_set_low'
308
309 INTEGER :: handle, iatom, ikind, imat, iset, &
310 ishell, ispin, lmax, lshell, &
311 max_block, nao, natom, nmo, nset, &
312 nset_max, nshell_max, nspin
313 INTEGER, DIMENSION(:), POINTER :: nset_info, nshell
314 INTEGER, DIMENSION(:, :), POINTER :: l, nshell_info
315 INTEGER, DIMENSION(:, :, :), POINTER :: nso_info
316 REAL(kind=dp), DIMENSION(:), POINTER :: mo_eigenvalues, mo_occupation_numbers
317 TYPE(cp_fm_type), POINTER :: mo_coeff
318 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
319 TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
320
321 CALL timeset(routinen, handle)
322
323 NULLIFY (mo_coeff)
324 NULLIFY (mo_eigenvalues)
325 NULLIFY (mo_occupation_numbers)
326
327 nspin = SIZE(mo_array)
328 nao = mo_array(1)%nao
329
330 IF (ires > 0) THEN
331 ! Create some info about the basis set first
332 natom = SIZE(particle_set, 1)
333 nset_max = 0
334 nshell_max = 0
335
336 DO iatom = 1, natom
337 NULLIFY (orb_basis_set, dftb_parameter)
338 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
339 CALL get_qs_kind(qs_kind_set(ikind), &
340 basis_set=orb_basis_set, &
341 dftb_parameter=dftb_parameter)
342 IF (ASSOCIATED(orb_basis_set)) THEN
343 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
344 nset=nset, &
345 nshell=nshell, &
346 l=l)
347 nset_max = max(nset_max, nset)
348 DO iset = 1, nset
349 nshell_max = max(nshell_max, nshell(iset))
350 END DO
351 ELSE IF (ASSOCIATED(dftb_parameter)) THEN
352 CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
353 nset_max = max(nset_max, 1)
354 nshell_max = max(nshell_max, lmax + 1)
355 ELSE
356 ! We assume here an atom without a basis set
357 ! CPABORT("Unknown basis type. ")
358 END IF
359 END DO
360
361 ALLOCATE (nso_info(nshell_max, nset_max, natom))
362 nso_info(:, :, :) = 0
363
364 ALLOCATE (nshell_info(nset_max, natom))
365 nshell_info(:, :) = 0
366
367 ALLOCATE (nset_info(natom))
368 nset_info(:) = 0
369
370 DO iatom = 1, natom
371 NULLIFY (orb_basis_set, dftb_parameter)
372 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
373 CALL get_qs_kind(qs_kind_set(ikind), &
374 basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
375 IF (ASSOCIATED(orb_basis_set)) THEN
376 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
377 nset=nset, &
378 nshell=nshell, &
379 l=l)
380 nset_info(iatom) = nset
381 DO iset = 1, nset
382 nshell_info(iset, iatom) = nshell(iset)
383 DO ishell = 1, nshell(iset)
384 lshell = l(ishell, iset)
385 nso_info(ishell, iset, iatom) = nso(lshell)
386 END DO
387 END DO
388 ELSE IF (ASSOCIATED(dftb_parameter)) THEN
389 CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
390 nset_info(iatom) = 1
391 nshell_info(1, iatom) = lmax + 1
392 DO ishell = 1, lmax + 1
393 lshell = ishell - 1
394 nso_info(ishell, 1, iatom) = nso(lshell)
395 END DO
396 ELSE
397 ! We assume here an atom without a basis set
398 ! CPABORT("Unknown basis type. ")
399 END IF
400 END DO
401
402 WRITE (ires) natom, nspin, nao, nset_max, nshell_max
403 WRITE (ires) nset_info
404 WRITE (ires) nshell_info
405 WRITE (ires) nso_info
406
407 DEALLOCATE (nset_info)
408
409 DEALLOCATE (nshell_info)
410
411 DEALLOCATE (nso_info)
412 END IF
413
414 ! Use the ScaLAPACK block size as a default for buffering columns
415 CALL cp_fm_get_info(mo_array(1)%mo_coeff, ncol_block=max_block)
416 DO ispin = 1, nspin
417 mo_coeff => mo_array(ispin)%mo_coeff
418 nmo = mo_array(ispin)%nmo
419 IF (nmo > 0) THEN
420 mo_eigenvalues => mo_array(ispin)%eigenvalues
421 mo_occupation_numbers => mo_array(ispin)%occupation_numbers
422 IF (PRESENT(matrix_ks)) THEN
423 ! With OT: use the Kohn-Sham matrix for the update of the MO eigenvalues
424 CALL calculate_subspace_eigenvalues(orbitals=mo_coeff, &
425 ks_matrix=matrix_ks(ispin)%matrix, &
426 evals_arg=mo_eigenvalues)
427 END IF
428 IF (ires > 0) THEN
429 WRITE (ires) nmo, &
430 mo_array(ispin)%homo, &
431 mo_array(ispin)%lfomo, &
432 mo_array(ispin)%nelectron
433 WRITE (ires) mo_eigenvalues(1:nmo), mo_occupation_numbers(1:nmo)
434 END IF
435 END IF
436 IF (PRESENT(rt_mos)) THEN
437 DO imat = 2*ispin - 1, 2*ispin
438 CALL cp_fm_write_unformatted(rt_mos(imat), ires)
439 END DO
440 ELSE
441 CALL cp_fm_write_unformatted(mo_coeff, ires)
442 END IF
443 END DO
444
445 CALL timestop(handle)
446
447 END SUBROUTINE write_mo_set_low
448
449! **************************************************************************************************
450!> \brief ...
451!> \param filename ...
452!> \param exist ...
453!> \param section ...
454!> \param logger ...
455!> \param kp ...
456!> \param xas ...
457!> \param rtp ...
458! **************************************************************************************************
459 SUBROUTINE wfn_restart_file_name(filename, exist, section, logger, kp, xas, rtp)
460 CHARACTER(LEN=default_path_length), INTENT(OUT) :: filename
461 LOGICAL, INTENT(OUT) :: exist
462 TYPE(section_vals_type), POINTER :: section
463 TYPE(cp_logger_type), POINTER :: logger
464 LOGICAL, INTENT(IN), OPTIONAL :: kp, xas, rtp
465
466 INTEGER :: n_rep_val
467 LOGICAL :: my_kp, my_rtp, my_xas
468 TYPE(section_vals_type), POINTER :: print_key
469
470 my_kp = .false.
471 my_xas = .false.
472 my_rtp = .false.
473 IF (PRESENT(kp)) my_kp = kp
474 IF (PRESENT(xas)) my_xas = xas
475 IF (PRESENT(rtp)) my_rtp = rtp
476
477 exist = .false.
478 CALL section_vals_val_get(section, "WFN_RESTART_FILE_NAME", n_rep_val=n_rep_val)
479 IF (n_rep_val > 0) THEN
480 CALL section_vals_val_get(section, "WFN_RESTART_FILE_NAME", c_val=filename)
481 ELSE
482 IF (my_xas) THEN
483 ! try to read from the filename that is generated automatically from the printkey
484 print_key => section_vals_get_subs_vals(section, "PRINT%RESTART")
485 filename = cp_print_key_generate_filename(logger, print_key, &
486 extension="", my_local=.false.)
487 ELSE IF (my_rtp) THEN
488 ! try to read from the filename that is generated automatically from the printkey
489 print_key => section_vals_get_subs_vals(section, "REAL_TIME_PROPAGATION%PRINT%RESTART")
490 filename = cp_print_key_generate_filename(logger, print_key, &
491 extension=".rtpwfn", my_local=.false.)
492 ELSE IF (my_kp) THEN
493 ! try to read from the filename that is generated automatically from the printkey
494 print_key => section_vals_get_subs_vals(section, "SCF%PRINT%RESTART")
495 filename = cp_print_key_generate_filename(logger, print_key, &
496 extension=".kp", my_local=.false.)
497 ELSE
498 ! try to read from the filename that is generated automatically from the printkey
499 print_key => section_vals_get_subs_vals(section, "SCF%PRINT%RESTART")
500 filename = cp_print_key_generate_filename(logger, print_key, &
501 extension=".wfn", my_local=.false.)
502 END IF
503 END IF
504 IF (.NOT. my_xas) THEN
505 INQUIRE (file=filename, exist=exist)
506 END IF
507
508 END SUBROUTINE wfn_restart_file_name
509
510! **************************************************************************************************
511!> \brief ...
512!> \param mo_array ...
513!> \param qs_kind_set ...
514!> \param particle_set ...
515!> \param para_env ...
516!> \param id_nr ...
517!> \param multiplicity ...
518!> \param dft_section ...
519!> \param natom_mismatch ...
520!> \param cdft ...
521!> \param out_unit ...
522! **************************************************************************************************
523 SUBROUTINE read_mo_set_from_restart(mo_array, qs_kind_set, particle_set, &
524 para_env, id_nr, multiplicity, dft_section, natom_mismatch, &
525 cdft, out_unit)
526
527 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mo_array
528 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
529 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
530 TYPE(mp_para_env_type), POINTER :: para_env
531 INTEGER, INTENT(IN) :: id_nr, multiplicity
532 TYPE(section_vals_type), POINTER :: dft_section
533 LOGICAL, INTENT(OUT), OPTIONAL :: natom_mismatch
534 LOGICAL, INTENT(IN), OPTIONAL :: cdft
535 INTEGER, INTENT(IN), OPTIONAL :: out_unit
536
537 CHARACTER(LEN=*), PARAMETER :: routinen = 'read_mo_set_from_restart'
538
539 CHARACTER(LEN=default_path_length) :: file_name
540 INTEGER :: handle, ispin, my_out_unit, natom, &
541 nspin, restart_unit
542 LOGICAL :: exist, my_cdft
543 TYPE(cp_logger_type), POINTER :: logger
544
545 CALL timeset(routinen, handle)
546 logger => cp_get_default_logger()
547 my_cdft = .false.
548 IF (PRESENT(cdft)) my_cdft = cdft
549 my_out_unit = -1
550 IF (PRESENT(out_unit)) my_out_unit = out_unit
551
552 nspin = SIZE(mo_array)
553 restart_unit = -1
554
555 IF (para_env%is_source()) THEN
556
557 natom = SIZE(particle_set, 1)
558 CALL wfn_restart_file_name(file_name, exist, dft_section, logger)
559 IF (id_nr /= 0) THEN
560 ! Is it one of the backup files?
561 file_name = trim(file_name)//".bak-"//adjustl(cp_to_string(id_nr))
562 END IF
563
564 CALL open_file(file_name=file_name, &
565 file_action="READ", &
566 file_form="UNFORMATTED", &
567 file_status="OLD", &
568 unit_number=restart_unit)
569
570 END IF
571
572 CALL read_mos_restart_low(mo_array, para_env=para_env, qs_kind_set=qs_kind_set, &
573 particle_set=particle_set, natom=natom, &
574 rst_unit=restart_unit, multiplicity=multiplicity, natom_mismatch=natom_mismatch)
575
576 IF (PRESENT(natom_mismatch)) THEN
577 ! read_mos_restart_low only the io_node returns natom_mismatch, must broadcast it
578 CALL para_env%bcast(natom_mismatch)
579 IF (natom_mismatch) THEN
580 IF (para_env%is_source()) CALL close_file(unit_number=restart_unit)
581 CALL timestop(handle)
582 RETURN
583 END IF
584 END IF
585
586 ! Close restart file
587 IF (para_env%is_source()) THEN
588 IF (my_out_unit > 0) THEN
589 WRITE (unit=my_out_unit, fmt="(T2,A)") &
590 "WFN_RESTART| Restart file "//trim(file_name)//" read"
591 END IF
592 CALL close_file(unit_number=restart_unit)
593 END IF
594
595 ! CDFT has no real dft_section and does not need to print
596 IF (.NOT. my_cdft) THEN
597 DO ispin = 1, nspin
598 CALL write_mo_set_to_output_unit(mo_array(ispin), qs_kind_set, particle_set, &
599 dft_section, 4, 0, final_mos=.false.)
600 END DO
601 END IF
602
603 CALL timestop(handle)
604
605 END SUBROUTINE read_mo_set_from_restart
606
607! **************************************************************************************************
608!> \brief ...
609!> \param mo_array ...
610!> \param rt_mos ...
611!> \param qs_kind_set ...
612!> \param particle_set ...
613!> \param para_env ...
614!> \param id_nr ...
615!> \param multiplicity ...
616!> \param dft_section ...
617! **************************************************************************************************
618 SUBROUTINE read_rt_mos_from_restart(mo_array, rt_mos, qs_kind_set, particle_set, para_env, id_nr, multiplicity, dft_section)
619
620 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mo_array
621 TYPE(cp_fm_type), DIMENSION(:), POINTER :: rt_mos
622 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
623 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
624 TYPE(mp_para_env_type), POINTER :: para_env
625 INTEGER, INTENT(IN) :: id_nr, multiplicity
626 TYPE(section_vals_type), POINTER :: dft_section
627
628 CHARACTER(LEN=*), PARAMETER :: routinen = 'read_rt_mos_from_restart'
629
630 CHARACTER(LEN=default_path_length) :: file_name
631 INTEGER :: handle, ispin, natom, nspin, &
632 restart_unit, unit_nr
633 LOGICAL :: exist
634 TYPE(cp_logger_type), POINTER :: logger
635
636 CALL timeset(routinen, handle)
637 logger => cp_get_default_logger()
638
639 nspin = SIZE(mo_array)
640 restart_unit = -1
641
642 IF (para_env%is_source()) THEN
643
644 natom = SIZE(particle_set, 1)
645 CALL wfn_restart_file_name(file_name, exist, dft_section, logger, rtp=.true.)
646 IF (id_nr /= 0) THEN
647 ! Is it one of the backup files?
648 file_name = trim(file_name)//".bak-"//adjustl(cp_to_string(id_nr))
649 END IF
650
651 unit_nr = cp_logger_get_default_unit_nr(logger, local=.true.)
652 IF (unit_nr > 0) THEN
653 WRITE (unit_nr, '(T2,A)') "Read RTP restart from the file: "//trim(file_name)
654 END IF
655
656 CALL open_file(file_name=file_name, &
657 file_action="READ", &
658 file_form="UNFORMATTED", &
659 file_status="OLD", &
660 unit_number=restart_unit)
661
662 END IF
663
664 CALL read_mos_restart_low(mo_array, rt_mos=rt_mos, para_env=para_env, &
665 particle_set=particle_set, qs_kind_set=qs_kind_set, natom=natom, &
666 rst_unit=restart_unit, multiplicity=multiplicity)
667
668 ! Close restart file
669 IF (para_env%is_source()) CALL close_file(unit_number=restart_unit)
670
671 DO ispin = 1, nspin
672 CALL write_mo_set_to_output_unit(mo_array(ispin), qs_kind_set, particle_set, &
673 dft_section, 4, 0, final_mos=.false.)
674 END DO
675
676 CALL timestop(handle)
677
678 END SUBROUTINE read_rt_mos_from_restart
679
680! **************************************************************************************************
681!> \brief Reading the mos from apreviously defined restart file
682!> \param mos ...
683!> \param para_env ...
684!> \param qs_kind_set ...
685!> \param particle_set ...
686!> \param natom ...
687!> \param rst_unit ...
688!> \param multiplicity ...
689!> \param rt_mos ...
690!> \param natom_mismatch ...
691!> \par History
692!> 12.2007 created [MI]
693!> \author MI
694! **************************************************************************************************
695 SUBROUTINE read_mos_restart_low(mos, para_env, qs_kind_set, particle_set, natom, rst_unit, &
696 multiplicity, rt_mos, natom_mismatch)
697
698 TYPE(mo_set_type), DIMENSION(:), INTENT(INOUT) :: mos
699 TYPE(mp_para_env_type), POINTER :: para_env
700 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
701 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
702 INTEGER, INTENT(IN) :: natom, rst_unit
703 INTEGER, INTENT(in), OPTIONAL :: multiplicity
704 TYPE(cp_fm_type), DIMENSION(:), OPTIONAL, POINTER :: rt_mos
705 LOGICAL, INTENT(OUT), OPTIONAL :: natom_mismatch
706
707 INTEGER :: homo, homo_read, i, iatom, ikind, imat, irow, iset, iset_read, ishell, &
708 ishell_read, iso, ispin, lfomo_read, lmax, lshell, my_mult, nao, nao_read, natom_read, &
709 nelectron, nelectron_read, nmo, nmo_read, nnshell, nset, nset_max, nshell_max, nspin, &
710 nspin_read, offset_read
711 INTEGER, DIMENSION(:), POINTER :: nset_info, nshell
712 INTEGER, DIMENSION(:, :), POINTER :: l, nshell_info
713 INTEGER, DIMENSION(:, :, :), POINTER :: nso_info, offset_info
714 LOGICAL :: minbas, natom_match, use_this
715 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: eig_read, occ_read
716 REAL(kind=dp), DIMENSION(:, :), POINTER :: vecbuffer, vecbuffer_read
717 TYPE(cp_logger_type), POINTER :: logger
718 TYPE(gto_basis_set_type), POINTER :: orb_basis_set
719 TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
720
721 logger => cp_get_default_logger()
722
723 nspin = SIZE(mos)
724 nao = mos(1)%nao
725 my_mult = 0
726 IF (PRESENT(multiplicity)) my_mult = multiplicity
727
728 IF (para_env%is_source()) THEN
729 READ (rst_unit) natom_read, nspin_read, nao_read, nset_max, nshell_max
730 IF (PRESENT(rt_mos)) THEN
731 IF (nspin_read /= nspin) THEN
732 cpabort("To change nspin is not possible. ")
733 END IF
734 ELSE
735 ! we should allow for restarting with different spin settings
736 IF (nspin_read /= nspin) THEN
737 WRITE (cp_logger_get_default_unit_nr(logger), *) &
738 "READ RESTART : WARNING : nspin is not equal "
739 END IF
740 ! this case needs fixing of homo/lfomo/nelec/occupations ...
741 IF (nspin_read > nspin) THEN
742 cpabort("Reducing nspin is not possible. ")
743 END IF
744 END IF
745
746 natom_match = (natom_read == natom)
747
748 IF (natom_match) THEN ! actually do the read read
749
750 ! Let's make it possible to change the basis set
751 ALLOCATE (nso_info(nshell_max, nset_max, natom_read))
752 ALLOCATE (nshell_info(nset_max, natom_read))
753 ALLOCATE (nset_info(natom_read))
754 ALLOCATE (offset_info(nshell_max, nset_max, natom_read))
755
756 IF (nao_read /= nao) THEN
757 WRITE (cp_logger_get_default_unit_nr(logger), *) &
758 " READ RESTART : WARNING : DIFFERENT # AOs ", nao, nao_read
759 IF (PRESENT(rt_mos)) THEN
760 cpabort("To change basis is not possible. ")
761 END IF
762 END IF
763
764 READ (rst_unit) nset_info
765 READ (rst_unit) nshell_info
766 READ (rst_unit) nso_info
767
768 i = 1
769 DO iatom = 1, natom
770 DO iset = 1, nset_info(iatom)
771 DO ishell = 1, nshell_info(iset, iatom)
772 offset_info(ishell, iset, iatom) = i
773 i = i + nso_info(ishell, iset, iatom)
774 END DO
775 END DO
776 END DO
777
778 ALLOCATE (vecbuffer_read(1, nao_read))
779
780 END IF ! natom_match
781 END IF ! ionode
782
783 ! make natom_match and natom_mismatch uniform across all nodes
784 CALL para_env%bcast(natom_match)
785 IF (PRESENT(natom_mismatch)) natom_mismatch = .NOT. natom_match
786 ! handle natom_match false
787 IF (.NOT. natom_match) THEN
788 IF (PRESENT(natom_mismatch)) THEN
789 WRITE (cp_logger_get_default_unit_nr(logger), *) &
790 " READ RESTART : WARNING : DIFFERENT natom, returning ", natom, natom_read
791 RETURN
792 ELSE
793 cpabort("Incorrect number of atoms in restart file. ")
794 END IF
795 END IF
796
797 CALL para_env%bcast(nspin_read)
798
799 ALLOCATE (vecbuffer(1, nao))
800
801 DO ispin = 1, nspin
802
803 nmo = mos(ispin)%nmo
804 homo = mos(ispin)%homo
805 mos(ispin)%eigenvalues(:) = 0.0_dp
806 mos(ispin)%occupation_numbers(:) = 0.0_dp
807 CALL cp_fm_set_all(mos(ispin)%mo_coeff, 0.0_dp)
808
809 IF (para_env%is_source() .AND. (nmo > 0)) THEN
810 READ (rst_unit) nmo_read, homo_read, lfomo_read, nelectron_read
811 ALLOCATE (eig_read(nmo_read), occ_read(nmo_read))
812 eig_read = 0.0_dp
813 occ_read = 0.0_dp
814
815 nmo = min(nmo, nmo_read)
816 IF (nmo_read < nmo) THEN
817 CALL cp_warn(__location__, &
818 "The number of MOs on the restart unit is smaller than the number of "// &
819 "the allocated MOs. The MO set will be padded with zeros!")
820 END IF
821 IF (nmo_read > nmo) THEN
822 CALL cp_warn(__location__, &
823 "The number of MOs on the restart unit is greater than the number of "// &
824 "the allocated MOs. The read MO set will be truncated!")
825 END IF
826
827 READ (rst_unit) eig_read(1:nmo_read), occ_read(1:nmo_read)
828 mos(ispin)%eigenvalues(1:nmo) = eig_read(1:nmo)
829 mos(ispin)%occupation_numbers(1:nmo) = occ_read(1:nmo)
830 DEALLOCATE (eig_read, occ_read)
831
832 mos(ispin)%homo = homo_read
833 mos(ispin)%lfomo = lfomo_read
834 IF (min(homo_read, homo) > nmo) THEN
835 IF (nelectron_read == mos(ispin)%nelectron) THEN
836 CALL cp_warn(__location__, &
837 "The number of occupied MOs on the restart unit is larger than "// &
838 "the allocated MOs. The read MO set will be truncated and the occupation numbers recalculated!")
839 CALL set_mo_occupation(mo_set=mos(ispin))
840 ELSE
841 ! can not make this a warning i.e. homo must be smaller than nmo
842 ! otherwise e.g. set_mo_occupation will go out of bounds
843 cpabort("Number of occupied MOs on restart unit larger than allocated MOs. ")
844 END IF
845 END IF
846 END IF
847
848 CALL para_env%bcast(nmo)
849 CALL para_env%bcast(mos(ispin)%homo)
850 CALL para_env%bcast(mos(ispin)%lfomo)
851 CALL para_env%bcast(mos(ispin)%nelectron)
852 CALL para_env%bcast(mos(ispin)%eigenvalues)
853 CALL para_env%bcast(mos(ispin)%occupation_numbers)
854
855 IF (PRESENT(rt_mos)) THEN
856 DO imat = 2*ispin - 1, 2*ispin
857 DO i = 1, nmo
858 IF (para_env%is_source()) THEN
859 READ (rst_unit) vecbuffer
860 ELSE
861 vecbuffer(1, :) = 0.0_dp
862 END IF
863 CALL para_env%bcast(vecbuffer)
864 CALL cp_fm_set_submatrix(rt_mos(imat), &
865 vecbuffer, 1, i, nao, 1, transpose=.true.)
866 END DO
867 END DO
868 ELSE
869 DO i = 1, nmo
870 IF (para_env%is_source()) THEN
871 READ (rst_unit) vecbuffer_read
872 ! now, try to assign the read to the real vector
873 ! in case the basis set changed this involves some guessing
874 irow = 1
875 DO iatom = 1, natom
876 NULLIFY (orb_basis_set, dftb_parameter, l, nshell)
877 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
878 CALL get_qs_kind(qs_kind_set(ikind), &
879 basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
880 IF (ASSOCIATED(orb_basis_set)) THEN
881 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
882 nset=nset, &
883 nshell=nshell, &
884 l=l)
885 minbas = .false.
886 ELSE IF (ASSOCIATED(dftb_parameter)) THEN
887 CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
888 nset = 1
889 minbas = .true.
890 ELSE
891 ! assume an atom without basis set
892 ! CPABORT("Unknown basis set type. ")
893 nset = 0
894 END IF
895
896 use_this = .true.
897 iset_read = 1
898 DO iset = 1, nset
899 ishell_read = 1
900 IF (minbas) THEN
901 nnshell = lmax + 1
902 ELSE
903 nnshell = nshell(iset)
904 END IF
905 DO ishell = 1, nnshell
906 IF (minbas) THEN
907 lshell = ishell - 1
908 ELSE
909 lshell = l(ishell, iset)
910 END IF
911 IF (iset_read > nset_info(iatom)) use_this = .false.
912 IF (use_this) THEN ! avoids out of bound access of the lower line if false
913 IF (nso(lshell) == nso_info(ishell_read, iset_read, iatom)) THEN
914 offset_read = offset_info(ishell_read, iset_read, iatom)
915 ishell_read = ishell_read + 1
916 IF (ishell_read > nshell_info(iset, iatom)) THEN
917 ishell_read = 1
918 iset_read = iset_read + 1
919 END IF
920 ELSE
921 use_this = .false.
922 END IF
923 END IF
924 DO iso = 1, nso(lshell)
925 IF (use_this) THEN
926 IF (offset_read - 1 + iso < 1 .OR. offset_read - 1 + iso > nao_read) THEN
927 vecbuffer(1, irow) = 0.0_dp
928 ELSE
929 vecbuffer(1, irow) = vecbuffer_read(1, offset_read - 1 + iso)
930 END IF
931 ELSE
932 vecbuffer(1, irow) = 0.0_dp
933 END IF
934 irow = irow + 1
935 END DO
936 use_this = .true.
937 END DO
938 END DO
939 END DO
940
941 ELSE
942
943 vecbuffer(1, :) = 0.0_dp
944
945 END IF
946
947 CALL para_env%bcast(vecbuffer)
948 CALL cp_fm_set_submatrix(mos(ispin)%mo_coeff, &
949 vecbuffer, 1, i, nao, 1, transpose=.true.)
950 END DO
951 END IF
952 ! Skip extra MOs if there any
953 IF (para_env%is_source()) THEN
954 !ignore nmo = 0
955 IF (nmo > 0) THEN
956 DO i = nmo + 1, nmo_read
957 READ (rst_unit) vecbuffer_read
958 END DO
959 END IF
960 END IF
961
962 IF (.NOT. PRESENT(rt_mos)) THEN
963 IF (ispin == 1 .AND. nspin_read < nspin) THEN
964
965 mos(ispin + 1)%homo = mos(ispin)%homo
966 mos(ispin + 1)%lfomo = mos(ispin)%lfomo
967 nelectron = mos(ispin)%nelectron
968 IF (my_mult /= 1) THEN
969 CALL cp_abort(__location__, &
970 "Restarting an LSD calculation from an LDA wfn only works for multiplicity=1 (singlets).")
971 END IF
972 IF (mos(ispin + 1)%nelectron < 0) THEN
973 cpabort("LSD: too few electrons for this multiplisity. ")
974 END IF
975 mos(ispin + 1)%eigenvalues = mos(ispin)%eigenvalues
976 mos(ispin)%occupation_numbers = mos(ispin)%occupation_numbers/2.0_dp
977 mos(ispin + 1)%occupation_numbers = mos(ispin)%occupation_numbers
978 CALL cp_fm_to_fm(mos(ispin)%mo_coeff, mos(ispin + 1)%mo_coeff)
979 EXIT
980 END IF
981 END IF
982 END DO ! ispin
983
984 DEALLOCATE (vecbuffer)
985
986 IF (para_env%is_source()) THEN
987 DEALLOCATE (vecbuffer_read)
988 DEALLOCATE (offset_info)
989 DEALLOCATE (nso_info)
990 DEALLOCATE (nshell_info)
991 DEALLOCATE (nset_info)
992 END IF
993
994 END SUBROUTINE read_mos_restart_low
995
996! **************************************************************************************************
997!> \brief Write MO information to output file (eigenvalues, occupation numbers, coefficients)
998!> \param mo_set ...
999!> \param qs_kind_set ...
1000!> \param particle_set ...
1001!> \param dft_section ...
1002!> \param before Digits before the dot
1003!> \param kpoint An integer that labels the current k point, e.g. its index
1004!> \param final_mos ...
1005!> \param spin ...
1006!> \param solver_method ...
1007!> \param rtp ...
1008!> \param cpart ...
1009!> \param sim_step ...
1010!> \param umo_set ...
1011!> \param qs_env ...
1012!> \param para_env_inter_kp ...
1013!> \date 15.05.2001
1014!> \par History:
1015!> - Optionally print Cartesian MOs (20.04.2005, MK)
1016!> - Revise printout of MO information (05.05.2021, MK)
1017!> \par Variables
1018!> - after : Number of digits after point.
1019!> - before: Number of digits before point.
1020!> \author Matthias Krack (MK)
1021!> \version 1.1
1022! **************************************************************************************************
1023 SUBROUTINE write_mo_set_to_output_unit(mo_set, qs_kind_set, particle_set, &
1024 dft_section, before, kpoint, final_mos, spin, &
1025 solver_method, rtp, cpart, sim_step, umo_set, qs_env, &
1026 para_env_inter_kp)
1027
1028 TYPE(mo_set_type), INTENT(IN), OPTIONAL :: mo_set
1029 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
1030 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
1031 TYPE(section_vals_type), POINTER :: dft_section
1032 INTEGER, INTENT(IN) :: before, kpoint
1033 LOGICAL, INTENT(IN), OPTIONAL :: final_mos
1034 CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: spin
1035 CHARACTER(LEN=2), INTENT(IN), OPTIONAL :: solver_method
1036 LOGICAL, INTENT(IN), OPTIONAL :: rtp
1037 INTEGER, INTENT(IN), OPTIONAL :: cpart, sim_step
1038 TYPE(mo_set_type), INTENT(IN), OPTIONAL :: umo_set
1039 TYPE(qs_environment_type), OPTIONAL, POINTER :: qs_env
1040 TYPE(mp_para_env_type), OPTIONAL, POINTER :: para_env_inter_kp
1041
1042 CHARACTER(LEN=12) :: symbol
1043 CHARACTER(LEN=12), DIMENSION(:), POINTER :: bcgf_symbol
1044 CHARACTER(LEN=14) :: fmtstr5
1045 CHARACTER(LEN=15) :: energy_str, orbital_str, step_string
1046 CHARACTER(LEN=2) :: element_symbol, my_solver_method
1047 CHARACTER(LEN=2*default_string_length) :: name
1048 CHARACTER(LEN=21) :: vector_str
1049 CHARACTER(LEN=22) :: fmtstr4
1050 CHARACTER(LEN=24) :: fmtstr2
1051 CHARACTER(LEN=25) :: fmtstr1
1052 CHARACTER(LEN=29) :: fmtstr6
1053 CHARACTER(LEN=4) :: reim
1054 CHARACTER(LEN=40) :: fmtstr3
1055 CHARACTER(LEN=6), DIMENSION(:), POINTER :: bsgf_symbol
1056 COMPLEX(KIND=dp), ALLOCATABLE, DIMENSION(:, :) :: cmo_coeff_buf
1057 INTEGER :: after, first_mo, from, homo, iatom, icgf, ico, icol, ikind, imo, irow, iset, &
1058 isgf, ishell, iso, iw, jcol, last_mo, left, lmax, lshell, nao, natom, ncgf, ncol, nkind, &
1059 nmo, nmo_local, nset, nsgf, numo, right, scf_step, to, width
1060 INTEGER, DIMENSION(:), POINTER :: mo_index_range, nshell
1061 INTEGER, DIMENSION(:, :), POINTER :: l
1062 LOGICAL :: ionode, my_final, my_rtp, omit_headers, print_cartesian, print_cartesian_overlap, &
1063 print_eigvals, print_eigvecs, print_occup, should_output
1064 REAL(kind=dp) :: chemical_potential, gap, maxocc
1065 REAL(kind=dp), ALLOCATABLE, DIMENSION(:) :: mo_eigenvalues, mo_occupation_numbers
1066 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: cmatrix, smatrix
1067 REAL(kind=dp), DIMENSION(:), POINTER :: eigenvalues, occupation_numbers
1068 TYPE(cp_cfm_type), POINTER :: cmo_coeff
1069 TYPE(cp_fm_type), POINTER :: mo_coeff, umo_coeff
1070 TYPE(cp_logger_type), POINTER :: logger
1071 TYPE(dbcsr_p_type), DIMENSION(:), POINTER :: sro
1072 TYPE(gto_basis_set_p_type), DIMENSION(:), POINTER :: orb_basis_set_list
1073 TYPE(gto_basis_set_type), POINTER :: orb_basis_set, orbbasis
1074 TYPE(mp_para_env_type), POINTER :: para_env
1075 TYPE(neighbor_list_set_p_type), DIMENSION(:), &
1076 POINTER :: sro_list
1077 TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
1078 TYPE(qs_kind_type), POINTER :: qs_kind
1079 TYPE(qs_ks_env_type), POINTER :: ks_env
1080
1081 NULLIFY (bcgf_symbol)
1082 NULLIFY (bsgf_symbol)
1083 NULLIFY (logger)
1084 NULLIFY (mo_index_range)
1085 NULLIFY (nshell)
1086 NULLIFY (mo_coeff)
1087 NULLIFY (cmo_coeff)
1088
1089 logger => cp_get_default_logger()
1090 ionode = logger%para_env%is_source()
1091 CALL section_vals_val_get(dft_section, "PRINT%MO%EIGENVALUES", l_val=print_eigvals)
1092 CALL section_vals_val_get(dft_section, "PRINT%MO%EIGENVECTORS", l_val=print_eigvecs)
1093 CALL section_vals_val_get(dft_section, "PRINT%MO%OCCUPATION_NUMBERS", l_val=print_occup)
1094 CALL section_vals_val_get(dft_section, "PRINT%MO%CARTESIAN", l_val=print_cartesian)
1095 CALL section_vals_val_get(dft_section, "PRINT%MO%MO_INDEX_RANGE", i_vals=mo_index_range)
1096 CALL section_vals_val_get(dft_section, "PRINT%MO%NDIGITS", i_val=after)
1097 CALL section_vals_val_get(dft_section, "PRINT%MO%CARTESIAN_OVERLAP", l_val=print_cartesian_overlap)
1098 after = min(max(after, 1), 16)
1099
1100 ! Do we print the final MO information after SCF convergence is reached (default: no)
1101 IF (PRESENT(final_mos)) THEN
1102 my_final = final_mos
1103 ELSE
1104 my_final = .false.
1105 END IF
1106
1107 ! complex MOS for RTP, no eigenvalues
1108 my_rtp = .false.
1109 IF (PRESENT(rtp)) THEN
1110 my_rtp = rtp
1111 ! print the first time step if MO print required
1112 should_output = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
1113 "PRINT%MO"), cp_p_file) &
1114 .OR. (sim_step == 1)
1115 ELSE
1116 should_output = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
1117 "PRINT%MO"), cp_p_file) .OR. my_final
1118 END IF
1119
1120 IF ((.NOT. should_output) .OR. (.NOT. (print_eigvals .OR. print_eigvecs .OR. print_occup))) RETURN
1121
1122 IF (my_rtp) THEN
1123 cpassert(PRESENT(sim_step))
1124 cpassert(PRESENT(cpart))
1125 scf_step = sim_step
1126 IF (cpart == 0) THEN
1127 reim = "IMAG"
1128 ELSE
1129 reim = "REAL"
1130 END IF
1131 print_eigvals = .false.
1132 ELSE
1133 scf_step = max(0, logger%iter_info%iteration(logger%iter_info%n_rlevel) - 1)
1134 END IF
1135
1136 IF (.NOT. my_final) THEN
1137 IF (.NOT. my_rtp) THEN
1138 step_string = " AFTER SCF STEP"
1139 ELSE
1140 step_string = " AFTER RTP STEP"
1141 END IF
1142 END IF
1143
1144 IF (PRESENT(solver_method)) THEN
1145 my_solver_method = solver_method
1146 ELSE
1147 ! Traditional diagonalization is assumed as default solver method
1148 my_solver_method = "TD"
1149 END IF
1150
1151 ! Retrieve MO information
1152 IF (PRESENT(para_env_inter_kp)) THEN
1153 cpassert(ASSOCIATED(para_env_inter_kp))
1154 cpassert(.NOT. PRESENT(umo_set))
1155 nmo_local = 0
1156 homo = 0
1157 maxocc = 0.0_dp
1158 chemical_potential = 0.0_dp
1159 NULLIFY (cmo_coeff, eigenvalues, occupation_numbers, mo_coeff)
1160 IF (PRESENT(mo_set)) THEN
1161 CALL get_mo_set(mo_set=mo_set, eigenvalues=eigenvalues, &
1162 occupation_numbers=occupation_numbers, mo_coeff=mo_coeff, &
1163 homo=homo, maxocc=maxocc, mu=chemical_potential, nmo=nmo_local, &
1164 cmo_coeff=cmo_coeff)
1165 END IF
1166 nmo = nmo_local
1167 CALL para_env_inter_kp%max(nmo)
1168 CALL para_env_inter_kp%sum(homo)
1169 CALL para_env_inter_kp%sum(maxocc)
1170 CALL para_env_inter_kp%sum(chemical_potential)
1171 ALLOCATE (mo_eigenvalues(nmo), mo_occupation_numbers(nmo))
1172 mo_eigenvalues = 0.0_dp
1173 mo_occupation_numbers = 0.0_dp
1174 IF (nmo_local > 0 .AND. PRESENT(mo_set)) THEN
1175 mo_eigenvalues(1:nmo_local) = eigenvalues(1:nmo_local)
1176 mo_occupation_numbers(1:nmo_local) = occupation_numbers(1:nmo_local)
1177 END IF
1178 CALL para_env_inter_kp%sum(mo_eigenvalues)
1179 CALL para_env_inter_kp%sum(mo_occupation_numbers)
1180 IF (print_eigvecs) THEN
1181 CALL get_qs_kind_set(qs_kind_set, nsgf=nao)
1182 ALLOCATE (smatrix(nao, nmo))
1183 smatrix = 0.0_dp
1184 IF (nmo_local > 0 .AND. PRESENT(mo_set)) THEN
1185 IF (ASSOCIATED(cmo_coeff)) THEN
1186 ALLOCATE (cmo_coeff_buf(nao, nmo_local))
1187 CALL cp_cfm_get_submatrix(cmo_coeff, cmo_coeff_buf)
1188 ! Preserve the legacy text output by printing only the real part.
1189 smatrix(:, 1:nmo_local) = real(cmo_coeff_buf, kind=dp)
1190 DEALLOCATE (cmo_coeff_buf)
1191 ELSE
1192 CALL cp_fm_get_submatrix(mo_coeff, smatrix(:, 1:nmo_local))
1193 END IF
1194 END IF
1195 CALL para_env_inter_kp%sum(smatrix)
1196 END IF
1197 numo = 0
1198 ELSE
1199 cpassert(PRESENT(mo_set))
1200 CALL get_mo_set(mo_set=mo_set, &
1201 mo_coeff=mo_coeff, &
1202 eigenvalues=eigenvalues, &
1203 occupation_numbers=occupation_numbers, &
1204 homo=homo, &
1205 maxocc=maxocc, &
1206 nao=nao, &
1207 nmo=nmo, &
1208 mu=chemical_potential)
1209 IF (PRESENT(umo_set)) THEN
1210 CALL get_mo_set(mo_set=umo_set, &
1211 mo_coeff=umo_coeff, &
1212 nmo=numo)
1213 nmo = nmo + numo
1214 ELSE
1215 numo = 0
1216 END IF
1217 ALLOCATE (mo_eigenvalues(nmo), mo_occupation_numbers(nmo))
1218 mo_eigenvalues(1:nmo - numo) = eigenvalues(1:nmo - numo)
1219 mo_occupation_numbers = 0.0_dp
1220 mo_occupation_numbers(1:nmo - numo) = occupation_numbers(1:nmo - numo)
1221 IF (numo > 0) THEN
1222 CALL get_mo_set(mo_set=umo_set, eigenvalues=eigenvalues)
1223 mo_eigenvalues(nmo - numo + 1:nmo) = eigenvalues(1:numo)
1224 END IF
1225 END IF
1226
1227 IF (print_eigvecs) THEN
1228 IF (.NOT. ALLOCATED(smatrix)) THEN
1229 ALLOCATE (smatrix(nao, nmo))
1230 CALL cp_fm_get_submatrix(mo_coeff, smatrix(1:nao, 1:nmo - numo))
1231 IF (numo > 0) THEN
1232 CALL cp_fm_get_submatrix(umo_coeff, smatrix(1:nao, nmo - numo + 1:nmo))
1233 END IF
1234 END IF
1235 IF (.NOT. ionode) THEN
1236 DEALLOCATE (smatrix)
1237 END IF
1238 END IF
1239
1240 IF (PRESENT(qs_env)) THEN
1241 IF (ASSOCIATED(qs_env) .AND. my_final .AND. print_cartesian_overlap) THEN
1242 NULLIFY (qs_kind_set)
1243 CALL get_qs_env(qs_env=qs_env, qs_kind_set=qs_kind_set)
1244 nkind = SIZE(qs_kind_set)
1245
1246 IF (btest(cp_print_key_should_output(logger%iter_info, &
1247 qs_env%input, "DFT%PRINT%AO_MATRICES/OVERLAP"), cp_p_file)) THEN
1248 ALLOCATE (orb_basis_set_list(nkind))
1249 DO ikind = 1, nkind
1250 qs_kind => qs_kind_set(ikind)
1251 NULLIFY (orb_basis_set_list(ikind)%gto_basis_set)
1252 NULLIFY (orbbasis)
1253 CALL get_qs_kind(qs_kind=qs_kind, basis_set=orbbasis, basis_type="ORB")
1254 IF (ASSOCIATED(orbbasis)) orb_basis_set_list(ikind)%gto_basis_set => orbbasis
1255 END DO
1256 NULLIFY (sro_list)
1257 CALL setup_neighbor_list(sro_list, orb_basis_set_list, qs_env=qs_env)
1258 NULLIFY (sro)
1259 NULLIFY (para_env)
1260 CALL get_qs_env(qs_env, ks_env=ks_env, para_env=para_env)
1261 CALL build_overlap_matrix_simple(ks_env, sro, &
1262 orb_basis_set_list, orb_basis_set_list, sro_list, .true.)
1263 CALL release_neighbor_list_sets(sro_list)
1264
1265 iw = cp_print_key_unit_nr(logger, qs_env%input, "DFT%PRINT%AO_MATRICES/OVERLAP", &
1266 extension=".Log")
1267 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%OMIT_HEADERS", l_val=omit_headers)
1268 CALL section_vals_val_get(qs_env%input, "DFT%PRINT%AO_MATRICES%NDIGITS", i_val=after)
1269 after = min(max(after, 1), 16)
1270 IF (ASSOCIATED(sro)) THEN
1271 CALL cp_dbcsr_write_sparse_matrix(sro(1)%matrix, 4, after, qs_env, para_env, &
1272 output_unit=iw, omit_headers=omit_headers, &
1273 cartesian_basis=.true.)
1274 END IF
1275 CALL cp_print_key_finished_output(iw, logger, qs_env%input, &
1276 "DFT%PRINT%AO_MATRICES/OVERLAP")
1277 IF (ASSOCIATED(sro)) CALL dbcsr_deallocate_matrix_set(sro)
1278 DEALLOCATE (orb_basis_set_list)
1279 END IF
1280 END IF
1281 END IF
1282
1283 iw = cp_print_key_unit_nr(logger, dft_section, "PRINT%MO", &
1284 ignore_should_output=should_output, &
1285 extension=".MOLog")
1286
1287 IF (iw > 0) THEN
1288
1289 natom = SIZE(particle_set)
1290 CALL get_qs_kind_set(qs_kind_set, ncgf=ncgf, nsgf=nsgf)
1291
1292 ! Definition of the variable formats
1293
1294 fmtstr1 = "(T2,A,21X, ( X,I5, X))"
1295 fmtstr2 = "(T2,A,21X, (1X,F . ))"
1296 fmtstr3 = "(T2,A,I5,1X,I5,1X,A,1X,A6, (1X,F . ))"
1297
1298 width = before + after + 3
1299 ncol = int(56/width)
1300
1301 right = max((after - 2), 1)
1302 left = width - right - 5
1303
1304 WRITE (unit=fmtstr1(11:12), fmt="(I2)") ncol
1305 WRITE (unit=fmtstr1(14:15), fmt="(I2)") left
1306 WRITE (unit=fmtstr1(21:22), fmt="(I2)") right
1307
1308 WRITE (unit=fmtstr2(11:12), fmt="(I2)") ncol
1309 WRITE (unit=fmtstr2(18:19), fmt="(I2)") width - 1
1310 WRITE (unit=fmtstr2(21:22), fmt="(I2)") after
1311
1312 WRITE (unit=fmtstr3(27:28), fmt="(I2)") ncol
1313 WRITE (unit=fmtstr3(34:35), fmt="(I2)") width - 1
1314 WRITE (unit=fmtstr3(37:38), fmt="(I2)") after
1315
1316 IF (my_final .OR. (my_solver_method == "TD")) THEN
1317 energy_str = "EIGENVALUES"
1318 vector_str = "EIGENVECTORS"
1319 ELSE
1320 energy_str = "ENERGIES"
1321 vector_str = "COEFFICIENTS"
1322 END IF
1323
1324 IF (my_rtp) THEN
1325 energy_str = "ZEROS"
1326 vector_str = trim(reim)//" RTP COEFFICIENTS"
1327 END IF
1328
1329 IF (print_eigvecs) THEN
1330
1331 IF (print_cartesian) THEN
1332
1333 orbital_str = "CARTESIAN"
1334
1335 ALLOCATE (cmatrix(ncgf, ncgf))
1336 cmatrix = 0.0_dp
1337
1338 ! Transform spherical MOs to Cartesian MOs
1339 icgf = 1
1340 isgf = 1
1341 DO iatom = 1, natom
1342 NULLIFY (orb_basis_set, dftb_parameter)
1343 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
1344 CALL get_qs_kind(qs_kind_set(ikind), &
1345 basis_set=orb_basis_set, &
1346 dftb_parameter=dftb_parameter)
1347 IF (ASSOCIATED(orb_basis_set)) THEN
1348 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1349 nset=nset, &
1350 nshell=nshell, &
1351 l=l)
1352 DO iset = 1, nset
1353 DO ishell = 1, nshell(iset)
1354 lshell = l(ishell, iset)
1355 CALL dgemm("T", "N", nco(lshell), nmo, nso(lshell), 1.0_dp, &
1356 orbtramat(lshell)%c2s, nso(lshell), &
1357 smatrix(isgf, 1), nsgf, 0.0_dp, &
1358 cmatrix(icgf, 1), ncgf)
1359 icgf = icgf + nco(lshell)
1360 isgf = isgf + nso(lshell)
1361 END DO
1362 END DO
1363 ELSE IF (ASSOCIATED(dftb_parameter)) THEN
1364 CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
1365 DO ishell = 1, lmax + 1
1366 lshell = ishell - 1
1367 CALL dgemm("T", "N", nco(lshell), nsgf, nso(lshell), 1.0_dp, &
1368 orbtramat(lshell)%c2s, nso(lshell), &
1369 smatrix(isgf, 1), nsgf, 0.0_dp, &
1370 cmatrix(icgf, 1), ncgf)
1371 icgf = icgf + nco(lshell)
1372 isgf = isgf + nso(lshell)
1373 END DO
1374 ELSE
1375 ! assume atom without basis set
1376 ! CPABORT("Unknown basis set type")
1377 END IF
1378 END DO ! iatom
1379
1380 ELSE
1381
1382 orbital_str = "SPHERICAL"
1383
1384 END IF ! print_cartesian
1385
1386 name = trim(energy_str)//", OCCUPATION NUMBERS, AND "// &
1387 trim(orbital_str)//" "//trim(vector_str)
1388
1389 IF (.NOT. my_final) THEN
1390 WRITE (unit=name, fmt="(A,1X,I0)") trim(name)//step_string, scf_step
1391 END IF
1392
1393 ELSE IF (print_occup .OR. print_eigvals) THEN
1394 name = trim(energy_str)//" AND OCCUPATION NUMBERS"
1395
1396 IF (.NOT. my_final) THEN
1397 WRITE (unit=name, fmt="(A,1X,I0)") trim(name)//step_string, scf_step
1398 END IF
1399 END IF ! print_eigvecs
1400
1401 ! Print headline
1402 IF (PRESENT(spin) .AND. (kpoint > 0)) THEN
1403 WRITE (unit=iw, fmt="(/,T2,A,I0)") &
1404 "MO| "//trim(spin)//" "//trim(name)//" FOR K POINT ", kpoint
1405 ELSE IF (PRESENT(spin)) THEN
1406 WRITE (unit=iw, fmt="(/,T2,A)") &
1407 "MO| "//trim(spin)//" "//trim(name)
1408 ELSE IF (kpoint > 0) THEN
1409 WRITE (unit=iw, fmt="(/,T2,A,I0)") &
1410 "MO| "//trim(name)//" FOR K POINT ", kpoint
1411 ELSE
1412 WRITE (unit=iw, fmt="(/,T2,A)") &
1413 "MO| "//trim(name)
1414 END IF
1415
1416 ! Check if only a subset of the MOs has to be printed
1417 IF (all(mo_index_range > 0)) THEN
1418 IF (mo_index_range(2) > nmo) THEN
1419 CALL cp_warn(__location__, &
1420 "The last orbital index is larger than the number of orbitals.")
1421 END IF
1422 IF (mo_index_range(1) > mo_index_range(2)) THEN
1423 CALL cp_warn(__location__, &
1424 "The first orbital index is larger than the last orbital index.")
1425 END IF
1426 first_mo = min(max(1, mo_index_range(1)), nmo)
1427 last_mo = min(max(first_mo, mo_index_range(2)), nmo)
1428 ELSE IF (mo_index_range(2) < 0) THEN
1429 IF (mo_index_range(1) > nmo) THEN
1430 CALL cp_warn(__location__, &
1431 "The first orbital index is larger than the number of orbitals.")
1432 END IF
1433 first_mo = min(max(1, mo_index_range(1)), nmo)
1434 last_mo = nmo
1435 ELSE
1436 first_mo = 1
1437 last_mo = nmo
1438 END IF
1439
1440 IF (print_eigvecs) THEN
1441
1442 ! Print full MO information
1443
1444 DO icol = first_mo, last_mo, ncol
1445
1446 from = icol
1447 to = min((from + ncol - 1), last_mo)
1448
1449 WRITE (unit=iw, fmt="(T2,A)") "MO|"
1450 WRITE (unit=iw, fmt=fmtstr1) &
1451 "MO|", (jcol, jcol=from, to)
1452 WRITE (unit=iw, fmt=fmtstr2) &
1453 "MO|", (mo_eigenvalues(jcol), jcol=from, to)
1454 WRITE (unit=iw, fmt="(T2,A)") "MO|"
1455 WRITE (unit=iw, fmt=fmtstr2) &
1456 "MO|", (mo_occupation_numbers(jcol), jcol=from, to)
1457 WRITE (unit=iw, fmt="(T2,A)") "MO|"
1458
1459 irow = 1
1460
1461 DO iatom = 1, natom
1462
1463 IF (iatom /= 1) WRITE (unit=iw, fmt="(T2,A)") "MO|"
1464
1465 NULLIFY (orb_basis_set, dftb_parameter)
1466 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, &
1467 element_symbol=element_symbol, kind_number=ikind)
1468 CALL get_qs_kind(qs_kind_set(ikind), &
1469 basis_set=orb_basis_set, &
1470 dftb_parameter=dftb_parameter)
1471
1472 IF (print_cartesian) THEN
1473
1474 IF (ASSOCIATED(orb_basis_set)) THEN
1475 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1476 nset=nset, &
1477 nshell=nshell, &
1478 l=l, &
1479 cgf_symbol=bcgf_symbol)
1480
1481 icgf = 1
1482 DO iset = 1, nset
1483 DO ishell = 1, nshell(iset)
1484 lshell = l(ishell, iset)
1485 DO ico = 1, nco(lshell)
1486 WRITE (unit=iw, fmt=fmtstr3) &
1487 "MO|", irow, iatom, adjustr(element_symbol), bcgf_symbol(icgf), &
1488 (cmatrix(irow, jcol), jcol=from, to)
1489 icgf = icgf + 1
1490 irow = irow + 1
1491 END DO
1492 END DO
1493 END DO
1494 ELSE IF (ASSOCIATED(dftb_parameter)) THEN
1495 CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
1496 icgf = 1
1497 DO ishell = 1, lmax + 1
1498 lshell = ishell - 1
1499 DO ico = 1, nco(lshell)
1500 symbol = cgf_symbol(1, indco(1:3, icgf))
1501 symbol(1:2) = " "
1502 WRITE (unit=iw, fmt=fmtstr3) &
1503 "MO|", irow, iatom, adjustr(element_symbol), symbol, &
1504 (cmatrix(irow, jcol), jcol=from, to)
1505 icgf = icgf + 1
1506 irow = irow + 1
1507 END DO
1508 END DO
1509 ELSE
1510 ! assume atom without basis set
1511 ! CPABORT("Unknown basis set type")
1512 END IF
1513
1514 ELSE
1515
1516 IF (ASSOCIATED(orb_basis_set)) THEN
1517 CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
1518 nset=nset, &
1519 nshell=nshell, &
1520 l=l, &
1521 sgf_symbol=bsgf_symbol)
1522 isgf = 1
1523 DO iset = 1, nset
1524 DO ishell = 1, nshell(iset)
1525 lshell = l(ishell, iset)
1526 DO iso = 1, nso(lshell)
1527 WRITE (unit=iw, fmt=fmtstr3) &
1528 "MO|", irow, iatom, adjustr(element_symbol), bsgf_symbol(isgf), &
1529 (smatrix(irow, jcol), jcol=from, to)
1530 isgf = isgf + 1
1531 irow = irow + 1
1532 END DO
1533 END DO
1534 END DO
1535 ELSE IF (ASSOCIATED(dftb_parameter)) THEN
1536 CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
1537 isgf = 1
1538 DO ishell = 1, lmax + 1
1539 lshell = ishell - 1
1540 DO iso = 1, nso(lshell)
1541 symbol = sgf_symbol(1, lshell, -lshell + iso - 1)
1542 symbol(1:2) = " "
1543 WRITE (unit=iw, fmt=fmtstr3) &
1544 "MO|", irow, iatom, adjustr(element_symbol), symbol, &
1545 (smatrix(irow, jcol), jcol=from, to)
1546 isgf = isgf + 1
1547 irow = irow + 1
1548 END DO
1549 END DO
1550 ELSE
1551 ! assume atom without basis set
1552 ! CPABORT("Unknown basis set type")
1553 END IF
1554
1555 END IF ! print_cartesian
1556
1557 END DO ! iatom
1558
1559 END DO ! icol
1560
1561 WRITE (unit=iw, fmt="(T2,A)") "MO|"
1562
1563 ! Release work storage
1564
1565 IF (print_cartesian) THEN
1566 DEALLOCATE (cmatrix)
1567 END IF
1568 DEALLOCATE (smatrix)
1569
1570 ELSE IF (print_occup .OR. print_eigvals) THEN
1571
1572 WRITE (unit=iw, fmt="(T2,A)") "MO|"
1573 fmtstr4 = "(T2,A,I7,3(1X,F22. ))"
1574 WRITE (unit=fmtstr4(19:20), fmt="(I2)") after
1575 IF (my_final .OR. (my_solver_method == "TD")) THEN
1576 WRITE (unit=iw, fmt="(A)") &
1577 " MO| Index Eigenvalue [a.u.] Eigenvalue [eV] Occupation"
1578 ELSE
1579 WRITE (unit=iw, fmt="(A)") &
1580 " MO| Index Energy [a.u.] Energy [eV] Occupation"
1581 END IF
1582 DO imo = first_mo, last_mo
1583 WRITE (unit=iw, fmt=fmtstr4) &
1584 "MO|", imo, mo_eigenvalues(imo), &
1585 mo_eigenvalues(imo)*evolt, &
1586 mo_occupation_numbers(imo)
1587 END DO
1588 fmtstr5 = "(A,T59,F22. )"
1589 WRITE (unit=fmtstr5(12:13), fmt="(I2)") after
1590 WRITE (unit=iw, fmt=fmtstr5) &
1591 " MO| Sum:", accurate_sum(mo_occupation_numbers(:))
1592
1593 END IF ! print_eigvecs
1594
1595 IF (.NOT. my_rtp) THEN
1596 fmtstr6 = "(A,T18,F17. ,A,T41,F17. ,A)"
1597 WRITE (unit=fmtstr6(12:13), fmt="(I2)") after
1598 WRITE (unit=fmtstr6(25:26), fmt="(I2)") after
1599 WRITE (unit=iw, fmt=fmtstr6) &
1600 " MO| E(Fermi):", chemical_potential, " a.u.", chemical_potential*evolt, " eV"
1601 END IF
1602 IF ((homo > 0) .AND. .NOT. my_rtp) THEN
1603 IF ((mo_occupation_numbers(homo) == maxocc) .AND. (last_mo > homo)) THEN
1604 gap = mo_eigenvalues(homo + 1) - &
1605 mo_eigenvalues(homo)
1606 WRITE (unit=iw, fmt=fmtstr6) &
1607 " MO| Band gap:", gap, " a.u.", gap*evolt, " eV"
1608 END IF
1609 END IF
1610 WRITE (unit=iw, fmt="(A)") ""
1611
1612 END IF ! iw
1613
1614 IF (ALLOCATED(mo_eigenvalues)) DEALLOCATE (mo_eigenvalues)
1615 IF (ALLOCATED(mo_occupation_numbers)) DEALLOCATE (mo_occupation_numbers)
1616
1617 CALL cp_print_key_finished_output(iw, logger, dft_section, "PRINT%MO", &
1618 ignore_should_output=should_output)
1619
1620 END SUBROUTINE write_mo_set_to_output_unit
1621
1622END MODULE qs_mo_io
dgemm
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.
Definition grid_cpu_task_list.c:225
cp_dbcsr_api::dbcsr_create
Definition cp_dbcsr_api.F:194
cp_dbcsr_operations::dbcsr_deallocate_matrix_set
Definition cp_dbcsr_operations.F:104
cp_fm_types::cp_fm_to_fm
Definition cp_fm_types.F:85
cp_log_handling::cp_to_string
Definition cp_log_handling.F:90
kahan_sum::accurate_sum
Definition kahan_sum.F:41
qs_density_matrices::calculate_density_matrix
Definition qs_density_matrices.F:58
qs_mo_methods::calculate_subspace_eigenvalues
Definition qs_mo_methods.F:70
qs_mo_occupation::set_mo_occupation
Definition qs_mo_occupation.F:50
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind
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.
Definition atomic_kind_types.F:131
basis_set_types
Definition basis_set_types.F:15
basis_set_types::get_gto_basis_set
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, ccon)
...
Definition basis_set_types.F:643
cp_cfm_types
Represents a complex full matrix distributed on many processors.
Definition cp_cfm_types.F:12
cp_cfm_types::cp_cfm_get_submatrix
subroutine, public cp_cfm_get_submatrix(fm, target_m, start_row, start_col, n_rows, n_cols, transpose)
Extract a sub-matrix from the full matrix: op(target_m)(1:n_rows,1:n_cols) = fm(start_row:start_row+n...
Definition cp_cfm_types.F:357
cp_dbcsr_api
Definition cp_dbcsr_api.F:8
cp_dbcsr_api::dbcsr_binary_write
subroutine, public dbcsr_binary_write(matrix, filepath)
...
Definition cp_dbcsr_api.F:290
cp_dbcsr_api::dbcsr_release
subroutine, public dbcsr_release(matrix)
...
Definition cp_dbcsr_api.F:1142
cp_dbcsr_contrib
Definition cp_dbcsr_contrib.F:8
cp_dbcsr_contrib::dbcsr_checksum
real(kind=dp) function, public dbcsr_checksum(matrix, pos)
Calculates the checksum of a DBCSR matrix.
Definition cp_dbcsr_contrib.F:677
cp_dbcsr_operations
DBCSR operations in CP2K.
Definition cp_dbcsr_operations.F:18
cp_dbcsr_operations::copy_dbcsr_to_fm
subroutine, public copy_dbcsr_to_fm(matrix, fm, plan, check_plan)
Copy a DBCSR matrix to a BLACS matrix.
Definition cp_dbcsr_operations.F:344
cp_dbcsr_operations::copy_fm_to_dbcsr
subroutine, public copy_fm_to_dbcsr(fm, matrix, keep_sparsity)
Copy a BLACS matrix to a dbcsr matrix.
Definition cp_dbcsr_operations.F:240
cp_dbcsr_output
DBCSR output in CP2K.
Definition cp_dbcsr_output.F:17
cp_dbcsr_output::cp_dbcsr_write_sparse_matrix
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, cartesian_basis)
...
Definition cp_dbcsr_output.F:165
cp_files
Utility routines to open and close files.
Definition cp_files.F:22
cp_files::open_file
subroutine, public open_file(file_name, unit_number, file_status, file_form, file_action, file_position, file_pad, file_access, debug)
Opens a file and returns its logical unit number.
Definition cp_files.F:197
cp_files::close_file
subroutine, public close_file(unit_number, file_status)
Closes an open file given by its logical unit number.
Definition cp_files.F:64
cp_fm_types
represent a full matrix distributed on many processors
Definition cp_fm_types.F:15
cp_fm_types::cp_fm_get_info
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
Definition cp_fm_types.F:1046
cp_fm_types::cp_fm_write_unformatted
subroutine, public cp_fm_write_unformatted(fm, unit)
...
Definition cp_fm_types.F:2195
cp_fm_types::cp_fm_set_submatrix
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...
Definition cp_fm_types.F:758
cp_fm_types::cp_fm_set_all
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
Definition cp_fm_types.F:536
cp_fm_types::cp_fm_get_submatrix
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,...
Definition cp_fm_types.F:933
cp_log_handling
various routines to log and control the output. The idea is that decisions about where to log should ...
Definition cp_log_handling.F:41
cp_log_handling::cp_logger_get_default_unit_nr
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
Definition cp_log_handling.F:583
cp_log_handling::cp_get_default_logger
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
Definition cp_log_handling.F:234
cp_output_handling
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
Definition cp_output_handling.F:25
cp_output_handling::cp_print_key_generate_filename
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...
Definition cp_output_handling.F:835
cp_output_handling::cp_print_key_finished_output
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,...
Definition cp_output_handling.F:1156
cp_output_handling::cp_p_file
integer, parameter, public cp_p_file
Definition cp_output_handling.F:107
cp_output_handling::cp_print_key_should_output
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...
Definition cp_output_handling.F:345
cp_output_handling::cp_print_key_unit_nr
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, file_access)
...
Definition cp_output_handling.F:935
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_get_subs_vals
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Definition input_section_types.F:735
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1052
kahan_sum
sums arrays of real/complex numbers with much reduced round-off as compared to a naive implementation...
Definition kahan_sum.F:29
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_string_length
integer, parameter, public default_string_length
Definition kinds.F:57
kinds::default_path_length
integer, parameter, public default_path_length
Definition kinds.F:58
message_passing
Interface to the message passing library MPI.
Definition message_passing.F:27
orbital_pointers
Provides Cartesian and spherical orbital pointers and indices.
Definition orbital_pointers.F:15
orbital_pointers::nco
integer, dimension(:), allocatable, public nco
Definition orbital_pointers.F:42
orbital_pointers::indco
integer, dimension(:, :), allocatable, public indco
Definition orbital_pointers.F:43
orbital_pointers::nso
integer, dimension(:), allocatable, public nso
Definition orbital_pointers.F:42
orbital_symbols
orbital_symbols
Definition orbital_symbols.F:14
orbital_symbols::cgf_symbol
character(len=12) function, public cgf_symbol(n, lxyz)
Build a Cartesian orbital symbol (orbital labels for printing).
Definition orbital_symbols.F:44
orbital_symbols::sgf_symbol
character(len=6) function, public sgf_symbol(n, l, m)
Build a spherical orbital symbol (orbital labels for printing).
Definition orbital_symbols.F:107
orbital_transformation_matrices
Calculation of the spherical harmonics and the corresponding orbital transformation matrices.
Definition orbital_transformation_matrices.F:17
orbital_transformation_matrices::orbtramat
type(orbtramat_type), dimension(:), pointer, public orbtramat
Definition orbital_transformation_matrices.F:47
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
physcon
Definition of physical constants:
Definition physcon.F:68
physcon::evolt
real(kind=dp), parameter, public evolt
Definition physcon.F:183
qs_density_matrices
collects routines that calculate density matrices
Definition qs_density_matrices.F:14
qs_dftb_types
Definition of the DFTB parameter types.
Definition qs_dftb_types.F:12
qs_dftb_utils
Working with the DFTB parameter types.
Definition qs_dftb_utils.F:12
qs_dftb_utils::get_dftb_atom_param
subroutine, public get_dftb_atom_param(dftb_parameter, name, typ, defined, z, zeff, natorb, lmax, skself, occupation, eta, energy, cutoff, xi, di, rcdisp, dudq)
...
Definition qs_dftb_utils.F:124
qs_environment_types
Definition qs_environment_types.F:14
qs_environment_types::get_qs_env
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, mimic, 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, matrix_vhxc, 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, xcint_weights, 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.
Definition qs_environment_types.F:544
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, cneo_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, hund_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, proj_shell_charge, lr_atom, do_mtlr, u_j_loop, ao_coef, init_u_ramping_each_scf, reltmat, ghost, monovalent, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition qs_kind_types.F:483
qs_kind_types::get_qs_kind_set
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, do_mtlr_present, basis_type, total_zeff_corr, npgf_seg, cneo_potential_present, nkind_q, natom_q)
Get attributes of an atomic kind set.
Definition qs_kind_types.F:976
qs_ks_types
Definition qs_ks_types.F:15
qs_mo_io
Definition and initialisation of the mo data type.
Definition qs_mo_io.F:21
qs_mo_io::wfn_restart_file_name
subroutine, public wfn_restart_file_name(filename, exist, section, logger, kp, xas, rtp)
...
Definition qs_mo_io.F:460
qs_mo_io::read_mo_set_from_restart
subroutine, public read_mo_set_from_restart(mo_array, qs_kind_set, particle_set, para_env, id_nr, multiplicity, dft_section, natom_mismatch, cdft, out_unit)
...
Definition qs_mo_io.F:526
qs_mo_io::read_rt_mos_from_restart
subroutine, public read_rt_mos_from_restart(mo_array, rt_mos, qs_kind_set, particle_set, para_env, id_nr, multiplicity, dft_section)
...
Definition qs_mo_io.F:619
qs_mo_io::write_dm_binary_restart
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
Definition qs_mo_io.F:184
qs_mo_io::write_mo_set_low
subroutine, public write_mo_set_low(mo_array, qs_kind_set, particle_set, ires, rt_mos, matrix_ks)
...
Definition qs_mo_io.F:297
qs_mo_io::read_mos_restart_low
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.
Definition qs_mo_io.F:697
qs_mo_io::write_mo_set_to_output_unit
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, qs_env, para_env_inter_kp)
Write MO information to output file (eigenvalues, occupation numbers, coefficients).
Definition qs_mo_io.F:1027
qs_mo_io::write_mo_set_to_restart
subroutine, public write_mo_set_to_restart(mo_array, particle_set, dft_section, qs_kind_set, matrix_ks)
...
Definition qs_mo_io.F:119
qs_mo_io::write_rt_mos_to_restart
subroutine, public write_rt_mos_to_restart(mo_array, rt_mos, particle_set, dft_section, qs_kind_set)
...
Definition qs_mo_io.F:246
qs_mo_methods
collects routines that perform operations directly related to MOs
Definition qs_mo_methods.F:14
qs_mo_occupation
Set occupation of molecular orbitals.
Definition qs_mo_occupation.F:15
qs_mo_types
Definition and initialisation of the mo data type.
Definition qs_mo_types.F:22
qs_mo_types::get_mo_set
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, cmo_coeff)
Get the components of a MO set data structure.
Definition qs_mo_types.F:440
qs_neighbor_list_types
Define the neighbor list data types and the corresponding functionality.
Definition qs_neighbor_list_types.F:17
qs_neighbor_list_types::release_neighbor_list_sets
subroutine, public release_neighbor_list_sets(nlists)
releases an array of neighbor_list_sets
Definition qs_neighbor_list_types.F:1096
qs_neighbor_lists
Generate the atomic neighbor lists.
Definition qs_neighbor_lists.F:19
qs_neighbor_lists::setup_neighbor_list
subroutine, public setup_neighbor_list(ab_list, basis_set_a, basis_set_b, qs_env, mic, symmetric, molecular, operator_type)
Build a neighborlist.
Definition qs_neighbor_lists.F:1534
qs_overlap
Calculation of overlap matrix, its derivatives and forces.
Definition qs_overlap.F:19
qs_overlap::build_overlap_matrix_simple
subroutine, public build_overlap_matrix_simple(ks_env, matrix_s, basis_set_list_a, basis_set_list_b, sab_nl, lcart)
Calculation of the overlap matrix over Cartesian Gaussian functions.
Definition qs_overlap.F:575
basis_set_types::gto_basis_set_p_type
Definition basis_set_types.F:94
basis_set_types::gto_basis_set_type
Definition basis_set_types.F:72
cp_cfm_types::cp_cfm_type
Represent a complex full matrix.
Definition cp_cfm_types.F:70
cp_dbcsr_api::dbcsr_p_type
Definition cp_dbcsr_api.F:172
cp_dbcsr_api::dbcsr_type
Definition cp_dbcsr_api.F:176
cp_fm_types::cp_fm_type
represent a full matrix
Definition cp_fm_types.F:116
cp_log_handling::cp_logger_type
type of a logger, at the moment it contains just a print level starting at which level it should be l...
Definition cp_log_handling.F:140
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
message_passing::mp_para_env_type
stores all the informations relevant to an mpi environment
Definition message_passing.F:774
particle_types::particle_type
Definition particle_types.F:35
qs_dftb_types::qs_dftb_atom_type
Definition qs_dftb_types.F:27
qs_environment_types::qs_environment_type
Definition qs_environment_types.F:223
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:190
qs_ks_types::qs_ks_env_type
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...
Definition qs_ks_types.F:137
qs_mo_types::mo_set_type
Definition qs_mo_types.F:53
qs_neighbor_list_types::neighbor_list_set_p_type
Definition qs_neighbor_list_types.F:67