LCOV - code coverage report
Current view: top level - src/motion - pint_io.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:b195825) Lines: 223 271 82.3 %
Date: 2024-04-20 06:29:22 Functions: 8 8 100.0 %

          Line data    Source code
       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  I/O subroutines for pint_env
      10             : !> \author Lukasz Walewski
      11             : !> \date   2009-06-04
      12             : ! **************************************************************************************************
      13             : MODULE pint_io
      14             : 
      15             :    USE cell_types,                      ONLY: cell_type
      16             :    USE cp_log_handling,                 ONLY: cp_get_default_logger,&
      17             :                                               cp_logger_get_default_unit_nr,&
      18             :                                               cp_logger_type
      19             :    USE cp_output_handling,              ONLY: cp_p_file,&
      20             :                                               cp_print_key_finished_output,&
      21             :                                               cp_print_key_should_output,&
      22             :                                               cp_print_key_unit_nr
      23             :    USE cp_subsys_types,                 ONLY: cp_subsys_get,&
      24             :                                               cp_subsys_type
      25             :    USE cp_units,                        ONLY: cp_unit_from_cp2k
      26             :    USE f77_interface,                   ONLY: f_env_add_defaults,&
      27             :                                               f_env_rm_defaults,&
      28             :                                               f_env_type
      29             :    USE force_env_types,                 ONLY: force_env_get
      30             :    USE input_constants,                 ONLY: dump_atomic,&
      31             :                                               dump_dcd,&
      32             :                                               dump_dcd_aligned_cell,&
      33             :                                               dump_xmol
      34             :    USE input_section_types,             ONLY: section_vals_get_subs_vals,&
      35             :                                               section_vals_type,&
      36             :                                               section_vals_val_get
      37             :    USE kinds,                           ONLY: default_string_length,&
      38             :                                               dp
      39             :    USE machine,                         ONLY: m_flush
      40             :    USE particle_list_types,             ONLY: particle_list_type
      41             :    USE particle_methods,                ONLY: write_particle_coordinates
      42             :    USE pint_public,                     ONLY: pint_com_pos
      43             :    USE pint_transformations,            ONLY: pint_u2x
      44             :    USE pint_types,                      ONLY: e_conserved_id,&
      45             :                                               e_kin_thermo_id,&
      46             :                                               e_kin_virial_id,&
      47             :                                               e_potential_id,&
      48             :                                               pint_env_type
      49             : #include "../base/base_uses.f90"
      50             : 
      51             :    IMPLICIT NONE
      52             : 
      53             :    PRIVATE
      54             : 
      55             :    LOGICAL, PRIVATE, PARAMETER :: debug_this_module = .TRUE.
      56             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pint_io'
      57             : 
      58             :    PUBLIC :: pint_write_line
      59             :    PUBLIC :: pint_write_centroids
      60             :    PUBLIC :: pint_write_trajectory
      61             :    PUBLIC :: pint_write_com
      62             :    PUBLIC :: pint_write_ener
      63             :    PUBLIC :: pint_write_action
      64             :    PUBLIC :: pint_write_step_info
      65             :    PUBLIC :: pint_write_rgyr
      66             : 
      67             : CONTAINS
      68             : 
      69             : ! ***************************************************************************
      70             : !> \brief  Writes out a line of text to the default output unit.
      71             : !> \param line ...
      72             : !> \date   2009-07-10
      73             : !> \author Lukasz Walewski
      74             : ! **************************************************************************************************
      75         138 :    SUBROUTINE pint_write_line(line)
      76             : 
      77             :       CHARACTER(len=*), INTENT(IN)                       :: line
      78             : 
      79             :       CHARACTER(len=default_string_length)               :: my_label
      80             :       INTEGER                                            :: unit_nr
      81             :       TYPE(cp_logger_type), POINTER                      :: logger
      82             : 
      83         138 :       NULLIFY (logger)
      84         138 :       logger => cp_get_default_logger()
      85         138 :       my_label = "PINT|"
      86             : 
      87         138 :       IF (logger%para_env%is_source()) THEN
      88          79 :          unit_nr = cp_logger_get_default_unit_nr(logger)
      89          79 :          WRITE (unit_nr, '(T2,A)') TRIM(my_label)//" "//TRIM(line)
      90             :       END IF
      91             : 
      92         138 :    END SUBROUTINE pint_write_line
      93             : 
      94             : ! ***************************************************************************
      95             : !> \brief Write out the trajectory of the centroid (positions and velocities)
      96             : !> \param pint_env ...
      97             : !> \par History
      98             : !>      various bug fixes - hforbert
      99             : !>      2010-11-25 rewritten, added support for velocity printing,
     100             : !>                 calc of the stddev of the beads turned off [lwalewski]
     101             : !> \author fawzi
     102             : ! **************************************************************************************************
     103         494 :    SUBROUTINE pint_write_centroids(pint_env)
     104             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     105             : 
     106             :       CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_centroids'
     107             :       INTEGER, PARAMETER                                 :: n_ids = 2, pos_id = 1, vel_id = 2
     108             : 
     109             :       CHARACTER(len=default_string_length)               :: ext, form, my_middle_name, unit_str
     110             :       CHARACTER(len=default_string_length), DIMENSION(2) :: content_id, middle_name, sect_path, title
     111             :       INTEGER                                            :: handle, handle1, iat, ib, id, idim, &
     112             :                                                             idir, ierr, outformat, should_output, &
     113             :                                                             unit_nr
     114             :       LOGICAL                                            :: new_file, print_kind
     115             :       REAL(kind=dp)                                      :: nb, ss, unit_conv, vv
     116             :       TYPE(cell_type), POINTER                           :: cell
     117             :       TYPE(cp_logger_type), POINTER                      :: logger
     118             :       TYPE(cp_subsys_type), POINTER                      :: subsys
     119             :       TYPE(f_env_type), POINTER                          :: f_env
     120             :       TYPE(particle_list_type), POINTER                  :: particles
     121             :       TYPE(section_vals_type), POINTER                   :: print_key
     122             : 
     123         494 :       CALL timeset(routineN, handle1)
     124             : 
     125         494 :       sect_path(pos_id) = "MOTION%PINT%PRINT%CENTROID_POS"
     126         494 :       sect_path(vel_id) = "MOTION%PINT%PRINT%CENTROID_VEL"
     127         494 :       middle_name(pos_id) = "centroid-pos"
     128         494 :       middle_name(vel_id) = "centroid-vel"
     129         494 :       content_id(pos_id) = "POS"
     130         494 :       content_id(vel_id) = "VEL"
     131             :       WRITE (UNIT=title(pos_id), FMT="(A,I8,A,F20.10)") &
     132         494 :          " i =", pint_env%iter, &
     133        4476 :          ", E =", SUM(pint_env%e_pot_bead)*pint_env%propagator%physpotscale
     134             :       WRITE (UNIT=title(vel_id), FMT="(A,I8,A,F20.10,A,F20.10)") &
     135         494 :          " i =", pint_env%iter, &
     136         494 :          ", E_trm =", pint_env%energy(e_kin_thermo_id), &
     137         988 :          ", E_vir =", pint_env%energy(e_kin_virial_id)
     138             : 
     139         494 :       NULLIFY (logger)
     140         494 :       logger => cp_get_default_logger()
     141             : 
     142         494 :       CALL pint_u2x(pint_env, ux=pint_env%uv, x=pint_env%v)
     143             : 
     144             :       ! iterate over the properties that we know how to print
     145             :       ! (currently positions and velocities)
     146        1482 :       DO id = 1, n_ids
     147             : 
     148             :          print_key => section_vals_get_subs_vals(pint_env%input, &
     149         988 :                                                  TRIM(sect_path(id)))
     150             : 
     151             :          should_output = cp_print_key_should_output( &
     152             :                          iteration_info=logger%iter_info, &
     153         988 :                          basis_section=print_key)
     154             :          IF (.NOT. BTEST(should_output, cp_p_file)) CONTINUE
     155             : 
     156         988 :          print_kind = .FALSE.
     157             : 
     158             :          ! get units of measure for output (if available)
     159             :          CALL section_vals_val_get(print_key, "UNIT", &
     160         988 :                                    c_val=unit_str)
     161         988 :          unit_conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str))
     162             : 
     163             :          ! get the format for output
     164         988 :          CALL section_vals_val_get(print_key, "FORMAT", i_val=outformat)
     165             : 
     166           0 :          SELECT CASE (outformat)
     167             :          CASE (dump_dcd, dump_dcd_aligned_cell)
     168           0 :             form = "UNFORMATTED"
     169           0 :             ext = ".dcd"
     170             :          CASE (dump_atomic)
     171           0 :             form = "FORMATTED"
     172           0 :             ext = ""
     173             :          CASE (dump_xmol)
     174             :             CALL section_vals_val_get(print_key, "PRINT_ATOM_KIND", &
     175         988 :                                       l_val=print_kind)
     176         988 :             form = "FORMATTED"
     177         988 :             ext = ".xyz"
     178             :          CASE default
     179         988 :             CPABORT("")
     180             :          END SELECT
     181             : 
     182         988 :          NULLIFY (f_env, cell, subsys)
     183             :          CALL f_env_add_defaults(f_env_id=pint_env%replicas%f_env_id, &
     184         988 :                                  f_env=f_env, handle=handle)
     185             :          CALL force_env_get(force_env=f_env%force_env, &
     186         988 :                             cell=cell, subsys=subsys)
     187         988 :          CALL cp_subsys_get(subsys, particles=particles)
     188             : 
     189             :          ! calculate and copy the requested property
     190             :          ! to the particles structure
     191         988 :          nb = REAL(pint_env%p, dp)
     192         988 :          idim = 0
     193      224836 :          DO iat = 1, pint_env%ndim/3
     194      896380 :             DO idir = 1, 3
     195      671544 :                idim = idim + 1
     196      671544 :                ss = 0.0_dp
     197      671544 :                vv = 0.0_dp
     198             : !          ss2=0.0_dp
     199     4048776 :                DO ib = 1, pint_env%p
     200     3377232 :                   ss = ss + pint_env%x(ib, idim)
     201     4048776 :                   vv = vv + pint_env%v(ib, idim)
     202             : !            ss2=ss2+pint_env%x(ib,idim)**2
     203             :                END DO
     204      671544 :                particles%els(iat)%r(idir) = ss/nb
     205      895392 :                particles%els(iat)%v(idir) = vv/nb
     206             : !          particles%els(iat)%v(idir)=SQRT(ss2/nb-(ss/nb)**2)
     207             :             END DO
     208             :          END DO
     209             : 
     210             :          ! set up the output unit number and file name
     211             :          ! for the current property
     212         988 :          my_middle_name = TRIM(middle_name(id))
     213             :          unit_nr = cp_print_key_unit_nr(logger=logger, &
     214             :                                         basis_section=print_key, print_key_path="", &
     215             :                                         extension=TRIM(ext), middle_name=TRIM(my_middle_name), &
     216         988 :                                         local=.FALSE., file_form=form, is_new_file=new_file)
     217             : 
     218             :          ! don't write the 0-th frame if the file already exists
     219         988 :          IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN
     220             :             CALL cp_print_key_finished_output(unit_nr, logger, &
     221          82 :                                               print_key)
     222             :             CONTINUE
     223             :          END IF
     224             : 
     225             :          ! actually perform the i/o - on the ionode only
     226         988 :          IF (unit_nr > 0) THEN
     227             : 
     228             :             CALL write_particle_coordinates( &
     229             :                particles%els, &
     230             :                iunit=unit_nr, &
     231             :                output_format=outformat, &
     232             :                content=content_id(id), &
     233             :                title=title(id), &
     234             :                cell=cell, &
     235             :                unit_conv=unit_conv, &
     236         410 :                print_kind=print_kind)
     237             : 
     238             :             CALL cp_print_key_finished_output(unit_nr, logger, &
     239         410 :                                               print_key, "", local=.FALSE.)
     240             : 
     241             :          END IF
     242             : 
     243         988 :          CALL f_env_rm_defaults(f_env, ierr, handle)
     244        3458 :          CPASSERT(ierr == 0)
     245             : 
     246             :       END DO
     247             : 
     248         494 :       CALL timestop(handle1)
     249         494 :    END SUBROUTINE pint_write_centroids
     250             : 
     251             : ! ***************************************************************************
     252             : !> \brief  Write out the trajectory of the beads (positions and velocities)
     253             : !> \param pint_env ...
     254             : !> \par    History
     255             : !>         2010-11-25 added support for velocity printing [lwalewski]
     256             : !> \author hforbert
     257             : ! **************************************************************************************************
     258         494 :    SUBROUTINE pint_write_trajectory(pint_env)
     259             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     260             : 
     261             :       CHARACTER(len=*), PARAMETER :: routineN = 'pint_write_trajectory'
     262             :       INTEGER, PARAMETER                                 :: force_id = 3, n_ids = 3, pos_id = 1, &
     263             :                                                             vel_id = 2
     264             : 
     265             :       CHARACTER(len=default_string_length)               :: ext, form, ib_str, my_middle_name, &
     266             :                                                             title, unit_str
     267             :       CHARACTER(len=default_string_length), DIMENSION(3) :: content_id, middle_name, sect_path
     268             :       INTEGER                                            :: handle, handle1, iat, ib, id, idim, &
     269             :                                                             idir, ierr, imag_stride, outformat, &
     270             :                                                             should_output, unit_nr
     271             :       LOGICAL                                            :: new_file
     272             :       REAL(kind=dp)                                      :: unit_conv
     273             :       TYPE(cell_type), POINTER                           :: cell
     274             :       TYPE(cp_logger_type), POINTER                      :: logger
     275             :       TYPE(cp_subsys_type), POINTER                      :: subsys
     276             :       TYPE(f_env_type), POINTER                          :: f_env
     277             :       TYPE(particle_list_type), POINTER                  :: particles
     278             :       TYPE(section_vals_type), POINTER                   :: print_key
     279             : 
     280         494 :       CALL timeset(routineN, handle1)
     281             : 
     282         494 :       sect_path(pos_id) = "MOTION%PRINT%TRAJECTORY"
     283         494 :       sect_path(vel_id) = "MOTION%PRINT%VELOCITIES"
     284         494 :       sect_path(force_id) = "MOTION%PRINT%FORCES"
     285         494 :       middle_name(pos_id) = "pos-"
     286         494 :       middle_name(vel_id) = "vel-"
     287         494 :       middle_name(force_id) = "force-"
     288         494 :       content_id(pos_id) = "POS"
     289         494 :       content_id(vel_id) = "VEL"
     290         494 :       content_id(force_id) = "FORCE"
     291             : 
     292         494 :       NULLIFY (logger)
     293         494 :       logger => cp_get_default_logger()
     294             : 
     295         494 :       CALL pint_u2x(pint_env, ux=pint_env%uv, x=pint_env%v)
     296             : 
     297             :       ! iterate over the properties that we know how to print
     298             :       ! (currently positions and velocities)
     299        1976 :       DO id = 1, n_ids
     300             : 
     301             :          print_key => section_vals_get_subs_vals(pint_env%input, &
     302        1482 :                                                  TRIM(sect_path(id)))
     303             : 
     304             :          should_output = cp_print_key_should_output( &
     305             :                          iteration_info=logger%iter_info, &
     306        1482 :                          basis_section=print_key)
     307             :          IF (.NOT. BTEST(should_output, cp_p_file)) CONTINUE
     308             : 
     309             :          ! get units of measure for output (if available)
     310             :          CALL section_vals_val_get(print_key, "UNIT", &
     311        1482 :                                    c_val=unit_str)
     312        1482 :          unit_conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str))
     313             : 
     314             :          ! get the format for output
     315        1482 :          CALL section_vals_val_get(print_key, "FORMAT", i_val=outformat)
     316             : 
     317           0 :          SELECT CASE (outformat)
     318             :          CASE (dump_dcd, dump_dcd_aligned_cell)
     319           0 :             form = "UNFORMATTED"
     320           0 :             ext = ".dcd"
     321             :          CASE (dump_atomic)
     322           0 :             form = "FORMATTED"
     323           0 :             ext = ""
     324             :          CASE (dump_xmol)
     325        1482 :             form = "FORMATTED"
     326        1482 :             ext = ".xyz"
     327             :          CASE default
     328        1482 :             CPABORT("")
     329             :          END SELECT
     330             : 
     331        1482 :          NULLIFY (f_env, cell, subsys)
     332             :          CALL f_env_add_defaults(f_env_id=pint_env%replicas%f_env_id, &
     333        1482 :                                  f_env=f_env, handle=handle)
     334             :          CALL force_env_get(force_env=f_env%force_env, &
     335        1482 :                             cell=cell, subsys=subsys)
     336        1482 :          CALL cp_subsys_get(subsys, particles=particles)
     337             : 
     338             :          !Get print stride for bead trajectories
     339             :          CALL section_vals_val_get(pint_env%input, &
     340             :                                    "MOTION%PINT%PRINT%IMAGINARY_TIME_STRIDE", &
     341        1482 :                                    i_val=imag_stride)
     342             : 
     343             :          ! iterate over beads
     344       11946 :          DO ib = 1, pint_env%p, imag_stride
     345             : 
     346             :             ! copy the requested property of the current bead
     347             :             ! to the particles structure
     348       10464 :             idim = 0
     349     1699080 :             DO iat = 1, pint_env%ndim/3
     350     6764928 :                DO idir = 1, 3
     351     5065848 :                   idim = idim + 1
     352     5065848 :                   particles%els(iat)%r(idir) = pint_env%x(ib, idim)
     353     5065848 :                   particles%els(iat)%v(idir) = pint_env%v(ib, idim)
     354     6754464 :                   particles%els(iat)%f(idir) = pint_env%f(ib, idim)
     355             :                END DO
     356             :             END DO
     357             : 
     358             :             ! set up the output unit number and file name
     359             :             ! for the current property and bead
     360       10464 :             ib_str = ""
     361       10464 :             WRITE (ib_str, *) ib
     362       10464 :             my_middle_name = TRIM(middle_name(id))//TRIM(ADJUSTL(ib_str))
     363             :             unit_nr = cp_print_key_unit_nr(logger=logger, &
     364             :                                            basis_section=print_key, print_key_path="", &
     365             :                                            extension=TRIM(ext), middle_name=TRIM(my_middle_name), &
     366       10464 :                                            local=.FALSE., file_form=form, is_new_file=new_file)
     367             : 
     368             :             ! don't write the 0-th frame if the file already exists
     369       10464 :             IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN
     370             :                CALL cp_print_key_finished_output(unit_nr, logger, &
     371         916 :                                                  print_key)
     372             :                CONTINUE
     373             :             END IF
     374             : 
     375             :             ! actually perform the i/o - on the ionode only
     376       11946 :             IF (unit_nr > 0) THEN
     377             : 
     378        1540 :                IF (outformat == dump_xmol) THEN
     379             :                   WRITE (UNIT=title, FMT="(A,I8,A,F20.10)") &
     380        1540 :                      " i =", pint_env%iter, &
     381        3080 :                      ", E =", pint_env%e_pot_bead(ib)
     382             :                END IF
     383             : 
     384             :                CALL write_particle_coordinates( &
     385             :                   particles%els, &
     386             :                   iunit=unit_nr, &
     387             :                   output_format=outformat, &
     388             :                   content=content_id(id), &
     389             :                   title=title, &
     390             :                   cell=cell, &
     391        1540 :                   unit_conv=unit_conv)
     392             : 
     393             :                CALL cp_print_key_finished_output(unit_nr, logger, &
     394        1540 :                                                  print_key, "", local=.FALSE.)
     395             : 
     396             :             END IF
     397             : 
     398             :          END DO
     399             : 
     400        1482 :          CALL f_env_rm_defaults(f_env, ierr, handle)
     401        6422 :          CPASSERT(ierr == 0)
     402             : 
     403             :       END DO
     404             : 
     405         494 :       CALL timestop(handle1)
     406         494 :    END SUBROUTINE pint_write_trajectory
     407             : 
     408             : ! ***************************************************************************
     409             : !> \brief  Write center of mass (COM) position according to PINT%PRINT%COM
     410             : !> \param pint_env ...
     411             : !> \date   2010-02-17
     412             : !> \author Lukasz Walewski
     413             : ! **************************************************************************************************
     414         494 :    SUBROUTINE pint_write_com(pint_env)
     415             : 
     416             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     417             : 
     418             :       CHARACTER(len=default_string_length)               :: stmp1, stmp2
     419             :       INTEGER                                            :: ic, unit_nr
     420             :       LOGICAL                                            :: new_file, should_output
     421             :       REAL(kind=dp), DIMENSION(3)                        :: com_r
     422             :       TYPE(cp_logger_type), POINTER                      :: logger
     423             :       TYPE(section_vals_type), POINTER                   :: print_key
     424             : 
     425         494 :       NULLIFY (logger)
     426         494 :       logger => cp_get_default_logger()
     427             : 
     428             :       ! decide whether to write anything or not
     429         494 :       NULLIFY (print_key)
     430             :       print_key => section_vals_get_subs_vals(pint_env%input, &
     431         494 :                                               "MOTION%PINT%PRINT%COM")
     432             :       should_output = BTEST(cp_print_key_should_output( &
     433             :                             iteration_info=logger%iter_info, &
     434         494 :                             basis_section=print_key), cp_p_file)
     435         494 :       IF (.NOT. should_output) THEN
     436         494 :          RETURN
     437             :       END IF
     438             : 
     439           0 :       com_r = pint_com_pos(pint_env)
     440           0 :       DO ic = 1, 3
     441           0 :          com_r(ic) = cp_unit_from_cp2k(com_r(ic), "angstrom")
     442             :       END DO
     443             : 
     444             :       unit_nr = cp_print_key_unit_nr(logger, print_key, is_new_file=new_file, &
     445           0 :                                      middle_name="com-pos", extension=".xyz")
     446             : 
     447             :       ! don't write the 0-th frame if the file already exists
     448           0 :       IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN
     449             :          CALL cp_print_key_finished_output(unit_nr, logger, &
     450           0 :                                            print_key)
     451           0 :          RETURN
     452             :       END IF
     453             : 
     454             :       ! actually perform the i/o - on the ionode only
     455           0 :       IF (unit_nr > 0) THEN
     456             : 
     457           0 :          WRITE (unit_nr, '(I2)') 1
     458           0 :          WRITE (stmp1, *) pint_env%iter
     459           0 :          WRITE (stmp2, '(F20.10)') pint_env%energy(e_conserved_id)
     460           0 :          WRITE (unit_nr, '(4A)') " Iteration = ", TRIM(ADJUSTL(stmp1)), &
     461           0 :             ", E_conserved = ", TRIM(ADJUSTL(stmp2))
     462           0 :          WRITE (unit_nr, '(A2,3(1X,F20.10))') "X ", (com_r(ic), ic=1, 3)
     463             : 
     464           0 :          CALL m_flush(unit_nr)
     465             : 
     466             :       END IF
     467             : 
     468           0 :       CALL cp_print_key_finished_output(unit_nr, logger, print_key)
     469             : 
     470             :    END SUBROUTINE pint_write_com
     471             : 
     472             : ! ***************************************************************************
     473             : !> \brief  Writes out the energies according to PINT%PRINT%ENERGY
     474             : !> \param  pint_env path integral environment
     475             : !> \par    History
     476             : !>           various bug fixes [hforbert]
     477             : !>           2009-11-16 energy components calc moved out of here [lwalewski]
     478             : !> \author fawzi
     479             : ! **************************************************************************************************
     480         494 :    SUBROUTINE pint_write_ener(pint_env)
     481             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     482             : 
     483             :       INTEGER                                            :: ndof, unit_nr
     484             :       LOGICAL                                            :: file_is_new
     485             :       REAL(kind=dp)                                      :: t, temp
     486             :       TYPE(cp_logger_type), POINTER                      :: logger
     487             :       TYPE(section_vals_type), POINTER                   :: print_key
     488             : 
     489         494 :       NULLIFY (print_key, logger)
     490             :       print_key => section_vals_get_subs_vals(pint_env%input, &
     491         494 :                                               "MOTION%PINT%PRINT%ENERGY")
     492         494 :       logger => cp_get_default_logger()
     493         494 :       IF (BTEST(cp_print_key_should_output(iteration_info=logger%iter_info, &
     494             :                                            basis_section=print_key), cp_p_file)) THEN
     495             : 
     496             :          unit_nr = cp_print_key_unit_nr(logger, print_key, middle_name="energy", &
     497         494 :                                         extension=".dat", is_new_file=file_is_new)
     498             : 
     499             :          ! don't write the 0-th frame if the file already exists
     500         494 :          IF (.NOT. file_is_new .AND. (pint_env%iter .LE. pint_env%first_step)) THEN
     501             :             CALL cp_print_key_finished_output(unit_nr, logger, &
     502          35 :                                               print_key)
     503          35 :             RETURN
     504             :          END IF
     505             : 
     506             :          ! cp_print_key_unit_nr returns -1 on nodes other than logger%para_env%is_source()
     507         459 :          IF (unit_nr > 0) THEN
     508             : 
     509             :             ! please keep the format explanation up to date
     510             :             ! keep the constant of motion the true constant of motion !
     511         240 :             IF (file_is_new) THEN
     512             :                WRITE (unit_nr, "(A8,1X,A12,1X,5(A20,1X),A12)") &
     513          21 :                   "# StepNr", &
     514          21 :                   "   Time [fs]", &
     515          21 :                   "      Kinetic [a.u.]", &
     516          21 :                   "    VirialKin [a.u.]", &
     517          21 :                   "     Temperature [K]", &
     518          21 :                   "    Potential [a.u.]", &
     519          21 :                   "      ConsQty [a.u.]", &
     520          42 :                   "     CPU [s]"
     521             :             END IF
     522             : 
     523         240 :             t = cp_unit_from_cp2k(pint_env%t, "fs")
     524             : 
     525         240 :             ndof = pint_env%p
     526         240 :             IF (pint_env%first_propagated_mode .EQ. 2) THEN
     527           0 :                ndof = ndof - 1
     528             :             END IF
     529             :             temp = cp_unit_from_cp2k(2.0_dp*pint_env%e_kin_beads/ &
     530             :                                      REAL(ndof, dp)/REAL(pint_env%ndim, dp), &
     531         240 :                                      "K")*pint_env%propagator%temp_sim2phys
     532             : 
     533             :             WRITE (unit_nr, "(I8,1X,F12.3,1X,5(F20.9,1X),F12.1)") &
     534         240 :                pint_env%iter, &
     535         240 :                t, &
     536         240 :                pint_env%energy(e_kin_thermo_id), &
     537         240 :                pint_env%energy(e_kin_virial_id), &
     538         240 :                temp, &
     539         240 :                pint_env%energy(e_potential_id), &
     540         240 :                pint_env%energy(e_conserved_id), &
     541         480 :                pint_env%time_per_step
     542         240 :             CALL m_flush(unit_nr)
     543             : 
     544             :          END IF
     545             : 
     546         459 :          CALL cp_print_key_finished_output(unit_nr, logger, print_key)
     547             :       END IF
     548             : 
     549             :    END SUBROUTINE pint_write_ener
     550             : 
     551             : ! ***************************************************************************
     552             : !> \brief  Writes out the actions according to PINT%PRINT%ACTION
     553             : !> \param  pint_env path integral environment
     554             : !> \author Felix Uhl
     555             : ! **************************************************************************************************
     556         494 :    SUBROUTINE pint_write_action(pint_env)
     557             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     558             : 
     559             :       INTEGER                                            :: unit_nr
     560             :       LOGICAL                                            :: file_is_new
     561             :       REAL(kind=dp)                                      :: t
     562             :       TYPE(cp_logger_type), POINTER                      :: logger
     563             :       TYPE(section_vals_type), POINTER                   :: print_key
     564             : 
     565         494 :       NULLIFY (print_key, logger)
     566             :       print_key => section_vals_get_subs_vals(pint_env%input, &
     567         494 :                                               "MOTION%PINT%PRINT%ACTION")
     568         494 :       logger => cp_get_default_logger()
     569         494 :       IF (BTEST(cp_print_key_should_output(iteration_info=logger%iter_info, &
     570             :                                            basis_section=print_key), cp_p_file)) THEN
     571             : 
     572             :          unit_nr = cp_print_key_unit_nr(logger, print_key, middle_name="action", &
     573           0 :                                         extension=".dat", is_new_file=file_is_new)
     574             : 
     575             :          ! don't write the 0-th frame if the file already exists
     576           0 :          IF (.NOT. file_is_new .AND. (pint_env%iter .LE. pint_env%first_step)) THEN
     577             :             CALL cp_print_key_finished_output(unit_nr, logger, &
     578           0 :                                               print_key)
     579           0 :             RETURN
     580             :          END IF
     581             : 
     582             :          ! cp_print_key_unit_nr returns -1 on nodes other than logger%para_env%is_source()
     583           0 :          IF (unit_nr > 0) THEN
     584             : 
     585             :             ! please keep the format explanation up to date
     586             :             ! keep the constant of motion the true constant of motion !
     587           0 :             IF (file_is_new) THEN
     588             :                WRITE (unit_nr, "(A8,1X,A12,1X,2(A25,1X))") &
     589           0 :                   "# StepNr", &
     590           0 :                   "   Time [fs]", &
     591           0 :                   "       Link Action [a.u.]", &
     592           0 :                   "  Potential Action [a.u.]"
     593             :             END IF
     594             : 
     595           0 :             t = cp_unit_from_cp2k(pint_env%t, "fs")
     596             : 
     597             :             WRITE (unit_nr, "(I8,1X,F12.3,1X,5(F20.9,1X),F12.1)") &
     598           0 :                pint_env%iter, &
     599           0 :                t, &
     600           0 :                pint_env%link_action, &
     601           0 :                pint_env%pot_action
     602           0 :             CALL m_flush(unit_nr)
     603             : 
     604             :          END IF
     605             : 
     606           0 :          CALL cp_print_key_finished_output(unit_nr, logger, print_key)
     607             :       END IF
     608             : 
     609             :    END SUBROUTINE pint_write_action
     610             : 
     611             : ! ***************************************************************************
     612             : !> \brief  Write step info to the output file.
     613             : !> \param pint_env ...
     614             : !> \date   2009-11-16
     615             : !> \par History
     616             : !>      2010-01-27 getting default unit nr now only on ionode [lwalewski]
     617             : !> \author Lukasz Walewski
     618             : ! **************************************************************************************************
     619         438 :    SUBROUTINE pint_write_step_info(pint_env)
     620             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     621             : 
     622             :       CHARACTER(len=default_string_length)               :: msgstr, stmp, time_unit
     623             :       INTEGER                                            :: unit_nr
     624             :       REAL(kind=dp)                                      :: time_used
     625             :       TYPE(cp_logger_type), POINTER                      :: logger
     626             : 
     627         438 :       unit_nr = 0
     628         438 :       NULLIFY (logger)
     629         438 :       logger => cp_get_default_logger()
     630             : 
     631         438 :       time_used = pint_env%time_per_step
     632         438 :       time_unit = "sec"
     633         438 :       IF (time_used .GE. 60.0_dp) THEN
     634           0 :          time_used = time_used/60.0_dp
     635           0 :          time_unit = "min"
     636             :       END IF
     637         438 :       IF (time_used .GE. 60.0_dp) THEN
     638           0 :          time_used = time_used/60.0_dp
     639           0 :          time_unit = "hours"
     640             :       END IF
     641         438 :       msgstr = "PINT step"
     642         438 :       stmp = ""
     643         438 :       WRITE (stmp, *) pint_env%iter
     644         438 :       msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(stmp))//" of"
     645         438 :       stmp = ""
     646         438 :       WRITE (stmp, *) pint_env%last_step
     647         438 :       msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(stmp))//" in"
     648         438 :       stmp = ""
     649         438 :       WRITE (stmp, '(F20.1)') time_used
     650         438 :       msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(stmp))
     651         438 :       msgstr = TRIM(ADJUSTL(msgstr))//" "//TRIM(ADJUSTL(time_unit))//"."
     652             : 
     653         438 :       IF (logger%para_env%is_source()) THEN
     654         219 :          unit_nr = cp_logger_get_default_unit_nr(logger)
     655         219 :          WRITE (unit_nr, '(T2,A)') "PINT| "//TRIM(ADJUSTL(msgstr))
     656             :       END IF
     657             : 
     658             :       ! print out the total energy - for regtest evaluation
     659         438 :       stmp = ""
     660         438 :       WRITE (stmp, *) pint_env%energy(e_conserved_id)
     661         438 :       msgstr = "Total energy = "//TRIM(ADJUSTL(stmp))
     662         438 :       IF (logger%para_env%is_source()) THEN
     663         219 :          WRITE (unit_nr, '(T2,A)') "PINT| "//TRIM(ADJUSTL(msgstr))
     664             :       END IF
     665             : 
     666         438 :    END SUBROUTINE pint_write_step_info
     667             : 
     668             : ! ***************************************************************************
     669             : !> \brief  Write radii of gyration according to PINT%PRINT%CENTROID_GYR
     670             : !> \param pint_env ...
     671             : !> \date   2011-01-07
     672             : !> \author Lukasz Walewski
     673             : ! **************************************************************************************************
     674         494 :    SUBROUTINE pint_write_rgyr(pint_env)
     675             : 
     676             :       TYPE(pint_env_type), INTENT(IN)                    :: pint_env
     677             : 
     678             :       CHARACTER(len=default_string_length)               :: unit_str
     679             :       INTEGER                                            :: ia, ib, ic, idim, unit_nr
     680             :       LOGICAL                                            :: new_file, should_output
     681             :       REAL(kind=dp)                                      :: nb, ss, unit_conv
     682             :       TYPE(cp_logger_type), POINTER                      :: logger
     683             :       TYPE(section_vals_type), POINTER                   :: print_key
     684             : 
     685         494 :       NULLIFY (logger)
     686         494 :       logger => cp_get_default_logger()
     687             : 
     688             :       ! decide whether to write anything or not
     689         494 :       NULLIFY (print_key)
     690             :       print_key => section_vals_get_subs_vals(pint_env%input, &
     691         494 :                                               "MOTION%PINT%PRINT%CENTROID_GYR")
     692             :       should_output = BTEST(cp_print_key_should_output( &
     693             :                             iteration_info=logger%iter_info, &
     694         494 :                             basis_section=print_key), cp_p_file)
     695         494 :       IF (.NOT. should_output) THEN
     696          99 :          RETURN
     697             :       END IF
     698             : 
     699             :       ! get the units conversion factor
     700         422 :       CALL section_vals_val_get(print_key, "UNIT", c_val=unit_str)
     701         422 :       unit_conv = cp_unit_from_cp2k(1.0_dp, TRIM(unit_str))
     702             : 
     703             :       ! calculate the centroid positions
     704         422 :       nb = REAL(pint_env%p, dp)
     705         422 :       idim = 0
     706        1754 :       DO ia = 1, pint_env%ndim/3
     707        5750 :          DO ic = 1, 3
     708        3996 :             idim = idim + 1
     709        3996 :             ss = 0.0_dp
     710       33732 :             DO ib = 1, pint_env%p
     711       33732 :                ss = ss + pint_env%x(ib, idim)
     712             :             END DO
     713        5328 :             pint_env%rtmp_ndim(idim) = ss/nb
     714             :          END DO
     715             :       END DO
     716             : 
     717             :       ! calculate the radii of gyration
     718             :       idim = 0
     719        1754 :       DO ia = 1, pint_env%ndim/3
     720             :          ss = 0.0_dp
     721        5328 :          DO ic = 1, 3
     722        3996 :             idim = idim + 1
     723       35064 :             DO ib = 1, pint_env%p
     724       33732 :                ss = ss + (pint_env%x(ib, idim) - pint_env%rtmp_ndim(idim))**2
     725             :             END DO
     726             :          END DO
     727        1754 :          pint_env%rtmp_natom(ia) = SQRT(ss/nb)*unit_conv
     728             :       END DO
     729             : 
     730             :       unit_nr = cp_print_key_unit_nr(logger, print_key, is_new_file=new_file, &
     731         422 :                                      middle_name="centroid-gyr", extension=".dat")
     732             : 
     733             :       ! don't write the 0-th frame if the file already exists
     734         422 :       IF (.NOT. new_file .AND. (pint_env%iter .LE. pint_env%first_step)) THEN
     735             :          CALL cp_print_key_finished_output(unit_nr, logger, &
     736          27 :                                            print_key)
     737          27 :          RETURN
     738             :       END IF
     739             : 
     740             :       ! actually perform the i/o - on the ionode only
     741         395 :       IF (unit_nr > 0) THEN
     742             : 
     743         853 :          DO ia = 1, pint_env%ndim/3
     744         853 :             WRITE (unit_nr, '(F20.10,1X)', ADVANCE='NO') pint_env%rtmp_natom(ia)
     745             :          END DO
     746         205 :          WRITE (unit_nr, '(A)') ""
     747             : 
     748         205 :          CALL m_flush(unit_nr)
     749             : 
     750             :       END IF
     751             : 
     752         395 :       CALL cp_print_key_finished_output(unit_nr, logger, print_key)
     753             : 
     754             :    END SUBROUTINE pint_write_rgyr
     755             : 
     756             : END MODULE pint_io

Generated by: LCOV version 1.15