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