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 Handles the type to compute averages during an MD
      10             : !> \author Teodoro Laino [tlaino] - 03.2008 - University of Zurich
      11             : ! **************************************************************************************************
      12             : MODULE averages_types
      13             :    USE cell_types,                      ONLY: cell_type
      14             :    USE colvar_utils,                    ONLY: get_clv_force,&
      15             :                                               number_of_colvar
      16             :    USE cp_log_handling,                 ONLY: cp_get_default_logger,&
      17             :                                               cp_logger_type,&
      18             :                                               cp_to_string
      19             :    USE cp_output_handling,              ONLY: cp_print_key_finished_output,&
      20             :                                               cp_print_key_unit_nr
      21             :    USE force_env_types,                 ONLY: force_env_type
      22             :    USE input_section_types,             ONLY: section_vals_get,&
      23             :                                               section_vals_get_subs_vals,&
      24             :                                               section_vals_remove_values,&
      25             :                                               section_vals_type,&
      26             :                                               section_vals_val_get
      27             :    USE kinds,                           ONLY: default_string_length,&
      28             :                                               dp
      29             :    USE md_ener_types,                   ONLY: md_ener_type
      30             :    USE virial_types,                    ONLY: virial_type
      31             : #include "../base/base_uses.f90"
      32             : 
      33             :    IMPLICIT NONE
      34             : 
      35             :    PRIVATE
      36             : 
      37             : ! **************************************************************************************************
      38             :    TYPE average_quantities_type
      39             :       INTEGER                                       :: ref_count = 0, itimes_start = -1
      40             :       LOGICAL                                       :: do_averages = .FALSE.
      41             :       TYPE(section_vals_type), POINTER              :: averages_section => NULL()
      42             :       ! Real Scalar Quantities
      43             :       REAL(KIND=dp)                                 :: avetemp = 0.0_dp, avepot = 0.0_dp, avekin = 0.0_dp, &
      44             :                                                        avevol = 0.0_dp, aveca = 0.0_dp, avecb = 0.0_dp, avecc = 0.0_dp
      45             :       REAL(KIND=dp)                                 :: avetemp_baro = 0.0_dp, avehugoniot = 0.0_dp, avecpu = 0.0_dp
      46             :       REAL(KIND=dp)                                 :: aveal = 0.0_dp, avebe = 0.0_dp, avega = 0.0_dp, avepress = 0.0_dp, &
      47             :                                                        avekinc = 0.0_dp, avetempc = 0.0_dp, avepxx = 0.0_dp
      48             :       REAL(KIND=dp)                                 :: avetemp_qm = 0.0_dp, avekin_qm = 0.0_dp, econs = 0.0_dp
      49             :       ! Virial
      50             :       TYPE(virial_type), POINTER                    :: virial => NULL()
      51             :       ! Colvar
      52             :       REAL(KIND=dp), POINTER, DIMENSION(:)          :: avecolvar => NULL()
      53             :       REAL(KIND=dp), POINTER, DIMENSION(:)          :: aveMmatrix => NULL()
      54             :    END TYPE average_quantities_type
      55             : 
      56             : ! **************************************************************************************************
      57             :    INTERFACE get_averages
      58             :       MODULE PROCEDURE get_averages_rs, get_averages_rv, get_averages_rm
      59             :    END INTERFACE get_averages
      60             : 
      61             : ! *** Public subroutines and data types ***
      62             :    PUBLIC :: average_quantities_type, create_averages, release_averages, &
      63             :              retain_averages, compute_averages
      64             : 
      65             : ! *** Global parameters ***
      66             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'averages_types'
      67             : 
      68             : CONTAINS
      69             : 
      70             : ! **************************************************************************************************
      71             : !> \brief Creates averages environment
      72             : !> \param averages ...
      73             : !> \param averages_section ...
      74             : !> \param virial_avg ...
      75             : !> \param force_env ...
      76             : !> \author Teodoro Laino [tlaino] - 03.2008 - University of Zurich
      77             : ! **************************************************************************************************
      78        1806 :    SUBROUTINE create_averages(averages, averages_section, virial_avg, force_env)
      79             :       TYPE(average_quantities_type), POINTER             :: averages
      80             :       TYPE(section_vals_type), POINTER                   :: averages_section
      81             :       LOGICAL, INTENT(IN), OPTIONAL                      :: virial_avg
      82             :       TYPE(force_env_type), POINTER                      :: force_env
      83             : 
      84             :       INTEGER                                            :: i, nint
      85             :       LOGICAL                                            :: do_colvars
      86             : 
      87        1806 :       ALLOCATE (averages)
      88             :       ! Point to the averages section
      89        1806 :       averages%averages_section => averages_section
      90        1806 :       averages%ref_count = 1
      91        1806 :       CALL section_vals_val_get(averages_section, "_SECTION_PARAMETERS_", l_val=averages%do_averages)
      92        1806 :       IF (averages%do_averages) THEN
      93             :          ! Setup Virial if requested
      94        1804 :          IF (PRESENT(virial_avg)) THEN
      95          20 :             IF (virial_avg) THEN
      96        4960 :                ALLOCATE (averages%virial)
      97             :             END IF
      98             :          END IF
      99        1804 :          CALL section_vals_val_get(averages_section, "AVERAGE_COLVAR", l_val=do_colvars)
     100             :          ! Total number of COLVARs
     101        1804 :          nint = 0
     102        1804 :          IF (do_colvars) nint = number_of_colvar(force_env)
     103        3610 :          ALLOCATE (averages%avecolvar(nint))
     104        3610 :          ALLOCATE (averages%aveMmatrix(nint*nint))
     105        1900 :          DO i = 1, nint
     106        1900 :             averages%avecolvar(i) = 0.0_dp
     107             :          END DO
     108        8216 :          DO i = 1, nint*nint
     109        6412 :             averages%aveMmatrix(i) = 0.0_dp
     110             :          END DO
     111             :       END IF
     112        1806 :    END SUBROUTINE create_averages
     113             : 
     114             : ! **************************************************************************************************
     115             : !> \brief retains the given averages env
     116             : !> \param averages ...
     117             : !> \author Teodoro Laino [tlaino] - 03.2008 - University of Zurich
     118             : ! **************************************************************************************************
     119          20 :    SUBROUTINE retain_averages(averages)
     120             :       TYPE(average_quantities_type), POINTER             :: averages
     121             : 
     122          20 :       CPASSERT(ASSOCIATED(averages))
     123          20 :       CPASSERT(averages%ref_count > 0)
     124          20 :       averages%ref_count = averages%ref_count + 1
     125          20 :    END SUBROUTINE retain_averages
     126             : 
     127             : ! **************************************************************************************************
     128             : !> \brief releases the given averages env
     129             : !> \param averages ...
     130             : !> \author Teodoro Laino [tlaino] - 03.2008 - University of Zurich
     131             : ! **************************************************************************************************
     132        1826 :    SUBROUTINE release_averages(averages)
     133             :       TYPE(average_quantities_type), POINTER             :: averages
     134             : 
     135             :       TYPE(section_vals_type), POINTER                   :: work_section
     136             : 
     137        1826 :       IF (ASSOCIATED(averages)) THEN
     138        1826 :          CPASSERT(averages%ref_count > 0)
     139        1826 :          averages%ref_count = averages%ref_count - 1
     140        1826 :          IF (averages%ref_count == 0) THEN
     141        1806 :             IF (ASSOCIATED(averages%virial)) DEALLOCATE (averages%virial)
     142        1806 :             IF (ASSOCIATED(averages%avecolvar)) THEN
     143        1804 :                DEALLOCATE (averages%avecolvar)
     144             :             END IF
     145        1806 :             IF (ASSOCIATED(averages%aveMmatrix)) THEN
     146        1804 :                DEALLOCATE (averages%aveMmatrix)
     147             :             END IF
     148             :             ! Removes restart values from the corresponding restart section..
     149        1806 :             work_section => section_vals_get_subs_vals(averages%averages_section, "RESTART_AVERAGES")
     150        1806 :             CALL section_vals_remove_values(work_section)
     151        1806 :             NULLIFY (averages%averages_section)
     152             :             !
     153        1806 :             DEALLOCATE (averages)
     154             :          END IF
     155             :       END IF
     156             : 
     157        1826 :    END SUBROUTINE release_averages
     158             : 
     159             : ! **************************************************************************************************
     160             : !> \brief computes the averages
     161             : !> \param averages ...
     162             : !> \param force_env ...
     163             : !> \param md_ener ...
     164             : !> \param cell ...
     165             : !> \param virial ...
     166             : !> \param pv_scalar ...
     167             : !> \param pv_xx ...
     168             : !> \param used_time ...
     169             : !> \param hugoniot ...
     170             : !> \param abc ...
     171             : !> \param cell_angle ...
     172             : !> \param nat ...
     173             : !> \param itimes ...
     174             : !> \param time ...
     175             : !> \param my_pos ...
     176             : !> \param my_act ...
     177             : !> \author Teodoro Laino [tlaino] - 03.2008 - University of Zurich
     178             : ! **************************************************************************************************
     179       83106 :    SUBROUTINE compute_averages(averages, force_env, md_ener, cell, virial, &
     180             :                                pv_scalar, pv_xx, used_time, hugoniot, abc, cell_angle, nat, itimes, &
     181             :                                time, my_pos, my_act)
     182             :       TYPE(average_quantities_type), POINTER             :: averages
     183             :       TYPE(force_env_type), POINTER                      :: force_env
     184             :       TYPE(md_ener_type), POINTER                        :: md_ener
     185             :       TYPE(cell_type), POINTER                           :: cell
     186             :       TYPE(virial_type), POINTER                         :: virial
     187             :       REAL(KIND=dp), INTENT(IN)                          :: pv_scalar, pv_xx
     188             :       REAL(KIND=dp), POINTER                             :: used_time
     189             :       REAL(KIND=dp), INTENT(IN)                          :: hugoniot
     190             :       REAL(KIND=dp), DIMENSION(3), INTENT(IN)            :: abc, cell_angle
     191             :       INTEGER, INTENT(IN)                                :: nat, itimes
     192             :       REAL(KIND=dp), INTENT(IN)                          :: time
     193             :       CHARACTER(LEN=default_string_length), INTENT(IN)   :: my_pos, my_act
     194             : 
     195             :       CHARACTER(len=*), PARAMETER                        :: routineN = 'compute_averages'
     196             : 
     197             :       CHARACTER(LEN=default_string_length)               :: ctmp
     198             :       INTEGER                                            :: delta_t, handle, i, nint, output_unit
     199             :       LOGICAL                                            :: restart_averages
     200             :       REAL(KIND=dp)                                      :: start_time
     201       41553 :       REAL(KIND=dp), DIMENSION(:), POINTER               :: cvalues, Mmatrix, tmp
     202             :       TYPE(cp_logger_type), POINTER                      :: logger
     203             :       TYPE(section_vals_type), POINTER                   :: restart_section
     204             : 
     205       41553 :       CALL timeset(routineN, handle)
     206             :       CALL section_vals_val_get(averages%averages_section, "ACQUISITION_START_TIME", &
     207       41553 :                                 r_val=start_time)
     208       41553 :       IF (averages%do_averages) THEN
     209       41533 :          NULLIFY (cvalues, Mmatrix, logger)
     210       41533 :          logger => cp_get_default_logger()
     211             :          ! Determine the nr. of internal colvar (if any/requested)
     212       41533 :          nint = 0
     213       41533 :          IF (ASSOCIATED(averages%avecolvar)) nint = SIZE(averages%avecolvar)
     214             : 
     215             :          ! Evaluate averages if we collected enough statistics (user defined)
     216       41533 :          IF (time >= start_time) THEN
     217             : 
     218             :             ! Handling properly the restart
     219       41519 :             IF (averages%itimes_start == -1) THEN
     220        1770 :                restart_section => section_vals_get_subs_vals(averages%averages_section, "RESTART_AVERAGES")
     221        1770 :                CALL section_vals_get(restart_section, explicit=restart_averages)
     222        1770 :                IF (restart_averages) THEN
     223         172 :                   CALL section_vals_val_get(restart_section, "ITIMES_START", i_val=averages%itimes_start)
     224         172 :                   CALL section_vals_val_get(restart_section, "AVECPU", r_val=averages%avecpu)
     225         172 :                   CALL section_vals_val_get(restart_section, "AVEHUGONIOT", r_val=averages%avehugoniot)
     226         172 :                   CALL section_vals_val_get(restart_section, "AVETEMP_BARO", r_val=averages%avetemp_baro)
     227         172 :                   CALL section_vals_val_get(restart_section, "AVEPOT", r_val=averages%avepot)
     228         172 :                   CALL section_vals_val_get(restart_section, "AVEKIN", r_val=averages%avekin)
     229         172 :                   CALL section_vals_val_get(restart_section, "AVETEMP", r_val=averages%avetemp)
     230         172 :                   CALL section_vals_val_get(restart_section, "AVEKIN_QM", r_val=averages%avekin_qm)
     231         172 :                   CALL section_vals_val_get(restart_section, "AVETEMP_QM", r_val=averages%avetemp_qm)
     232         172 :                   CALL section_vals_val_get(restart_section, "AVEVOL", r_val=averages%avevol)
     233         172 :                   CALL section_vals_val_get(restart_section, "AVECELL_A", r_val=averages%aveca)
     234         172 :                   CALL section_vals_val_get(restart_section, "AVECELL_B", r_val=averages%avecb)
     235         172 :                   CALL section_vals_val_get(restart_section, "AVECELL_C", r_val=averages%avecc)
     236         172 :                   CALL section_vals_val_get(restart_section, "AVEALPHA", r_val=averages%aveal)
     237         172 :                   CALL section_vals_val_get(restart_section, "AVEBETA", r_val=averages%avebe)
     238         172 :                   CALL section_vals_val_get(restart_section, "AVEGAMMA", r_val=averages%avega)
     239         172 :                   CALL section_vals_val_get(restart_section, "AVE_ECONS", r_val=averages%econs)
     240             :                   ! Virial
     241         172 :                   IF (virial%pv_availability) THEN
     242          76 :                      CALL section_vals_val_get(restart_section, "AVE_PRESS", r_val=averages%avepress)
     243          76 :                      CALL section_vals_val_get(restart_section, "AVE_PXX", r_val=averages%avepxx)
     244          76 :                      IF (ASSOCIATED(averages%virial)) THEN
     245           0 :                         CALL section_vals_val_get(restart_section, "AVE_PV_TOT", r_vals=tmp)
     246           0 :                         averages%virial%pv_total = RESHAPE(tmp, (/3, 3/))
     247           0 :                         CALL section_vals_val_get(restart_section, "AVE_PV_VIR", r_vals=tmp)
     248           0 :                         averages%virial%pv_virial = RESHAPE(tmp, (/3, 3/))
     249           0 :                         CALL section_vals_val_get(restart_section, "AVE_PV_KIN", r_vals=tmp)
     250           0 :                         averages%virial%pv_kinetic = RESHAPE(tmp, (/3, 3/))
     251           0 :                         CALL section_vals_val_get(restart_section, "AVE_PV_CNSTR", r_vals=tmp)
     252           0 :                         averages%virial%pv_constraint = RESHAPE(tmp, (/3, 3/))
     253           0 :                         CALL section_vals_val_get(restart_section, "AVE_PV_XC", r_vals=tmp)
     254           0 :                         averages%virial%pv_xc = RESHAPE(tmp, (/3, 3/))
     255           0 :                         CALL section_vals_val_get(restart_section, "AVE_PV_FOCK_4C", r_vals=tmp)
     256           0 :                         averages%virial%pv_fock_4c = RESHAPE(tmp, (/3, 3/))
     257             :                      END IF
     258             :                   END IF
     259             :                   ! Colvars
     260         172 :                   IF (nint > 0) THEN
     261           0 :                      CALL section_vals_val_get(restart_section, "AVE_COLVARS", r_vals=cvalues)
     262           0 :                      CALL section_vals_val_get(restart_section, "AVE_MMATRIX", r_vals=Mmatrix)
     263           0 :                      CPASSERT(nint == SIZE(cvalues))
     264           0 :                      CPASSERT(nint*nint == SIZE(Mmatrix))
     265           0 :                      averages%avecolvar = cvalues
     266           0 :                      averages%aveMmatrix = Mmatrix
     267             :                   END IF
     268             :                ELSE
     269        1598 :                   averages%itimes_start = itimes
     270             :                END IF
     271             :             END IF
     272       41519 :             delta_t = itimes - averages%itimes_start + 1
     273             : 
     274             :             ! Perform averages
     275        1604 :             SELECT CASE (delta_t)
     276             :             CASE (1)
     277        1604 :                averages%avecpu = used_time
     278        1604 :                averages%avehugoniot = hugoniot
     279        1604 :                averages%avetemp_baro = md_ener%temp_baro
     280        1604 :                averages%avepot = md_ener%epot
     281        1604 :                averages%avekin = md_ener%ekin
     282        1604 :                averages%avetemp = md_ener%temp_part
     283        1604 :                averages%avekin_qm = md_ener%ekin_qm
     284        1604 :                averages%avetemp_qm = md_ener%temp_qm
     285        1604 :                averages%avevol = cell%deth
     286        1604 :                averages%aveca = abc(1)
     287        1604 :                averages%avecb = abc(2)
     288        1604 :                averages%avecc = abc(3)
     289        1604 :                averages%aveal = cell_angle(3)
     290        1604 :                averages%avebe = cell_angle(2)
     291        1604 :                averages%avega = cell_angle(1)
     292        1604 :                averages%econs = 0._dp
     293             :                ! Virial
     294        1604 :                IF (virial%pv_availability) THEN
     295         224 :                   averages%avepress = pv_scalar
     296         224 :                   averages%avepxx = pv_xx
     297         224 :                   IF (ASSOCIATED(averages%virial)) THEN
     298         520 :                      averages%virial%pv_total = virial%pv_total
     299         520 :                      averages%virial%pv_virial = virial%pv_virial
     300         520 :                      averages%virial%pv_kinetic = virial%pv_kinetic
     301         520 :                      averages%virial%pv_constraint = virial%pv_constraint
     302         520 :                      averages%virial%pv_xc = virial%pv_xc
     303         520 :                      averages%virial%pv_fock_4c = virial%pv_fock_4c
     304             :                   END IF
     305             :                END IF
     306             :                ! Colvars
     307        1604 :                IF (nint > 0) THEN
     308             :                   CALL get_clv_force(force_env, nsize_xyz=nat*3, nsize_int=nint, &
     309           2 :                                      cvalues=averages%avecolvar, Mmatrix=averages%aveMmatrix)
     310             :                END IF
     311             :             CASE DEFAULT
     312       39915 :                CALL get_averages(averages%avecpu, used_time, delta_t)
     313       39915 :                CALL get_averages(averages%avehugoniot, hugoniot, delta_t)
     314       39915 :                CALL get_averages(averages%avetemp_baro, md_ener%temp_baro, delta_t)
     315       39915 :                CALL get_averages(averages%avepot, md_ener%epot, delta_t)
     316       39915 :                CALL get_averages(averages%avekin, md_ener%ekin, delta_t)
     317       39915 :                CALL get_averages(averages%avetemp, md_ener%temp_part, delta_t)
     318       39915 :                CALL get_averages(averages%avekin_qm, md_ener%ekin_qm, delta_t)
     319       39915 :                CALL get_averages(averages%avetemp_qm, md_ener%temp_qm, delta_t)
     320       39915 :                CALL get_averages(averages%avevol, cell%deth, delta_t)
     321       39915 :                CALL get_averages(averages%aveca, abc(1), delta_t)
     322       39915 :                CALL get_averages(averages%avecb, abc(2), delta_t)
     323       39915 :                CALL get_averages(averages%avecc, abc(3), delta_t)
     324       39915 :                CALL get_averages(averages%aveal, cell_angle(3), delta_t)
     325       39915 :                CALL get_averages(averages%avebe, cell_angle(2), delta_t)
     326       39915 :                CALL get_averages(averages%avega, cell_angle(1), delta_t)
     327       39915 :                CALL get_averages(averages%econs, md_ener%delta_cons, delta_t)
     328             :                ! Virial
     329       39915 :                IF (virial%pv_availability) THEN
     330        4220 :                   CALL get_averages(averages%avepress, pv_scalar, delta_t)
     331        4220 :                   CALL get_averages(averages%avepxx, pv_xx, delta_t)
     332        4220 :                   IF (ASSOCIATED(averages%virial)) THEN
     333         180 :                      CALL get_averages(averages%virial%pv_total, virial%pv_total, delta_t)
     334         180 :                      CALL get_averages(averages%virial%pv_virial, virial%pv_virial, delta_t)
     335         180 :                      CALL get_averages(averages%virial%pv_kinetic, virial%pv_kinetic, delta_t)
     336         180 :                      CALL get_averages(averages%virial%pv_constraint, virial%pv_constraint, delta_t)
     337         180 :                      CALL get_averages(averages%virial%pv_xc, virial%pv_xc, delta_t)
     338         180 :                      CALL get_averages(averages%virial%pv_fock_4c, virial%pv_fock_4c, delta_t)
     339             :                   END IF
     340             :                END IF
     341             :                ! Colvars
     342       81434 :                IF (nint > 0) THEN
     343          72 :                   ALLOCATE (cvalues(nint))
     344          72 :                   ALLOCATE (Mmatrix(nint*nint))
     345             :                   CALL get_clv_force(force_env, nsize_xyz=nat*3, nsize_int=nint, cvalues=cvalues, &
     346          24 :                                      Mmatrix=Mmatrix)
     347          24 :                   CALL get_averages(averages%avecolvar, cvalues, delta_t)
     348          24 :                   CALL get_averages(averages%aveMmatrix, Mmatrix, delta_t)
     349          24 :                   DEALLOCATE (cvalues)
     350          24 :                   DEALLOCATE (Mmatrix)
     351             :                END IF
     352             :             END SELECT
     353             :          END IF
     354             : 
     355             :          ! Possibly print averages
     356             :          output_unit = cp_print_key_unit_nr(logger, averages%averages_section, "PRINT_AVERAGES", &
     357       41533 :                                             extension=".avg", file_position=my_pos, file_action=my_act)
     358       41533 :          IF (output_unit > 0) THEN
     359             :             WRITE (output_unit, FMT='(A15,1X,"=",1X,G15.9,"   NSTEP #",I15)') &
     360           4 :                "AVECPU", averages%avecpu, itimes, &
     361           4 :                "AVEHUGONIOT", averages%avehugoniot, itimes, &
     362           4 :                "AVETEMP_BARO", averages%avetemp_baro, itimes, &
     363           4 :                "AVEPOT", averages%avepot, itimes, &
     364           4 :                "AVEKIN", averages%avekin, itimes, &
     365           4 :                "AVETEMP", averages%avetemp, itimes, &
     366           4 :                "AVEKIN_QM", averages%avekin_qm, itimes, &
     367           4 :                "AVETEMP_QM", averages%avetemp_qm, itimes, &
     368           4 :                "AVEVOL", averages%avevol, itimes, &
     369           4 :                "AVECELL_A", averages%aveca, itimes, &
     370           4 :                "AVECELL_B", averages%avecb, itimes, &
     371           4 :                "AVECELL_C", averages%avecc, itimes, &
     372           4 :                "AVEALPHA", averages%aveal, itimes, &
     373           4 :                "AVEBETA", averages%avebe, itimes, &
     374           4 :                "AVEGAMMA", averages%avega, itimes, &
     375           8 :                "AVE_ECONS", averages%econs, itimes
     376             :             ! Print the virial
     377           4 :             IF (virial%pv_availability) THEN
     378             :                WRITE (output_unit, FMT='(A15,1X,"=",1X,G15.9,"   NSTEP #",I15)') &
     379           0 :                   "AVE_PRESS", averages%avepress, itimes, &
     380           0 :                   "AVE_PXX", averages%avepxx, itimes
     381           0 :                IF (ASSOCIATED(averages%virial)) THEN
     382             :                   WRITE (output_unit, FMT='(A15,1X,"=",1X,G15.9,"   NSTEP #",I15)') &
     383           0 :                      "AVE_PV_TOT", averages%virial%pv_total, itimes, &
     384           0 :                      "AVE_PV_VIR", averages%virial%pv_virial, itimes, &
     385           0 :                      "AVE_PV_KIN", averages%virial%pv_kinetic, itimes, &
     386           0 :                      "AVE_PV_CNSTR", averages%virial%pv_constraint, itimes, &
     387           0 :                      "AVE_PV_XC", averages%virial%pv_xc, itimes, &
     388           0 :                      "AVE_PV_FOCK_4C", averages%virial%pv_fock_4c, itimes
     389             :                END IF
     390             :             END IF
     391         196 :             DO i = 1, nint
     392         192 :                ctmp = cp_to_string(i)
     393             :                WRITE (output_unit, FMT='(A15,1X,"=",1X,G15.9,"   NSTEP #",I15)') &
     394         196 :                   TRIM("AVE_CV-"//ADJUSTL(ctmp)), averages%avecolvar(i), itimes
     395             :             END DO
     396           4 :             WRITE (output_unit, FMT='(/)')
     397             :          END IF
     398             :          CALL cp_print_key_finished_output(output_unit, logger, averages%averages_section, &
     399       41533 :                                            "PRINT_AVERAGES")
     400             :       END IF
     401       41553 :       CALL timestop(handle)
     402       41553 :    END SUBROUTINE compute_averages
     403             : 
     404             : ! **************************************************************************************************
     405             : !> \brief computes the averages - low level for REAL
     406             : !> \param avg ...
     407             : !> \param add ...
     408             : !> \param delta_t ...
     409             : !> \author Teodoro Laino [tlaino] - 03.2008 - University of Zurich
     410             : ! **************************************************************************************************
     411      647080 :    SUBROUTINE get_averages_rs(avg, add, delta_t)
     412             :       REAL(KIND=dp), INTENT(INOUT)                       :: avg
     413             :       REAL(KIND=dp), INTENT(IN)                          :: add
     414             :       INTEGER, INTENT(IN)                                :: delta_t
     415             : 
     416      647080 :       avg = (avg*REAL(delta_t - 1, dp) + add)/REAL(delta_t, dp)
     417      647080 :    END SUBROUTINE get_averages_rs
     418             : 
     419             : ! **************************************************************************************************
     420             : !> \brief computes the averages - low level for REAL vector
     421             : !> \param avg ...
     422             : !> \param add ...
     423             : !> \param delta_t ...
     424             : !> \author Teodoro Laino [tlaino] - 10.2008 - University of Zurich
     425             : ! **************************************************************************************************
     426          48 :    SUBROUTINE get_averages_rv(avg, add, delta_t)
     427             :       REAL(KIND=dp), DIMENSION(:), INTENT(INOUT)         :: avg
     428             :       REAL(KIND=dp), DIMENSION(:), INTENT(IN)            :: add
     429             :       INTEGER, INTENT(IN)                                :: delta_t
     430             : 
     431             :       INTEGER                                            :: i
     432             :       LOGICAL                                            :: check
     433             : 
     434          48 :       check = SIZE(avg) == SIZE(add)
     435          48 :       CPASSERT(check)
     436       56496 :       DO i = 1, SIZE(avg)
     437       56496 :          avg(i) = (avg(i)*REAL(delta_t - 1, dp) + add(i))/REAL(delta_t, dp)
     438             :       END DO
     439          48 :    END SUBROUTINE get_averages_rv
     440             : 
     441             : ! **************************************************************************************************
     442             : !> \brief computes the averages - low level for REAL matrix
     443             : !> \param avg ...
     444             : !> \param add ...
     445             : !> \param delta_t ...
     446             : !> \author Teodoro Laino [tlaino] - 10.2008 - University of Zurich
     447             : ! **************************************************************************************************
     448        1080 :    SUBROUTINE get_averages_rm(avg, add, delta_t)
     449             :       REAL(KIND=dp), DIMENSION(:, :), INTENT(INOUT)      :: avg
     450             :       REAL(KIND=dp), DIMENSION(:, :), INTENT(IN)         :: add
     451             :       INTEGER, INTENT(IN)                                :: delta_t
     452             : 
     453             :       INTEGER                                            :: i, j
     454             :       LOGICAL                                            :: check
     455             : 
     456        1080 :       check = SIZE(avg, 1) == SIZE(add, 1)
     457        1080 :       CPASSERT(check)
     458        1080 :       check = SIZE(avg, 2) == SIZE(add, 2)
     459        1080 :       CPASSERT(check)
     460        4320 :       DO i = 1, SIZE(avg, 2)
     461       14040 :          DO j = 1, SIZE(avg, 1)
     462       12960 :             avg(j, i) = (avg(j, i)*REAL(delta_t - 1, dp) + add(j, i))/REAL(delta_t, dp)
     463             :          END DO
     464             :       END DO
     465        1080 :    END SUBROUTINE get_averages_rm
     466             : 
     467           0 : END MODULE averages_types