(git:34ef472)
pao_io.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief Routines for reading and writing restart files.
10 !> \author Ole Schuett
11 ! **************************************************************************************************
12 MODULE pao_io
13  USE atomic_kind_types, ONLY: atomic_kind_type,&
14  get_atomic_kind
15  USE basis_set_types, ONLY: gto_basis_set_type
16  USE cell_types, ONLY: cell_type
17  USE cp_files, ONLY: close_file,&
18  open_file
19  USE cp_log_handling, ONLY: cp_get_default_logger,&
20  cp_logger_get_default_io_unit,&
21  cp_logger_type
22  USE cp_output_handling, ONLY: cp_p_file,&
23  cp_print_key_finished_output,&
24  cp_print_key_should_output,&
25  cp_print_key_unit_nr
26  USE dbcsr_api, ONLY: &
27  dbcsr_convert_dbcsr_to_csr, dbcsr_copy, dbcsr_csr_create_from_dbcsr, &
28  dbcsr_csr_dbcsr_blkrow_dist, dbcsr_csr_destroy, dbcsr_csr_type, dbcsr_csr_write, &
29  dbcsr_desymmetrize, dbcsr_get_block_p, dbcsr_get_info, dbcsr_has_symmetry, dbcsr_release, &
30  dbcsr_type
31  USE dm_ls_scf_types, ONLY: ls_scf_env_type
32  USE input_section_types, ONLY: section_vals_get_subs_vals,&
33  section_vals_type,&
34  section_vals_val_get
35  USE kinds, ONLY: default_path_length,&
36  default_string_length,&
37  dp
38  USE message_passing, ONLY: mp_para_env_type
39  USE pao_input, ONLY: id2str
40  USE pao_param, ONLY: pao_param_count
41  USE pao_types, ONLY: pao_env_type
42  USE particle_types, ONLY: particle_type
43  USE physcon, ONLY: angstrom
44  USE qs_environment_types, ONLY: get_qs_env,&
45  qs_environment_type
46  USE qs_kind_types, ONLY: get_qs_kind,&
47  pao_potential_type,&
48  qs_kind_type
49 #include "./base/base_uses.f90"
50 
51  IMPLICIT NONE
52 
53  PRIVATE
54 
55  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pao_io'
56 
57  PUBLIC :: pao_read_restart, pao_write_restart
58  PUBLIC :: pao_read_raw, pao_kinds_ensure_equal
59  PUBLIC :: pao_ioblock_type, pao_iokind_type
60  PUBLIC :: pao_write_ks_matrix_csr, pao_write_s_matrix_csr
61 
62  ! data types used by pao_read_raw()
63  TYPE pao_ioblock_type
64  REAL(dp), DIMENSION(:, :), ALLOCATABLE :: p
65  END TYPE pao_ioblock_type
66 
67  TYPE pao_iokind_type
68  CHARACTER(LEN=default_string_length) :: name = ""
69  INTEGER :: z = -1
70  CHARACTER(LEN=default_string_length) :: prim_basis_name = ""
71  INTEGER :: prim_basis_size = -1
72  INTEGER :: pao_basis_size = -1
73  INTEGER :: nparams = -1
74  TYPE(pao_potential_type), ALLOCATABLE, DIMENSION(:) :: pao_potentials
75  END TYPE pao_iokind_type
76 
77  INTEGER, PARAMETER, PRIVATE :: file_format_version = 4
78 
79 CONTAINS
80 
81 ! **************************************************************************************************
82 !> \brief Reads restart file
83 !> \param pao ...
84 !> \param qs_env ...
85 ! **************************************************************************************************
86  SUBROUTINE pao_read_restart(pao, qs_env)
87  TYPE(pao_env_type), POINTER :: pao
88  TYPE(qs_environment_type), POINTER :: qs_env
89 
90  CHARACTER(LEN=default_string_length) :: param
91  INTEGER :: iatom, ikind, natoms
92  INTEGER, ALLOCATABLE, DIMENSION(:) :: atom2kind
93  INTEGER, DIMENSION(:), POINTER :: col_blk_sizes, row_blk_sizes
94  LOGICAL :: found
95  REAL(dp) :: diff
96  REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: hmat, positions
97  REAL(dp), DIMENSION(:, :), POINTER :: block_x, buffer
98  TYPE(cell_type), POINTER :: cell
99  TYPE(mp_para_env_type), POINTER :: para_env
100  TYPE(pao_ioblock_type), ALLOCATABLE, DIMENSION(:) :: xblocks
101  TYPE(pao_iokind_type), ALLOCATABLE, DIMENSION(:) :: kinds
102  TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
103 
104  cpassert(len_trim(pao%restart_file) > 0)
105  IF (pao%iw > 0) WRITE (pao%iw, '(A,A)') " PAO| Reading matrix_X from restart file: ", trim(pao%restart_file)
106 
107  CALL get_qs_env(qs_env, &
108  para_env=para_env, &
109  natom=natoms, &
110  cell=cell, &
111  particle_set=particle_set)
112 
113  ! read and check restart file on first rank only
114  IF (para_env%is_source()) THEN
115  CALL pao_read_raw(pao%restart_file, param, hmat, kinds, atom2kind, positions, xblocks)
116 
117  ! check cell
118  IF (maxval(abs(hmat - cell%hmat)) > 1e-10) &
119  cpwarn("Restarting from different cell")
120 
121  ! check parametrization
122  IF (trim(param) .NE. trim(adjustl(id2str(pao%parameterization)))) &
123  cpabort("Restart PAO parametrization does not match")
124 
125  ! check kinds
126  DO ikind = 1, SIZE(kinds)
127  CALL pao_kinds_ensure_equal(pao, qs_env, ikind, kinds(ikind))
128  END DO
129 
130  ! check number of atoms
131  IF (SIZE(positions, 1) /= natoms) &
132  cpabort("Number of atoms do not match")
133 
134  ! check atom2kind
135  DO iatom = 1, natoms
136  IF (atom2kind(iatom) /= particle_set(iatom)%atomic_kind%kind_number) &
137  cpabort("Restart atomic kinds do not match.")
138  END DO
139 
140  ! check positions, warning only
141  diff = 0.0_dp
142  DO iatom = 1, natoms
143  diff = max(diff, maxval(abs(positions(iatom, :) - particle_set(iatom)%r)))
144  END DO
145  IF (diff > 1e-10) &
146  cpwarn("Restarting from different atom positions")
147 
148  END IF
149 
150  ! scatter xblocks across ranks to fill pao%matrix_X
151  ! this could probably be done more efficiently
152  CALL dbcsr_get_info(pao%matrix_X, row_blk_size=row_blk_sizes, col_blk_size=col_blk_sizes)
153  DO iatom = 1, natoms
154  ALLOCATE (buffer(row_blk_sizes(iatom), col_blk_sizes(iatom)))
155  IF (para_env%is_source()) THEN
156  cpassert(row_blk_sizes(iatom) == SIZE(xblocks(iatom)%p, 1))
157  cpassert(col_blk_sizes(iatom) == SIZE(xblocks(iatom)%p, 2))
158  buffer = xblocks(iatom)%p
159  END IF
160  CALL para_env%bcast(buffer)
161  CALL dbcsr_get_block_p(matrix=pao%matrix_X, row=iatom, col=iatom, block=block_x, found=found)
162  IF (ASSOCIATED(block_x)) &
163  block_x = buffer
164  DEALLOCATE (buffer)
165  END DO
166 
167  ! ALLOCATABLEs deallocate themselves
168 
169  END SUBROUTINE pao_read_restart
170 
171 ! **************************************************************************************************
172 !> \brief Reads a restart file into temporary datastructures
173 !> \param filename ...
174 !> \param param ...
175 !> \param hmat ...
176 !> \param kinds ...
177 !> \param atom2kind ...
178 !> \param positions ...
179 !> \param xblocks ...
180 !> \param ml_range ...
181 ! **************************************************************************************************
182  SUBROUTINE pao_read_raw(filename, param, hmat, kinds, atom2kind, positions, xblocks, ml_range)
183  CHARACTER(LEN=default_path_length), INTENT(IN) :: filename
184  CHARACTER(LEN=default_string_length), INTENT(OUT) :: param
185  REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: hmat
186  TYPE(pao_iokind_type), ALLOCATABLE, DIMENSION(:) :: kinds
187  INTEGER, ALLOCATABLE, DIMENSION(:) :: atom2kind
188  REAL(dp), ALLOCATABLE, DIMENSION(:, :) :: positions
189  TYPE(pao_ioblock_type), ALLOCATABLE, DIMENSION(:) :: xblocks
190  INTEGER, DIMENSION(2), INTENT(OUT), OPTIONAL :: ml_range
191 
192  CHARACTER(LEN=default_string_length) :: label, str_in
193  INTEGER :: i1, i2, iatom, ikind, ipot, natoms, &
194  nkinds, nparams, unit_nr, xblocks_read
195  REAL(dp) :: r1, r2
196  REAL(dp), DIMENSION(3) :: pos_in
197  REAL(dp), DIMENSION(3, 3) :: hmat_angstrom
198 
199  cpassert(.NOT. ALLOCATED(hmat))
200  cpassert(.NOT. ALLOCATED(kinds))
201  cpassert(.NOT. ALLOCATED(atom2kind))
202  cpassert(.NOT. ALLOCATED(positions))
203  cpassert(.NOT. ALLOCATED(xblocks))
204 
205  natoms = -1
206  nkinds = -1
207  xblocks_read = 0
208 
209  CALL open_file(file_name=filename, file_status="OLD", file_form="FORMATTED", &
210  file_action="READ", unit_number=unit_nr)
211 
212  ! check if file starts with proper header !TODO: introduce a more unique header
213  READ (unit_nr, fmt=*) label, i1
214  IF (trim(label) /= "Version") &
215  cpabort("PAO restart file appears to be corrupted.")
216  IF (i1 /= file_format_version) cpabort("Restart PAO file format version is wrong")
217 
218  DO WHILE (.true.)
219  READ (unit_nr, fmt=*) label
220  backspace(unit_nr)
221 
222  IF (trim(label) == "Parametrization") THEN
223  READ (unit_nr, fmt=*) label, str_in
224  param = str_in
225 
226  ELSE IF (trim(label) == "Cell") THEN
227  READ (unit_nr, fmt=*) label, hmat_angstrom
228  ALLOCATE (hmat(3, 3))
229  hmat(:, :) = hmat_angstrom(:, :)/angstrom
230 
231  ELSE IF (trim(label) == "Nkinds") THEN
232  READ (unit_nr, fmt=*) label, nkinds
233  ALLOCATE (kinds(nkinds))
234 
235  ELSE IF (trim(label) == "Kind") THEN
236  READ (unit_nr, fmt=*) label, ikind, str_in, i1
237  cpassert(ALLOCATED(kinds))
238  kinds(ikind)%name = str_in
239  kinds(ikind)%z = i1
240 
241  ELSE IF (trim(label) == "PrimBasis") THEN
242  READ (unit_nr, fmt=*) label, ikind, i1, str_in
243  cpassert(ALLOCATED(kinds))
244  kinds(ikind)%prim_basis_size = i1
245  kinds(ikind)%prim_basis_name = str_in
246 
247  ELSE IF (trim(label) == "PaoBasis") THEN
248  READ (unit_nr, fmt=*) label, ikind, i1
249  cpassert(ALLOCATED(kinds))
250  kinds(ikind)%pao_basis_size = i1
251 
252  ELSE IF (trim(label) == "NPaoPotentials") THEN
253  READ (unit_nr, fmt=*) label, ikind, i1
254  cpassert(ALLOCATED(kinds))
255  ALLOCATE (kinds(ikind)%pao_potentials(i1))
256 
257  ELSE IF (trim(label) == "PaoPotential") THEN
258  READ (unit_nr, fmt=*) label, ikind, ipot, i1, i2, r1, r2
259  cpassert(ALLOCATED(kinds(ikind)%pao_potentials))
260  kinds(ikind)%pao_potentials(ipot)%maxl = i1
261  kinds(ikind)%pao_potentials(ipot)%max_projector = i2
262  kinds(ikind)%pao_potentials(ipot)%beta = r1
263  kinds(ikind)%pao_potentials(ipot)%weight = r2
264 
265  ELSE IF (trim(label) == "NParams") THEN
266  READ (unit_nr, fmt=*) label, ikind, i1
267  cpassert(ALLOCATED(kinds))
268  kinds(ikind)%nparams = i1
269 
270  ELSE IF (trim(label) == "Natoms") THEN
271  READ (unit_nr, fmt=*) label, natoms
272  ALLOCATE (positions(natoms, 3), atom2kind(natoms), xblocks(natoms))
273  positions = 0.0_dp; atom2kind = -1
274  IF (PRESENT(ml_range)) ml_range = (/1, natoms/)
275 
276  ELSE IF (trim(label) == "MLRange") THEN
277  ! Natoms entry has to come first
278  cpassert(natoms > 0)
279  ! range of atoms whose xblocks are used for machine learning
280  READ (unit_nr, fmt=*) label, i1, i2
281  IF (PRESENT(ml_range)) ml_range = (/i1, i2/)
282 
283  ELSE IF (trim(label) == "Atom") THEN
284  READ (unit_nr, fmt=*) label, iatom, str_in, pos_in
285  cpassert(ALLOCATED(kinds))
286  DO ikind = 1, nkinds
287  IF (trim(kinds(ikind)%name) .EQ. trim(str_in)) EXIT
288  END DO
289  cpassert(ALLOCATED(atom2kind) .AND. ALLOCATED(positions))
290  atom2kind(iatom) = ikind
291  positions(iatom, :) = pos_in/angstrom
292 
293  ELSE IF (trim(label) == "Xblock") THEN
294  READ (unit_nr, fmt=*) label, iatom
295  cpassert(ALLOCATED(kinds) .AND. ALLOCATED(atom2kind))
296  ikind = atom2kind(iatom)
297  nparams = kinds(ikind)%nparams
298  cpassert(nparams >= 0)
299  ALLOCATE (xblocks(iatom)%p(nparams, 1))
300  backspace(unit_nr)
301  READ (unit_nr, fmt=*) label, iatom, xblocks(iatom)%p
302  xblocks_read = xblocks_read + 1
303  cpassert(iatom == xblocks_read) ! ensure blocks are read in order
304 
305  ELSE IF (trim(label) == "THE_END") THEN
306  EXIT
307  ELSE
308  !CPWARN("Skipping restart header with label: "//TRIM(label))
309  READ (unit_nr, fmt=*) label ! just read again and ignore
310  END IF
311  END DO
312  CALL close_file(unit_number=unit_nr)
313 
314  cpassert(xblocks_read == natoms) ! ensure we read all blocks
315 
316  END SUBROUTINE pao_read_raw
317 
318 ! **************************************************************************************************
319 !> \brief Ensure that the kind read from the restart is equal to the kind curretly in use.
320 !> \param pao ...
321 !> \param qs_env ...
322 !> \param ikind ...
323 !> \param pao_kind ...
324 ! **************************************************************************************************
325  SUBROUTINE pao_kinds_ensure_equal(pao, qs_env, ikind, pao_kind)
326  TYPE(pao_env_type), POINTER :: pao
327  TYPE(qs_environment_type), POINTER :: qs_env
328  INTEGER, INTENT(IN) :: ikind
329  TYPE(pao_iokind_type), INTENT(IN) :: pao_kind
330 
331  CHARACTER(LEN=default_string_length) :: name
332  INTEGER :: ipot, nparams, pao_basis_size, z
333  TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
334  TYPE(gto_basis_set_type), POINTER :: basis_set
335  TYPE(pao_potential_type), DIMENSION(:), POINTER :: pao_potentials
336  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
337 
338  CALL get_qs_env(qs_env, &
339  atomic_kind_set=atomic_kind_set, &
340  qs_kind_set=qs_kind_set)
341 
342  IF (ikind > SIZE(atomic_kind_set) .OR. ikind > SIZE(qs_kind_set)) &
343  cpabort("Some kinds are missing.")
344 
345  CALL get_atomic_kind(atomic_kind_set(ikind), z=z, name=name)
346  CALL get_qs_kind(qs_kind_set(ikind), &
347  basis_set=basis_set, &
348  pao_basis_size=pao_basis_size, &
349  pao_potentials=pao_potentials)
350  CALL pao_param_count(pao, qs_env, ikind=ikind, nparams=nparams)
351 
352  IF (pao_kind%nparams /= nparams) &
353  cpabort("Number of parameters do not match")
354  IF (trim(pao_kind%name) .NE. trim(name)) &
355  cpabort("Kind names do not match")
356  IF (pao_kind%z /= z) &
357  cpabort("Atomic numbers do not match")
358  IF (trim(pao_kind%prim_basis_name) .NE. trim(basis_set%name)) &
359  cpabort("Primary Basis-set name does not match")
360  IF (pao_kind%prim_basis_size /= basis_set%nsgf) &
361  cpabort("Primary Basis-set size does not match")
362  IF (pao_kind%pao_basis_size /= pao_basis_size) &
363  cpabort("PAO basis size does not match")
364  IF (SIZE(pao_kind%pao_potentials) /= SIZE(pao_potentials)) &
365  cpabort("Number of PAO_POTENTIALS does not match")
366 
367  DO ipot = 1, SIZE(pao_potentials)
368  IF (pao_kind%pao_potentials(ipot)%maxl /= pao_potentials(ipot)%maxl) &
369  cpabort("PAO_POT_MAXL does not match")
370  IF (pao_kind%pao_potentials(ipot)%max_projector /= pao_potentials(ipot)%max_projector) &
371  cpabort("PAO_POT_MAX_PROJECTOR does not match")
372  IF (pao_kind%pao_potentials(ipot)%beta /= pao_potentials(ipot)%beta) &
373  cpwarn("PAO_POT_BETA does not match")
374  IF (pao_kind%pao_potentials(ipot)%weight /= pao_potentials(ipot)%weight) &
375  cpwarn("PAO_POT_WEIGHT does not match")
376  END DO
377 
378  END SUBROUTINE pao_kinds_ensure_equal
379 
380 ! **************************************************************************************************
381 !> \brief Writes restart file
382 !> \param pao ...
383 !> \param qs_env ...
384 !> \param energy ...
385 ! **************************************************************************************************
386  SUBROUTINE pao_write_restart(pao, qs_env, energy)
387  TYPE(pao_env_type), POINTER :: pao
388  TYPE(qs_environment_type), POINTER :: qs_env
389  REAL(dp) :: energy
390 
391  CHARACTER(len=*), PARAMETER :: printkey_section = 'DFT%LS_SCF%PAO%PRINT%RESTART', &
392  routinen = 'pao_write_restart'
393 
394  INTEGER :: handle, unit_max, unit_nr
395  TYPE(cp_logger_type), POINTER :: logger
396  TYPE(mp_para_env_type), POINTER :: para_env
397  TYPE(section_vals_type), POINTER :: input
398 
399  CALL timeset(routinen, handle)
400  logger => cp_get_default_logger()
401 
402  CALL get_qs_env(qs_env, input=input, para_env=para_env)
403 
404  ! open file
405  unit_nr = cp_print_key_unit_nr(logger, &
406  input, &
407  printkey_section, &
408  extension=".pao", &
409  file_action="WRITE", &
410  file_position="REWIND", &
411  file_status="UNKNOWN", &
412  do_backup=.true.)
413 
414  ! although just rank-0 writes the trajectory it requires collective MPI calls
415  unit_max = unit_nr
416  CALL para_env%max(unit_max)
417  IF (unit_max > 0) THEN
418  IF (pao%iw > 0) WRITE (pao%iw, '(A,A)') " PAO| Writing restart file."
419  IF (unit_nr > 0) &
420  CALL write_restart_header(pao, qs_env, energy, unit_nr)
421 
422  CALL pao_write_diagonal_blocks(para_env, pao%matrix_X, "Xblock", unit_nr)
423 
424  END IF
425 
426  ! close file
427  IF (unit_nr > 0) WRITE (unit_nr, '(A)') "THE_END"
428  CALL cp_print_key_finished_output(unit_nr, logger, input, printkey_section)
429 
430  CALL timestop(handle)
431  END SUBROUTINE pao_write_restart
432 
433 ! **************************************************************************************************
434 !> \brief Write the digonal blocks of given DBCSR matrix into the provided unit_nr
435 !> \param para_env ...
436 !> \param matrix ...
437 !> \param label ...
438 !> \param unit_nr ...
439 ! **************************************************************************************************
440  SUBROUTINE pao_write_diagonal_blocks(para_env, matrix, label, unit_nr)
441  TYPE(mp_para_env_type), POINTER :: para_env
442  TYPE(dbcsr_type) :: matrix
443  CHARACTER(LEN=*), INTENT(IN) :: label
444  INTEGER, INTENT(IN) :: unit_nr
445 
446  INTEGER :: iatom, natoms
447  INTEGER, DIMENSION(:), POINTER :: col_blk_sizes, row_blk_sizes
448  LOGICAL :: found
449  REAL(dp), DIMENSION(:, :), POINTER :: local_block, mpi_buffer
450 
451  !TODO: this is a serial algorithm
452  CALL dbcsr_get_info(matrix, row_blk_size=row_blk_sizes, col_blk_size=col_blk_sizes)
453  cpassert(SIZE(row_blk_sizes) == SIZE(col_blk_sizes))
454  natoms = SIZE(row_blk_sizes)
455 
456  DO iatom = 1, natoms
457  ALLOCATE (mpi_buffer(row_blk_sizes(iatom), col_blk_sizes(iatom)))
458  NULLIFY (local_block)
459  CALL dbcsr_get_block_p(matrix=matrix, row=iatom, col=iatom, block=local_block, found=found)
460  IF (ASSOCIATED(local_block)) THEN
461  IF (SIZE(local_block) > 0) & ! catch corner-case
462  mpi_buffer(:, :) = local_block(:, :)
463  ELSE
464  mpi_buffer(:, :) = 0.0_dp
465  END IF
466 
467  CALL para_env%sum(mpi_buffer)
468  IF (unit_nr > 0) THEN
469  WRITE (unit_nr, fmt="(A,1X,I10,1X)", advance='no') label, iatom
470  WRITE (unit_nr, *) mpi_buffer
471  END IF
472  DEALLOCATE (mpi_buffer)
473  END DO
474 
475  ! flush
476  IF (unit_nr > 0) FLUSH (unit_nr)
477 
478  END SUBROUTINE pao_write_diagonal_blocks
479 
480 ! **************************************************************************************************
481 !> \brief Writes header of restart file
482 !> \param pao ...
483 !> \param qs_env ...
484 !> \param energy ...
485 !> \param unit_nr ...
486 ! **************************************************************************************************
487  SUBROUTINE write_restart_header(pao, qs_env, energy, unit_nr)
488  TYPE(pao_env_type), POINTER :: pao
489  TYPE(qs_environment_type), POINTER :: qs_env
490  REAL(dp) :: energy
491  INTEGER, INTENT(IN) :: unit_nr
492 
493  CHARACTER(LEN=default_string_length) :: kindname
494  INTEGER :: iatom, ikind, ipot, nparams, &
495  pao_basis_size, z
496  TYPE(atomic_kind_type), DIMENSION(:), POINTER :: atomic_kind_set
497  TYPE(cell_type), POINTER :: cell
498  TYPE(gto_basis_set_type), POINTER :: basis_set
499  TYPE(pao_potential_type), DIMENSION(:), POINTER :: pao_potentials
500  TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
501  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
502 
503  CALL get_qs_env(qs_env, &
504  cell=cell, &
505  particle_set=particle_set, &
506  atomic_kind_set=atomic_kind_set, &
507  qs_kind_set=qs_kind_set)
508 
509  WRITE (unit_nr, "(A,5X,I0)") "Version", file_format_version
510  WRITE (unit_nr, "(A,5X,F20.10)") "Energy", energy
511  WRITE (unit_nr, "(A,5X,I0)") "Step", pao%istep
512  WRITE (unit_nr, "(A,5X,A)") "Parametrization", id2str(pao%parameterization)
513 
514  ! write kinds
515  WRITE (unit_nr, "(A,5X,I0)") "Nkinds", SIZE(atomic_kind_set)
516  DO ikind = 1, SIZE(atomic_kind_set)
517  CALL get_atomic_kind(atomic_kind_set(ikind), name=kindname, z=z)
518  CALL get_qs_kind(qs_kind_set(ikind), &
519  pao_basis_size=pao_basis_size, &
520  pao_potentials=pao_potentials, &
521  basis_set=basis_set)
522  CALL pao_param_count(pao, qs_env, ikind, nparams)
523  WRITE (unit_nr, "(A,5X,I10,1X,A,1X,I3)") "Kind", ikind, trim(kindname), z
524  WRITE (unit_nr, "(A,5X,I10,1X,I3)") "NParams", ikind, nparams
525  WRITE (unit_nr, "(A,5X,I10,1X,I10,1X,A)") "PrimBasis", ikind, basis_set%nsgf, trim(basis_set%name)
526  WRITE (unit_nr, "(A,5X,I10,1X,I3)") "PaoBasis", ikind, pao_basis_size
527  WRITE (unit_nr, "(A,5X,I10,1X,I3)") "NPaoPotentials", ikind, SIZE(pao_potentials)
528  DO ipot = 1, SIZE(pao_potentials)
529  WRITE (unit_nr, "(A,5X,I10,1X,I3)", advance='no') "PaoPotential", ikind, ipot
530  WRITE (unit_nr, "(1X,I3)", advance='no') pao_potentials(ipot)%maxl
531  WRITE (unit_nr, "(1X,I3)", advance='no') pao_potentials(ipot)%max_projector
532  WRITE (unit_nr, "(1X,F20.16)", advance='no') pao_potentials(ipot)%beta
533  WRITE (unit_nr, "(1X,F20.16)") pao_potentials(ipot)%weight
534  END DO
535  END DO
536 
537  ! write cell
538  WRITE (unit_nr, fmt="(A,5X)", advance='no') "Cell"
539  WRITE (unit_nr, *) cell%hmat*angstrom
540 
541  ! write atoms
542  WRITE (unit_nr, "(A,5X,I0)") "Natoms", SIZE(particle_set)
543  DO iatom = 1, SIZE(particle_set)
544  kindname = particle_set(iatom)%atomic_kind%name
545  WRITE (unit_nr, fmt="(A,5X,I10,5X,A,1X)", advance='no') "Atom ", iatom, trim(kindname)
546  WRITE (unit_nr, *) particle_set(iatom)%r*angstrom
547  END DO
548 
549  END SUBROUTINE write_restart_header
550 
551 !**************************************************************************************************
552 !> \brief writing the KS matrix (in terms of the PAO basis) in csr format into a file
553 !> \param qs_env qs environment
554 !> \param ls_scf_env ls environment
555 !> \author Mohammad Hossein Bani-Hashemian
556 ! **************************************************************************************************
557  SUBROUTINE pao_write_ks_matrix_csr(qs_env, ls_scf_env)
558  TYPE(qs_environment_type), POINTER :: qs_env
559  TYPE(ls_scf_env_type), TARGET :: ls_scf_env
560 
561  CHARACTER(len=*), PARAMETER :: routinen = 'pao_write_ks_matrix_csr'
562 
563  CHARACTER(LEN=default_path_length) :: file_name, fileformat
564  INTEGER :: handle, ispin, output_unit, unit_nr
565  LOGICAL :: bin, do_kpoints, do_ks_csr_write, uptr
566  REAL(kind=dp) :: thld
567  TYPE(cp_logger_type), POINTER :: logger
568  TYPE(dbcsr_csr_type) :: ks_mat_csr
569  TYPE(dbcsr_type) :: matrix_ks_nosym
570  TYPE(section_vals_type), POINTER :: dft_section, input
571 
572  CALL timeset(routinen, handle)
573 
574  NULLIFY (dft_section)
575 
576  logger => cp_get_default_logger()
577  output_unit = cp_logger_get_default_io_unit(logger)
578 
579  CALL get_qs_env(qs_env, input=input)
580  dft_section => section_vals_get_subs_vals(input, "DFT")
581  do_ks_csr_write = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
582  "PRINT%KS_CSR_WRITE"), cp_p_file)
583 
584  ! NOTE: k-points has to be treated differently later. k-points has KS matrix as double pointer.
585  CALL get_qs_env(qs_env=qs_env, do_kpoints=do_kpoints)
586 
587  IF (do_ks_csr_write .AND. (.NOT. do_kpoints)) THEN
588  CALL section_vals_val_get(dft_section, "PRINT%KS_CSR_WRITE%THRESHOLD", r_val=thld)
589  CALL section_vals_val_get(dft_section, "PRINT%KS_CSR_WRITE%UPPER_TRIANGULAR", l_val=uptr)
590  CALL section_vals_val_get(dft_section, "PRINT%KS_CSR_WRITE%BINARY", l_val=bin)
591 
592  IF (bin) THEN
593  fileformat = "UNFORMATTED"
594  ELSE
595  fileformat = "FORMATTED"
596  END IF
597 
598  DO ispin = 1, SIZE(ls_scf_env%matrix_ks)
599 
600  IF (dbcsr_has_symmetry(ls_scf_env%matrix_ks(ispin))) THEN
601  CALL dbcsr_desymmetrize(ls_scf_env%matrix_ks(ispin), matrix_ks_nosym)
602  ELSE
603  CALL dbcsr_copy(matrix_ks_nosym, ls_scf_env%matrix_ks(ispin))
604  END IF
605 
606  CALL dbcsr_csr_create_from_dbcsr(matrix_ks_nosym, ks_mat_csr, dbcsr_csr_dbcsr_blkrow_dist)
607  CALL dbcsr_convert_dbcsr_to_csr(matrix_ks_nosym, ks_mat_csr)
608 
609  WRITE (file_name, '(A,I0)') "PAO_KS_SPIN_", ispin
610  unit_nr = cp_print_key_unit_nr(logger, dft_section, "PRINT%KS_CSR_WRITE", &
611  extension=".csr", middle_name=trim(file_name), &
612  file_status="REPLACE", file_form=fileformat)
613  CALL dbcsr_csr_write(ks_mat_csr, unit_nr, upper_triangle=uptr, threshold=thld, binary=bin)
614 
615  CALL cp_print_key_finished_output(unit_nr, logger, dft_section, "PRINT%KS_CSR_WRITE")
616 
617  CALL dbcsr_csr_destroy(ks_mat_csr)
618  CALL dbcsr_release(matrix_ks_nosym)
619  END DO
620  END IF
621 
622  CALL timestop(handle)
623 
624  END SUBROUTINE pao_write_ks_matrix_csr
625 
626 !**************************************************************************************************
627 !> \brief writing the overlap matrix (in terms of the PAO basis) in csr format into a file
628 !> \param qs_env qs environment
629 !> \param ls_scf_env ls environment
630 !> \author Mohammad Hossein Bani-Hashemian
631 ! **************************************************************************************************
632  SUBROUTINE pao_write_s_matrix_csr(qs_env, ls_scf_env)
633  TYPE(qs_environment_type), POINTER :: qs_env
634  TYPE(ls_scf_env_type), TARGET :: ls_scf_env
635 
636  CHARACTER(len=*), PARAMETER :: routinen = 'pao_write_s_matrix_csr'
637 
638  CHARACTER(LEN=default_path_length) :: file_name, fileformat
639  INTEGER :: handle, output_unit, unit_nr
640  LOGICAL :: bin, do_kpoints, do_s_csr_write, uptr
641  REAL(kind=dp) :: thld
642  TYPE(cp_logger_type), POINTER :: logger
643  TYPE(dbcsr_csr_type) :: s_mat_csr
644  TYPE(dbcsr_type) :: matrix_s_nosym
645  TYPE(section_vals_type), POINTER :: dft_section, input
646 
647  CALL timeset(routinen, handle)
648 
649  NULLIFY (dft_section)
650 
651  logger => cp_get_default_logger()
652  output_unit = cp_logger_get_default_io_unit(logger)
653 
654  CALL get_qs_env(qs_env, input=input)
655  dft_section => section_vals_get_subs_vals(input, "DFT")
656  do_s_csr_write = btest(cp_print_key_should_output(logger%iter_info, dft_section, &
657  "PRINT%S_CSR_WRITE"), cp_p_file)
658 
659  ! NOTE: k-points has to be treated differently later. k-points has overlap matrix as double pointer.
660  CALL get_qs_env(qs_env=qs_env, do_kpoints=do_kpoints)
661 
662  IF (do_s_csr_write .AND. (.NOT. do_kpoints)) THEN
663  CALL section_vals_val_get(dft_section, "PRINT%S_CSR_WRITE%THRESHOLD", r_val=thld)
664  CALL section_vals_val_get(dft_section, "PRINT%S_CSR_WRITE%UPPER_TRIANGULAR", l_val=uptr)
665  CALL section_vals_val_get(dft_section, "PRINT%S_CSR_WRITE%BINARY", l_val=bin)
666 
667  IF (bin) THEN
668  fileformat = "UNFORMATTED"
669  ELSE
670  fileformat = "FORMATTED"
671  END IF
672 
673  IF (dbcsr_has_symmetry(ls_scf_env%matrix_s)) THEN
674  CALL dbcsr_desymmetrize(ls_scf_env%matrix_s, matrix_s_nosym)
675  ELSE
676  CALL dbcsr_copy(matrix_s_nosym, ls_scf_env%matrix_s)
677  END IF
678 
679  CALL dbcsr_csr_create_from_dbcsr(matrix_s_nosym, s_mat_csr, dbcsr_csr_dbcsr_blkrow_dist)
680  CALL dbcsr_convert_dbcsr_to_csr(matrix_s_nosym, s_mat_csr)
681 
682  WRITE (file_name, '(A,I0)') "PAO_S"
683  unit_nr = cp_print_key_unit_nr(logger, dft_section, "PRINT%S_CSR_WRITE", &
684  extension=".csr", middle_name=trim(file_name), &
685  file_status="REPLACE", file_form=fileformat)
686  CALL dbcsr_csr_write(s_mat_csr, unit_nr, upper_triangle=uptr, threshold=thld, binary=bin)
687 
688  CALL cp_print_key_finished_output(unit_nr, logger, dft_section, "PRINT%S_CSR_WRITE")
689 
690  CALL dbcsr_csr_destroy(s_mat_csr)
691  CALL dbcsr_release(matrix_s_nosym)
692  END IF
693 
694  CALL timestop(handle)
695 
696  END SUBROUTINE pao_write_s_matrix_csr
697 
698 END MODULE pao_io
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
cell_types
Handles all functions related to the CELL.
Definition: cell_types.F:15
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:308
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:119
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_io_unit
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
Definition: cp_log_handling.F:515
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:803
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:1013
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
dm_ls_scf_types
Types needed for a linear scaling quickstep SCF run based on the density matrix.
Definition: dm_ls_scf_types.F:15
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:724
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:1040
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:23
pao_input
Definition: pao_input.F:8
pao_input::id2str
character(len=11) function, public id2str(id)
Helper routine.
Definition: pao_input.F:216
pao_io
Routines for reading and writing restart files.
Definition: pao_io.F:12
pao_io::pao_kinds_ensure_equal
subroutine, public pao_kinds_ensure_equal(pao, qs_env, ikind, pao_kind)
Ensure that the kind read from the restart is equal to the kind curretly in use.
Definition: pao_io.F:326
pao_io::pao_write_restart
subroutine, public pao_write_restart(pao, qs_env, energy)
Writes restart file.
Definition: pao_io.F:387
pao_io::pao_write_ks_matrix_csr
subroutine, public pao_write_ks_matrix_csr(qs_env, ls_scf_env)
writing the KS matrix (in terms of the PAO basis) in csr format into a file
Definition: pao_io.F:558
pao_io::pao_read_raw
subroutine, public pao_read_raw(filename, param, hmat, kinds, atom2kind, positions, xblocks, ml_range)
Reads a restart file into temporary datastructures.
Definition: pao_io.F:183
pao_io::pao_write_s_matrix_csr
subroutine, public pao_write_s_matrix_csr(qs_env, ls_scf_env)
writing the overlap matrix (in terms of the PAO basis) in csr format into a file
Definition: pao_io.F:633
pao_io::pao_read_restart
subroutine, public pao_read_restart(pao, qs_env)
Reads restart file.
Definition: pao_io.F:87
pao_param
Front-End for any PAO parametrization.
Definition: pao_param.F:12
pao_param::pao_param_count
subroutine, public pao_param_count(pao, qs_env, ikind, nparams)
Returns the number of parameters for given atomic kind.
Definition: pao_param.F:173
pao_potentials
Factory routines for potentials used e.g. by pao_param_exp and pao_ml.
Definition: pao_potentials.F:12
pao_types
Types used by the PAO machinery.
Definition: pao_types.F:12
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::angstrom
real(kind=dp), parameter, public angstrom
Definition: physcon.F:144
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, 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_nonbond, sab_almo, sab_kp, sab_kp_nosym, 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, 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, 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, ecoul_1c, rho0_s_rs, rho0_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, 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, rhs)
Get the QUICKSTEP environment.
Definition: qs_environment_types.F:502
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, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, 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_r3d_rs_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, U_of_dft_plus_u, J_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, J0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition: qs_kind_types.F:451