Line data Source code
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 Restart file for k point calculations
10 : !> \par History
11 : !> \author JGH (30.09.2015)
12 : ! **************************************************************************************************
13 : MODULE kpoint_io
14 :
15 : USE atomic_kind_types, ONLY: get_atomic_kind
16 : USE basis_set_types, ONLY: get_gto_basis_set,&
17 : gto_basis_set_type
18 : USE cp_dbcsr_api, ONLY: dbcsr_get_info,&
19 : dbcsr_p_type
20 : USE cp_dbcsr_operations, ONLY: copy_dbcsr_to_fm,&
21 : copy_fm_to_dbcsr
22 : USE cp_files, ONLY: close_file,&
23 : open_file
24 : USE cp_fm_types, ONLY: cp_fm_read_unformatted,&
25 : cp_fm_type,&
26 : cp_fm_write_unformatted
27 : USE cp_log_handling, ONLY: cp_get_default_logger,&
28 : cp_logger_type,&
29 : cp_to_string
30 : USE cp_output_handling, ONLY: cp_p_file,&
31 : cp_print_key_finished_output,&
32 : cp_print_key_should_output,&
33 : cp_print_key_unit_nr
34 : USE input_section_types, ONLY: section_vals_type
35 : USE kinds, ONLY: default_path_length,&
36 : default_string_length
37 : USE kpoint_types, ONLY: get_kpoint_info,&
38 : kpoint_type
39 : USE message_passing, ONLY: mp_para_env_type
40 : USE orbital_pointers, ONLY: nso
41 : USE particle_types, ONLY: particle_type
42 : USE qs_dftb_types, ONLY: qs_dftb_atom_type
43 : USE qs_dftb_utils, ONLY: get_dftb_atom_param
44 : USE qs_kind_types, ONLY: get_qs_kind,&
45 : qs_kind_type
46 : USE qs_mo_io, ONLY: wfn_restart_file_name
47 : USE qs_scf_types, ONLY: qs_scf_env_type
48 : #include "./base/base_uses.f90"
49 :
50 : IMPLICIT NONE
51 :
52 : PRIVATE
53 :
54 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'kpoint_io'
55 :
56 : INTEGER, PARAMETER :: version = 1
57 :
58 : PUBLIC :: read_kpoints_restart, &
59 : write_kpoints_restart
60 : PUBLIC :: get_cell, write_kpoints_file_header
61 :
62 : ! **************************************************************************************************
63 :
64 : CONTAINS
65 :
66 : ! **************************************************************************************************
67 : !> \brief ...
68 : !> \param denmat ...
69 : !> \param kpoints ...
70 : !> \param scf_env ...
71 : !> \param dft_section ...
72 : !> \param particle_set ...
73 : !> \param qs_kind_set ...
74 : ! **************************************************************************************************
75 12088 : SUBROUTINE write_kpoints_restart(denmat, kpoints, scf_env, dft_section, &
76 : particle_set, qs_kind_set)
77 :
78 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: denmat
79 : TYPE(kpoint_type), POINTER :: kpoints
80 : TYPE(qs_scf_env_type), POINTER :: scf_env
81 : TYPE(section_vals_type), POINTER :: dft_section
82 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
83 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
84 :
85 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_kpoints_restart'
86 : CHARACTER(LEN=30), DIMENSION(2), PARAMETER :: &
87 : keys = ["SCF%PRINT%RESTART_HISTORY", "SCF%PRINT%RESTART "]
88 :
89 : INTEGER :: handle, ic, ikey, ires, ispin, nimages, &
90 : nspin
91 : INTEGER, DIMENSION(3) :: cell
92 12088 : INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
93 : TYPE(cp_fm_type), POINTER :: fmat
94 : TYPE(cp_logger_type), POINTER :: logger
95 :
96 12088 : CALL timeset(routineN, handle)
97 12088 : logger => cp_get_default_logger()
98 :
99 : IF (BTEST(cp_print_key_should_output(logger%iter_info, &
100 12088 : dft_section, keys(1)), cp_p_file) .OR. &
101 : BTEST(cp_print_key_should_output(logger%iter_info, &
102 : dft_section, keys(2)), cp_p_file)) THEN
103 :
104 1134 : DO ikey = 1, SIZE(keys)
105 :
106 756 : IF (BTEST(cp_print_key_should_output(logger%iter_info, &
107 12088 : dft_section, keys(ikey)), cp_p_file)) THEN
108 :
109 : ires = cp_print_key_unit_nr(logger, dft_section, keys(ikey), &
110 : extension=".kp", file_status="REPLACE", file_action="WRITE", &
111 378 : do_backup=.TRUE., file_form="UNFORMATTED")
112 :
113 378 : CALL write_kpoints_file_header(qs_kind_set, particle_set, ires)
114 :
115 378 : nspin = SIZE(denmat, 1)
116 378 : nimages = SIZE(denmat, 2)
117 378 : NULLIFY (cell_to_index)
118 378 : IF (nimages > 1) THEN
119 378 : CALL get_kpoint_info(kpoint=kpoints, cell_to_index=cell_to_index)
120 : END IF
121 378 : CPASSERT(ASSOCIATED(scf_env%scf_work1))
122 378 : fmat => scf_env%scf_work1(1)
123 :
124 880 : DO ispin = 1, nspin
125 502 : IF (ires > 0) WRITE (ires) ispin, nspin, nimages
126 40114 : DO ic = 1, nimages
127 39234 : IF (nimages > 1) THEN
128 39234 : cell = get_cell(ic, cell_to_index)
129 : ELSE
130 0 : cell = 0
131 : END IF
132 39234 : IF (ires > 0) WRITE (ires) ic, cell
133 39234 : CALL copy_dbcsr_to_fm(denmat(ispin, ic)%matrix, fmat)
134 39736 : CALL cp_fm_write_unformatted(fmat, ires)
135 : END DO
136 : END DO
137 :
138 378 : CALL cp_print_key_finished_output(ires, logger, dft_section, TRIM(keys(ikey)))
139 :
140 : END IF
141 :
142 : END DO
143 :
144 : END IF
145 :
146 12088 : CALL timestop(handle)
147 :
148 12088 : END SUBROUTINE write_kpoints_restart
149 :
150 : ! **************************************************************************************************
151 : !> \brief ...
152 : !> \param ic ...
153 : !> \param cell_to_index ...
154 : !> \return ...
155 : ! **************************************************************************************************
156 39806 : FUNCTION get_cell(ic, cell_to_index) RESULT(cell)
157 : INTEGER, INTENT(IN) :: ic
158 : INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
159 : INTEGER, DIMENSION(3) :: cell
160 :
161 : INTEGER :: i1, i2, i3
162 :
163 159224 : cell = 0
164 200746 : ALLCELL: DO i3 = LBOUND(cell_to_index, 3), UBOUND(cell_to_index, 3)
165 1026812 : DO i2 = LBOUND(cell_to_index, 2), UBOUND(cell_to_index, 2)
166 6045298 : DO i1 = LBOUND(cell_to_index, 1), UBOUND(cell_to_index, 1)
167 4557538 : IF (cell_to_index(i1, i2, i3) == ic) THEN
168 39806 : cell(1) = i1
169 39806 : cell(2) = i2
170 39806 : cell(3) = i3
171 39806 : EXIT ALLCELL
172 : END IF
173 : END DO
174 : END DO
175 : END DO ALLCELL
176 :
177 39806 : END FUNCTION get_cell
178 :
179 : ! **************************************************************************************************
180 : !> \brief ...
181 : !> \param qs_kind_set ...
182 : !> \param particle_set ...
183 : !> \param ires ...
184 : !> \param basis_type ...
185 : ! **************************************************************************************************
186 382 : SUBROUTINE write_kpoints_file_header(qs_kind_set, particle_set, ires, basis_type)
187 :
188 : TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
189 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
190 : INTEGER :: ires
191 : CHARACTER(LEN=*), INTENT(IN), OPTIONAL :: basis_type
192 :
193 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_kpoints_file_header'
194 :
195 : CHARACTER(LEN=default_string_length) :: basis_type_local
196 : INTEGER :: handle, iatom, ikind, iset, ishell, &
197 : lmax, lshell, nao, natom, nset, &
198 : nset_max, nsgf, nshell_max
199 382 : INTEGER, DIMENSION(:), POINTER :: nset_info, nshell
200 382 : INTEGER, DIMENSION(:, :), POINTER :: l, nshell_info
201 382 : INTEGER, DIMENSION(:, :, :), POINTER :: nso_info
202 : TYPE(gto_basis_set_type), POINTER :: orb_basis_set
203 : TYPE(qs_dftb_atom_type), POINTER :: dftb_parameter
204 :
205 382 : CALL timeset(routineN, handle)
206 :
207 382 : IF (ires > 0) THEN
208 :
209 191 : IF (PRESENT(basis_type)) THEN
210 2 : basis_type_local = basis_type
211 : ELSE
212 189 : basis_type_local = "ORB"
213 : END IF
214 :
215 : ! create some info about the basis set first
216 191 : natom = SIZE(particle_set, 1)
217 191 : nset_max = 0
218 191 : nshell_max = 0
219 :
220 914 : DO iatom = 1, natom
221 723 : NULLIFY (orb_basis_set, dftb_parameter)
222 723 : CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
223 : CALL get_qs_kind(qs_kind_set(ikind), basis_type=basis_type_local, &
224 723 : basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
225 1637 : IF (ASSOCIATED(orb_basis_set)) THEN
226 : CALL get_gto_basis_set(gto_basis_set=orb_basis_set, &
227 : nset=nset, &
228 : nshell=nshell, &
229 587 : l=l)
230 587 : nset_max = MAX(nset_max, nset)
231 1868 : DO iset = 1, nset
232 1868 : nshell_max = MAX(nshell_max, nshell(iset))
233 : END DO
234 136 : ELSEIF (ASSOCIATED(dftb_parameter)) THEN
235 136 : CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
236 136 : nset_max = MAX(nset_max, 1)
237 136 : nshell_max = MAX(nshell_max, lmax + 1)
238 : ELSE
239 0 : CPABORT("Unknown basis type.")
240 : END IF
241 : END DO
242 :
243 955 : ALLOCATE (nso_info(nshell_max, nset_max, natom))
244 4262 : nso_info(:, :, :) = 0
245 :
246 764 : ALLOCATE (nshell_info(nset_max, natom))
247 2340 : nshell_info(:, :) = 0
248 :
249 573 : ALLOCATE (nset_info(natom))
250 914 : nset_info(:) = 0
251 :
252 191 : nao = 0
253 914 : DO iatom = 1, natom
254 723 : NULLIFY (orb_basis_set, dftb_parameter)
255 723 : CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
256 : CALL get_qs_kind(qs_kind_set(ikind), basis_type=basis_type_local, &
257 723 : basis_set=orb_basis_set, dftb_parameter=dftb_parameter)
258 1637 : IF (ASSOCIATED(orb_basis_set)) THEN
259 587 : CALL get_gto_basis_set(gto_basis_set=orb_basis_set, nsgf=nsgf, nset=nset, nshell=nshell, l=l)
260 587 : nao = nao + nsgf
261 587 : nset_info(iatom) = nset
262 1868 : DO iset = 1, nset
263 1281 : nshell_info(iset, iatom) = nshell(iset)
264 3323 : DO ishell = 1, nshell(iset)
265 1455 : lshell = l(ishell, iset)
266 2736 : nso_info(ishell, iset, iatom) = nso(lshell)
267 : END DO
268 : END DO
269 136 : ELSEIF (ASSOCIATED(dftb_parameter)) THEN
270 136 : CALL get_dftb_atom_param(dftb_parameter, lmax=lmax)
271 136 : nset_info(iatom) = 1
272 136 : nshell_info(1, iatom) = lmax + 1
273 408 : DO ishell = 1, lmax + 1
274 272 : lshell = ishell - 1
275 408 : nso_info(ishell, 1, iatom) = nso(lshell)
276 : END DO
277 136 : nao = nao + (lmax + 1)**2
278 : ELSE
279 0 : CPABORT("Unknown basis set type.")
280 : END IF
281 : END DO
282 :
283 191 : WRITE (ires) version
284 191 : WRITE (ires) natom, nao, nset_max, nshell_max
285 914 : WRITE (ires) nset_info
286 2340 : WRITE (ires) nshell_info
287 4262 : WRITE (ires) nso_info
288 :
289 191 : DEALLOCATE (nset_info, nshell_info, nso_info)
290 : END IF
291 :
292 382 : CALL timestop(handle)
293 :
294 382 : END SUBROUTINE write_kpoints_file_header
295 :
296 : ! **************************************************************************************************
297 : !> \brief ...
298 : !> \param denmat ...
299 : !> \param kpoints ...
300 : !> \param fmwork ...
301 : !> \param natom ...
302 : !> \param para_env ...
303 : !> \param id_nr ...
304 : !> \param dft_section ...
305 : !> \param natom_mismatch ...
306 : ! **************************************************************************************************
307 24 : SUBROUTINE read_kpoints_restart(denmat, kpoints, fmwork, natom, &
308 : para_env, id_nr, dft_section, natom_mismatch)
309 :
310 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: denmat
311 : TYPE(kpoint_type), POINTER :: kpoints
312 : TYPE(cp_fm_type), DIMENSION(:), INTENT(INOUT) :: fmwork
313 : INTEGER, INTENT(IN) :: natom
314 : TYPE(mp_para_env_type), POINTER :: para_env
315 : INTEGER, INTENT(IN) :: id_nr
316 : TYPE(section_vals_type), POINTER :: dft_section
317 : LOGICAL, INTENT(OUT) :: natom_mismatch
318 :
319 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_kpoints_restart'
320 :
321 : CHARACTER(LEN=default_path_length) :: file_name
322 : INTEGER :: handle, restart_unit
323 : LOGICAL :: exist
324 : TYPE(cp_logger_type), POINTER :: logger
325 :
326 12 : CALL timeset(routineN, handle)
327 12 : logger => cp_get_default_logger()
328 :
329 12 : restart_unit = -1
330 :
331 12 : IF (para_env%is_source()) THEN
332 :
333 6 : CALL wfn_restart_file_name(file_name, exist, dft_section, logger, kp=.TRUE.)
334 6 : IF (id_nr /= 0) THEN
335 : ! Is it one of the backup files?
336 0 : file_name = TRIM(file_name)//".bak-"//ADJUSTL(cp_to_string(id_nr))
337 : END IF
338 :
339 : CALL open_file(file_name=file_name, &
340 : file_action="READ", &
341 : file_form="UNFORMATTED", &
342 : file_status="OLD", &
343 6 : unit_number=restart_unit)
344 :
345 : END IF
346 :
347 : CALL read_kpoints_restart_low(denmat, kpoints, fmwork(1), para_env, &
348 12 : natom, restart_unit, natom_mismatch)
349 :
350 : ! only the io_node returns natom_mismatch, must broadcast it
351 12 : CALL para_env%bcast(natom_mismatch)
352 :
353 : ! Close restart file
354 12 : IF (para_env%is_source()) CALL close_file(unit_number=restart_unit)
355 :
356 12 : CALL timestop(handle)
357 :
358 12 : END SUBROUTINE read_kpoints_restart
359 :
360 : ! **************************************************************************************************
361 : !> \brief Reading the mos from apreviously defined restart file
362 : !> \param denmat ...
363 : !> \param kpoints ...
364 : !> \param fmat ...
365 : !> \param para_env ...
366 : !> \param natom ...
367 : !> \param rst_unit ...
368 : !> \param natom_mismatch ...
369 : !> \par History
370 : !> 12.2007 created [MI]
371 : !> \author MI
372 : ! **************************************************************************************************
373 24 : SUBROUTINE read_kpoints_restart_low(denmat, kpoints, fmat, para_env, &
374 : natom, rst_unit, natom_mismatch)
375 :
376 : TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: denmat
377 : TYPE(kpoint_type), POINTER :: kpoints
378 : TYPE(cp_fm_type), INTENT(INOUT) :: fmat
379 : TYPE(mp_para_env_type), POINTER :: para_env
380 : INTEGER, INTENT(IN) :: natom, rst_unit
381 : LOGICAL, INTENT(OUT) :: natom_mismatch
382 :
383 : INTEGER :: ic, image, ispin, nao, nao_read, &
384 : natom_read, ni, nimages, nset_max, &
385 : nshell_max, nspin, nspin_read, &
386 : version_read
387 : INTEGER, DIMENSION(3) :: cell
388 12 : INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
389 : LOGICAL :: natom_match
390 :
391 0 : CPASSERT(ASSOCIATED(denmat))
392 12 : CALL dbcsr_get_info(denmat(1, 1)%matrix, nfullrows_total=nao)
393 :
394 12 : nspin = SIZE(denmat, 1)
395 12 : nimages = SIZE(denmat, 2)
396 12 : NULLIFY (cell_to_index)
397 12 : IF (nimages > 1) THEN
398 12 : CALL get_kpoint_info(kpoint=kpoints, cell_to_index=cell_to_index)
399 : END IF
400 :
401 12 : IF (para_env%is_source()) THEN
402 6 : READ (rst_unit) version_read
403 6 : IF (version_read /= version) THEN
404 : CALL cp_abort(__LOCATION__, &
405 0 : " READ RESTART : Version of restart file not supported")
406 : END IF
407 6 : READ (rst_unit) natom_read, nao_read, nset_max, nshell_max
408 6 : natom_match = (natom_read == natom)
409 6 : IF (natom_match) THEN
410 6 : IF (nao_read /= nao) THEN
411 : CALL cp_abort(__LOCATION__, &
412 0 : " READ RESTART : Different number of basis functions")
413 : END IF
414 6 : READ (rst_unit)
415 6 : READ (rst_unit)
416 6 : READ (rst_unit)
417 : END IF
418 : END IF
419 :
420 : ! make natom_match and natom_mismatch uniform across all nodes
421 12 : CALL para_env%bcast(natom_match)
422 12 : natom_mismatch = .NOT. natom_match
423 : ! handle natom_match false
424 12 : IF (natom_mismatch) THEN
425 : CALL cp_warn(__LOCATION__, &
426 0 : " READ RESTART : WARNING : DIFFERENT natom, returning ")
427 : ELSE
428 :
429 : DO
430 12 : ispin = 0
431 12 : IF (para_env%is_source()) THEN
432 6 : READ (rst_unit) ispin, nspin_read, ni
433 : END IF
434 12 : CALL para_env%bcast(ispin)
435 12 : CALL para_env%bcast(nspin_read)
436 12 : CALL para_env%bcast(ni)
437 12 : IF (nspin_read /= nspin) THEN
438 : CALL cp_abort(__LOCATION__, &
439 0 : " READ RESTART : Different spin polarisation not supported")
440 : END IF
441 1608 : DO ic = 1, ni
442 1596 : IF (para_env%is_source()) THEN
443 798 : READ (rst_unit) image, cell
444 : END IF
445 1596 : CALL para_env%bcast(image)
446 1596 : CALL para_env%bcast(cell)
447 : !
448 1596 : CALL cp_fm_read_unformatted(fmat, rst_unit)
449 : !
450 1596 : IF (nimages > 1) THEN
451 1596 : image = get_index(cell, cell_to_index)
452 0 : ELSE IF (SUM(ABS(cell(:))) == 0) THEN
453 0 : image = 1
454 : ELSE
455 0 : image = 0
456 : END IF
457 1608 : IF (image >= 1 .AND. image <= nimages) THEN
458 1596 : CALL copy_fm_to_dbcsr(fmat, denmat(ispin, image)%matrix)
459 : END IF
460 : END DO
461 12 : IF (ispin == nspin_read) EXIT
462 : END DO
463 :
464 : END IF
465 :
466 12 : END SUBROUTINE read_kpoints_restart_low
467 :
468 : ! **************************************************************************************************
469 : !> \brief ...
470 : !> \param cell ...
471 : !> \param cell_to_index ...
472 : !> \return ...
473 : ! **************************************************************************************************
474 1596 : FUNCTION get_index(cell, cell_to_index) RESULT(index)
475 : INTEGER, DIMENSION(3), INTENT(IN) :: cell
476 : INTEGER, DIMENSION(:, :, :), POINTER :: cell_to_index
477 : INTEGER :: index
478 :
479 : IF (CELL(1) >= LBOUND(cell_to_index, 1) .AND. CELL(1) <= UBOUND(cell_to_index, 1) .AND. &
480 : CELL(2) >= LBOUND(cell_to_index, 2) .AND. CELL(2) <= UBOUND(cell_to_index, 2) .AND. &
481 11172 : CELL(3) >= LBOUND(cell_to_index, 3) .AND. CELL(3) <= UBOUND(cell_to_index, 3)) THEN
482 :
483 1596 : index = cell_to_index(cell(1), cell(2), cell(3))
484 :
485 : ELSE
486 :
487 : index = 0
488 :
489 : END IF
490 :
491 1596 : END FUNCTION get_index
492 :
493 : END MODULE kpoint_io
|
