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