LCOV - code coverage report
Current view: top level - src/motion - pint_gle.F (source / functions) Coverage Total Hit
Test: CP2K Regtests (git:42dac4a) Lines: 100.0 % 61 61
Test Date: 2025-07-25 12:55:17 Functions: 100.0 % 3 3

            Line data    Source code
       1              : !--------------------------------------------------------------------------------------------------!
       2              : !   CP2K: A general program to perform molecular dynamics simulations                              !
       3              : !   Copyright 2000-2025 CP2K developers group <https://cp2k.org>                                   !
       4              : !                                                                                                  !
       5              : !   SPDX-License-Identifier: GPL-2.0-or-later                                                      !
       6              : !--------------------------------------------------------------------------------------------------!
       7              : 
       8              : ! **************************************************************************************************
       9              : !> \brief  Methods to apply GLE to PI runs.
      10              : !> \author michelec
      11              : !> \par History
      12              : !>      06.2010 created [michelec]
      13              : !> \note   trying to keep duplication at a minimum....
      14              : ! **************************************************************************************************
      15              : 
      16              : MODULE pint_gle
      17              :    USE gle_system_dynamics,             ONLY: gle_cholesky_stab
      18              :    USE gle_system_types,                ONLY: gle_type
      19              :    USE kinds,                           ONLY: dp
      20              :    USE pint_types,                      ONLY: pint_env_type
      21              : #include "../base/base_uses.f90"
      22              : 
      23              :    IMPLICIT NONE
      24              : 
      25              :    PRIVATE
      26              : 
      27              :    PUBLIC :: pint_gle_step, pint_gle_init, pint_calc_gle_energy
      28              : 
      29              :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'pint_gle'
      30              : 
      31              : CONTAINS
      32              : ! **************************************************************************************************
      33              : !> \brief ...
      34              : !> \param pint_env ...
      35              : ! **************************************************************************************************
      36            6 :    ELEMENTAL SUBROUTINE pint_calc_gle_energy(pint_env)
      37              :       TYPE(pint_env_type), INTENT(INOUT)                 :: pint_env
      38              : 
      39              :       INTEGER                                            :: i
      40              : 
      41            6 :       pint_env%e_gle = 0._dp
      42            6 :       IF (ASSOCIATED(pint_env%gle)) THEN
      43        55302 :          DO i = 1, pint_env%gle%loc_num_gle
      44        55302 :             pint_env%e_gle = pint_env%e_gle + pint_env%gle%nvt(i)%thermostat_energy
      45              :          END DO
      46              :       END IF
      47            6 :    END SUBROUTINE
      48              : 
      49              : ! **************************************************************************************************
      50              : !> \brief ...
      51              : !> \param pint_env ...
      52              : ! **************************************************************************************************
      53            2 :    SUBROUTINE pint_gle_init(pint_env)
      54              :       TYPE(pint_env_type), INTENT(INOUT)                 :: pint_env
      55              : 
      56              :       INTEGER                                            :: i, ib, idim, imap, j
      57            4 :       REAL(dp) :: mf, rr(pint_env%gle%ndim), cc(pint_env%gle%ndim, pint_env%gle%ndim)
      58              : 
      59            2 :       CALL gle_cholesky_stab(pint_env%gle%c_mat, cc, pint_env%gle%ndim)
      60        18434 :       DO i = 1, pint_env%gle%loc_num_gle
      61        18432 :          imap = pint_env%gle%map_info%index(i)
      62        18432 :          ib = 1 + (imap - 1)/pint_env%ndim
      63        18432 :          idim = 1 + MOD(imap - 1, pint_env%ndim)
      64        18432 :          mf = 1.0_dp/SQRT(pint_env%mass_fict(ib, idim))
      65       110592 :          DO j = 1, pint_env%gle%ndim
      66       110592 :             rr(j) = pint_env%gle%nvt(i)%gaussian_rng_stream%next()*mf
      67              :          END DO
      68       663554 :          pint_env%gle%nvt(i)%s = MATMUL(cc, rr)
      69              :       END DO
      70              : 
      71            2 :    END SUBROUTINE pint_gle_init
      72              : 
      73              : ! **************************************************************************************************
      74              : !> \brief ...
      75              : !> \param pint_env ...
      76              : ! **************************************************************************************************
      77            4 :    SUBROUTINE pint_gle_step(pint_env)
      78              :       TYPE(pint_env_type), INTENT(INOUT)                 :: pint_env
      79              : 
      80              :       CHARACTER(len=*), PARAMETER                        :: routineN = 'pint_gle_step'
      81              : 
      82              :       INTEGER                                            :: handle, iadd, ib, ideg, idim, imap, &
      83              :                                                             ndim, num
      84              :       REAL(dp)                                           :: alpha, beta, mf, rr
      85            4 :       REAL(dp), DIMENSION(:, :), POINTER                 :: a_mat, e_tmp, h_tmp, s_tmp
      86              :       TYPE(gle_type), POINTER                            :: gle
      87              : 
      88            4 :       CALL timeset(routineN, handle)
      89              : 
      90            4 :       gle => pint_env%gle
      91            4 :       ndim = gle%ndim
      92              : 
      93           16 :       ALLOCATE (s_tmp(ndim, gle%loc_num_gle))
      94       221188 :       s_tmp = 0.0_dp
      95           12 :       ALLOCATE (e_tmp(ndim, gle%loc_num_gle))
      96           12 :       ALLOCATE (h_tmp(ndim, gle%loc_num_gle))
      97              : 
      98        36868 :       DO ideg = 1, gle%loc_num_gle
      99        36864 :          imap = gle%map_info%index(ideg)
     100        36864 :          ib = 1 + (imap - 1)/pint_env%ndim
     101        36864 :          idim = 1 + MOD(imap - 1, pint_env%ndim)
     102              : 
     103        36864 :          gle%nvt(ideg)%s(1) = pint_env%uv_t(ib, idim)
     104              :          gle%nvt(ideg)%thermostat_energy = gle%nvt(ideg)%thermostat_energy &
     105        36864 :                                            + 0.5_dp*pint_env%mass_fict(ib, idim)*gle%nvt(ideg)%s(1)**2
     106        36864 :          s_tmp(1, imap) = gle%nvt(ideg)%s(1)
     107        36864 :          rr = gle%nvt(ideg)%gaussian_rng_stream%next()
     108        36864 :          mf = 1.0_dp/SQRT(pint_env%mass_fict(ib, idim))
     109        36864 :          e_tmp(1, imap) = rr*mf
     110       184324 :          DO iadd = 2, ndim
     111       147456 :             s_tmp(iadd, imap) = gle%nvt(ideg)%s(iadd)
     112       147456 :             rr = gle%nvt(ideg)%gaussian_rng_stream%next()
     113       184320 :             e_tmp(iadd, imap) = rr*mf
     114              :          END DO
     115              :       END DO
     116            4 :       num = gle%loc_num_gle
     117            4 :       a_mat => gle%gle_s
     118            4 :       alpha = 1.0_dp
     119            4 :       beta = 0.0_dp
     120              : 
     121            4 :       CALL DGEMM('N', 'N', ndim, num, ndim, alpha, a_mat(1, 1), ndim, e_tmp(1, 1), ndim, beta, h_tmp(1, 1), ndim)
     122              : 
     123            4 :       a_mat => gle%gle_t
     124            4 :       beta = 1.0_dp
     125            4 :       CALL DGEMM("N", "N", ndim, num, ndim, alpha, a_mat(1, 1), ndim, s_tmp(1, 1), ndim, beta, h_tmp(1, 1), ndim)
     126              : 
     127        36868 :       DO ideg = 1, gle%loc_num_gle
     128        36864 :          imap = gle%map_info%index(ideg)
     129              : 
     130       221184 :          DO iadd = 1, ndim
     131       221184 :             gle%nvt(ideg)%s(iadd) = h_tmp(iadd, imap)
     132              :          END DO
     133        36864 :          ib = 1 + (imap - 1)/pint_env%ndim
     134        36864 :          idim = 1 + MOD(imap - 1, pint_env%ndim)
     135        36864 :          pint_env%uv_t(ib, idim) = gle%nvt(ideg)%s(1)
     136              :          gle%nvt(ideg)%thermostat_energy = gle%nvt(ideg)%thermostat_energy &
     137        36868 :                                            - 0.5_dp*pint_env%mass_fict(ib, idim)*gle%nvt(ideg)%s(1)**2
     138              :       END DO
     139            4 :       pint_env%e_kin_t = 0.0_dp
     140            4 :       DEALLOCATE (e_tmp, s_tmp, h_tmp)
     141            4 :       CALL timestop(handle)
     142            4 :    END SUBROUTINE
     143              : END MODULE
        

Generated by: LCOV version 2.0-1