Small Voronoi Changes II (#1466)

* Add APPEND flag for &PRINT / &VORONOI. * Write unwrapped coordinates to .voronoi files. * Write additional comment line with cell vector to .voronoi files.
cp2k · Apr 13, 2021 · 3bec98f · 3bec98f
1 parent 264012e
commit 3bec98f
Show file tree

Hide file tree

Showing 3 changed files with 28 additions and 10 deletions.
diff --git a/src/input_cp2k_dft.F b/src/input_cp2k_dft.F
@@ -1626,6 +1626,12 @@ SUBROUTINE create_print_dft_section(section)
                                        print_level=debug_print_level + 1, filename="", &
                                        citations=(/Rycroft2009, Thomas2015, Brehm2018, Brehm2020, Brehm2021/))
 
+      CALL keyword_create(keyword, __LOCATION__, name="APPEND", &
+                          description="Appends frames to already existing .voronoi file.", &
+                          default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
+      CALL section_add_keyword(print_key, keyword)
+      CALL keyword_release(keyword)
+
       CALL keyword_create(keyword, __LOCATION__, name="SANITY_CHECK", &
                           description="Performs a sanity check before each Voronoi integration, i.e.,"// &
                           " checks if every grid point is located in exactly one Voronoi cell.", &
@@ -1635,7 +1641,7 @@ SUBROUTINE create_print_dft_section(section)
 
       CALL keyword_create(keyword, __LOCATION__, name="OVERWRITE", &
                           description="Specify this keyword to overwrite any existing ""properties.emp"" file if"// &
-                          " it already exists. By default, the data is appended to an existing file.", &
+                          " it already exists. By default, the data is appended to an existing .emp file.", &
                           usage="OVERWRITE T", default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
       CALL section_add_keyword(print_key, keyword)
       CALL keyword_release(keyword)
@@ -1653,13 +1659,15 @@ SUBROUTINE create_print_dft_section(section)
       CALL keyword_release(keyword)
 
       CALL keyword_create(keyword, __LOCATION__, name="OUTPUT_EMP", &
-                          description="Writes the resulting electromagnetic moments to a binary file ""properties.emp"".", &
+                          description="Writes the resulting electromagnetic moments to a binary file ""properties.emp""."// &
+                          " The file name cannot be changed.", &
                           usage="OUTPUT_EMP T", default_l_val=.FALSE., lone_keyword_l_val=.TRUE.)
       CALL section_add_keyword(print_key, keyword)
       CALL keyword_release(keyword)
 
       CALL keyword_create(keyword, __LOCATION__, name="OUTPUT_TEXT", &
-                          description="Writes the resulting electromagnetic moments to text files (*.voronoi).", &
+                          description="Writes the resulting electromagnetic moments to text files (*.voronoi)."// &
+                          " The file name is specified via FILENAME.", &
                           usage="OUTPUT_TEXT T", default_l_val=.TRUE., lone_keyword_l_val=.TRUE.)
       CALL section_add_keyword(print_key, keyword)
       CALL keyword_release(keyword)

diff --git a/src/qs_scf_post_gpw.F b/src/qs_scf_post_gpw.F
@@ -1820,13 +1820,14 @@ SUBROUTINE write_mo_free_results(qs_env)
       CHARACTER(len=1), DIMENSION(3), PARAMETER          :: cdir = (/"x", "y", "z"/)
 
       CHARACTER(LEN=2)                                   :: element_symbol
-      CHARACTER(LEN=default_path_length)                 :: filename, mpi_filename, my_pos_cube
+      CHARACTER(LEN=default_path_length)                 :: filename, mpi_filename, my_pos_cube, &
+                                                            my_pos_voro
       CHARACTER(LEN=default_string_length)               :: name
       INTEGER :: after, handle, i, iat, id, ikind, img, iso, ispin, iw, l, nd(3), ngto, niso, &
          nkind, np, nr, output_unit, print_level, should_print_bqb, should_print_voro, unit_nr, &
          unit_nr_voro
-      LOGICAL :: append_cube, do_kpoints, mpi_io, omit_headers, print_it, print_total_density, &
-         rho_r_valid, voro_print_txt, write_ks, write_xc, xrd_interface
+      LOGICAL :: append_cube, append_voro, do_kpoints, mpi_io, omit_headers, print_it, &
+         print_total_density, rho_r_valid, voro_print_txt, write_ks, write_xc, xrd_interface
       REAL(KIND=dp)                                      :: q_max, rho_hard, rho_soft, rho_total, &
                                                             rho_total_rspace, udvol, volume
       REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :)        :: bfun
@@ -2541,8 +2542,13 @@ SUBROUTINE write_mo_free_results(qs_env)
             IF (should_print_voro /= 0) THEN
                CALL section_vals_val_get(print_key_voro, "OUTPUT_TEXT", l_val=voro_print_txt)
                IF (voro_print_txt) THEN
+                  append_voro = section_get_lval(input, "DFT%PRINT%VORONOI%APPEND")
+                  my_pos_voro = "REWIND"
+                  IF (append_voro) THEN
+                     my_pos_voro = "APPEND"
+                  END IF
                   unit_nr_voro = cp_print_key_unit_nr(logger, input, "DFT%PRINT%VORONOI", extension=".voronoi", &
-                                                      log_filename=.FALSE.)
+                                                      file_position=my_pos_voro, log_filename=.FALSE.)
                ELSE
                   unit_nr_voro = 0
                END IF

diff --git a/src/voronoi_interface.F b/src/voronoi_interface.F
@@ -643,9 +643,13 @@ SUBROUTINE entry_voronoi_or_bqb(do_voro, do_bqb, input_voro, input_bqb, unit_vor
             END DO
 
             IF (unit_voro > 0) THEN
-               WRITE (unit_voro, FMT="(I0)") natom
-               WRITE (unit_voro, FMT="(A,I8,A,F12.4,A)") "# Simulation step ", sim_step, ",  Simulation time ", &
+               WRITE (unit_voro, FMT="(T2,I0)") natom
+               WRITE (unit_voro, FMT="(A,I8,A,F12.4,A)") "# Step ", sim_step, ",  Time ", &
                   sim_time*femtoseconds, " fs"
+               WRITE (unit_voro, FMT="(A,9F20.10)") "# Cell ", &
+                  cell%hmat(1, 1)*angstrom, cell%hmat(2, 1)*angstrom, cell%hmat(3, 1)*angstrom, &
+                  cell%hmat(1, 2)*angstrom, cell%hmat(2, 2)*angstrom, cell%hmat(3, 2)*angstrom, &
+                  cell%hmat(1, 3)*angstrom, cell%hmat(2, 3)*angstrom, cell%hmat(3, 3)*angstrom
                WRITE (unit_voro, FMT="(A,22A20)") "# Atom     Z", &
                   "Radius", "Position(X)", "Position(Y)", "Position(Z)", &
                   "Voronoi_Volume", "Z(eff)", "Charge", "Dipole(X)", "Dipole(Y)", "Dipole(Z)", &
@@ -658,7 +662,7 @@ SUBROUTINE entry_voronoi_or_bqb(do_voro, do_bqb, input_voro, input_bqb, unit_vor
                      i1, &
                      particles_z(i1), &
                      voro_radii(i1)/100.0_dp, &
-                     voro_wrapped_pos(i1, 1:3)/100.0_dp, &
+                     particles_r(1:3, i1)*angstrom, &
                      voro_volume(i1)/1000000.0_dp, &
                      particles_c(i1), &
                      voro_charge(i1), &