LCOV - code coverage report
Current view: top level - src/motion - wiener_process.F (source / functions) Hit Total Coverage
Test: CP2K Regtests (git:b195825) Lines: 60 60 100.0 %
Date: 2024-04-20 06:29:22 Functions: 3 3 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 Handling of the Wiener process currently employed in turn of the
      10             : !>      Langevin dynamics.
      11             : !> \par History
      12             : !>      none
      13             : !> \author Matthias Krack (05.07.2005)
      14             : ! **************************************************************************************************
      15             : MODULE wiener_process
      16             : 
      17             :    USE atomic_kind_list_types,          ONLY: atomic_kind_list_type
      18             :    USE cp_subsys_types,                 ONLY: cp_subsys_get,&
      19             :                                               cp_subsys_type
      20             :    USE distribution_1d_types,           ONLY: distribution_1d_type
      21             :    USE force_env_types,                 ONLY: force_env_get,&
      22             :                                               force_env_type
      23             :    USE input_section_types,             ONLY: section_vals_get,&
      24             :                                               section_vals_get_subs_vals,&
      25             :                                               section_vals_type,&
      26             :                                               section_vals_val_get
      27             :    USE kinds,                           ONLY: dp
      28             :    USE md_environment_types,            ONLY: get_md_env,&
      29             :                                               md_environment_type,&
      30             :                                               need_per_atom_wiener_process
      31             :    USE message_passing,                 ONLY: mp_para_env_type
      32             :    USE metadynamics_types,              ONLY: meta_env_type
      33             :    USE parallel_rng_types,              ONLY: GAUSSIAN,&
      34             :                                               next_rng_seed,&
      35             :                                               rng_record_length,&
      36             :                                               rng_stream_type,&
      37             :                                               rng_stream_type_from_record
      38             :    USE particle_list_types,             ONLY: particle_list_type
      39             :    USE simpar_types,                    ONLY: simpar_type
      40             :    USE string_utilities,                ONLY: compress
      41             : #include "../base/base_uses.f90"
      42             : 
      43             :    IMPLICIT NONE
      44             : 
      45             :    PRIVATE
      46             : 
      47             :    ! Global parameters in this module
      48             :    CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'wiener_process'
      49             : 
      50             :    ! Public subroutines
      51             :    PUBLIC :: create_wiener_process, create_wiener_process_cv
      52             : 
      53             : CONTAINS
      54             : 
      55             : ! **************************************************************************************************
      56             : !> \brief Create a Wiener process for Langevin dynamics and initialize an
      57             : !>      independent random number generator for each atom in all force
      58             : !>      environment and all the subsystems/fragments therein.
      59             : !> \param md_env ...
      60             : !> \par History
      61             : !>      Creation (06.07.2005,MK)
      62             : ! **************************************************************************************************
      63          46 :    SUBROUTINE create_wiener_process(md_env)
      64             : 
      65             :       TYPE(md_environment_type), POINTER                 :: md_env
      66             : 
      67             :       CHARACTER(LEN=40)                                  :: name
      68             :       INTEGER                                            :: iparticle, iparticle_kind, &
      69             :                                                             iparticle_local, nparticle, &
      70             :                                                             nparticle_kind, nparticle_local
      71          46 :       REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :)     :: seed
      72             :       TYPE(atomic_kind_list_type), POINTER               :: atomic_kinds
      73             :       TYPE(cp_subsys_type), POINTER                      :: subsys
      74             :       TYPE(distribution_1d_type), POINTER                :: local_particles
      75             :       TYPE(force_env_type), POINTER                      :: force_env
      76             :       TYPE(mp_para_env_type), POINTER                    :: para_env
      77             :       TYPE(particle_list_type), POINTER                  :: particles
      78             :       TYPE(section_vals_type), POINTER                   :: force_env_section, subsys_section, &
      79             :                                                             work_section
      80             :       TYPE(simpar_type), POINTER                         :: simpar
      81             : 
      82          46 :       NULLIFY (work_section, force_env)
      83          46 :       CPASSERT(ASSOCIATED(md_env))
      84             : 
      85             :       CALL get_md_env(md_env=md_env, force_env=force_env, para_env=para_env, &
      86          46 :                       simpar=simpar)
      87             : 
      88             :       ![NB] shouldn't the calling process know if it's needed
      89          46 :       IF (need_per_atom_wiener_process(md_env)) THEN
      90             : 
      91             :          CALL force_env_get(force_env, force_env_section=force_env_section, &
      92          46 :                             subsys=subsys)
      93             : 
      94          46 :          subsys_section => section_vals_get_subs_vals(force_env_section, "SUBSYS")
      95             : 
      96             :          CALL cp_subsys_get(subsys=subsys, atomic_kinds=atomic_kinds, local_particles=local_particles, &
      97          46 :                             particles=particles)
      98             : 
      99          46 :          nparticle_kind = atomic_kinds%n_els
     100          46 :          nparticle = particles%n_els
     101             : 
     102             :          ! Allocate the (local) data structures for the Wiener process
     103         862 :          ALLOCATE (local_particles%local_particle_set(nparticle_kind))
     104             : 
     105         770 :          DO iparticle_kind = 1, nparticle_kind
     106         724 :             nparticle_local = local_particles%n_el(iparticle_kind)
     107       19345 :             ALLOCATE (local_particles%local_particle_set(iparticle_kind)%rng(nparticle_local))
     108       18129 :             DO iparticle_local = 1, nparticle_local
     109      434699 :                ALLOCATE (local_particles%local_particle_set(iparticle_kind)%rng(iparticle_local)%stream)
     110             :             END DO
     111             :          END DO
     112             : 
     113             :          ! Each process generates all seeds. The seed generation should be
     114             :          ! quite fast and in this way a broadcast is avoided.
     115         138 :          ALLOCATE (seed(3, 2, nparticle))
     116             : 
     117             :          ! Load initial seed (not needed for a restart)
     118         414 :          seed(:, :, 1) = subsys%seed(:, :)
     119             : 
     120       34718 :          DO iparticle = 2, nparticle
     121       34718 :             seed(:, :, iparticle) = next_rng_seed(seed(:, :, iparticle - 1))
     122             :          END DO
     123             : 
     124             :          ! Update initial seed
     125         414 :          subsys%seed(:, :) = next_rng_seed(seed(:, :, nparticle))
     126             : 
     127             :          ! Create a random number stream (Wiener process) for each particle
     128         770 :          DO iparticle_kind = 1, nparticle_kind
     129         724 :             nparticle_local = local_particles%n_el(iparticle_kind)
     130       18129 :             DO iparticle_local = 1, nparticle_local
     131       17359 :                iparticle = local_particles%list(iparticle_kind)%array(iparticle_local)
     132       17359 :                WRITE (UNIT=name, FMT="(A,I8)") "Wiener process for particle", iparticle
     133       17359 :                CALL compress(name)
     134             :                local_particles%local_particle_set(iparticle_kind)%rng(iparticle_local)% &
     135             :                   stream = rng_stream_type(name=name, distribution_type=GAUSSIAN, &
     136       18083 :                                            extended_precision=.TRUE., seed=seed(:, :, iparticle))
     137             :             END DO
     138             :          END DO
     139             : 
     140          46 :          DEALLOCATE (seed)
     141             : 
     142             :          ! Possibly restart Wiener process
     143             :          NULLIFY (work_section)
     144             :          work_section => section_vals_get_subs_vals(section_vals=subsys_section, &
     145          46 :                                                     subsection_name="RNG_INIT")
     146             :          CALL init_local_particle_set(distribution_1d=local_particles, &
     147             :                                       nparticle_kind=nparticle_kind, &
     148          46 :                                       work_section=work_section)
     149             :       END IF
     150             : 
     151          46 :    END SUBROUTINE create_wiener_process
     152             : 
     153             : ! **************************************************************************************************
     154             : !> \brief Helper routine for create_wiener_process.
     155             : !> \param distribution_1d ...
     156             : !> \param nparticle_kind ...
     157             : !> \param work_section ...
     158             : !> \par History
     159             : !>      01.2014 moved from distribution_1d_types (Ole Schuett)
     160             : ! **************************************************************************************************
     161          46 :    SUBROUTINE init_local_particle_set(distribution_1d, nparticle_kind, &
     162             :                                       work_section)
     163             : 
     164             :       TYPE(distribution_1d_type), POINTER                :: distribution_1d
     165             :       INTEGER, INTENT(in)                                :: nparticle_kind
     166             :       TYPE(section_vals_type), POINTER                   :: work_section
     167             : 
     168             :       CHARACTER(LEN=rng_record_length)                   :: rng_record
     169             :       INTEGER                                            :: iparticle, iparticle_kind, &
     170             :                                                             iparticle_local, nparticle_local
     171             :       LOGICAL                                            :: explicit
     172             : 
     173             : ! -------------------------------------------------------------------------
     174             : 
     175          46 :       CPASSERT(ASSOCIATED(distribution_1d))
     176             : 
     177          46 :       IF (ASSOCIATED(work_section)) THEN
     178          46 :          CALL section_vals_get(work_section, explicit=explicit)
     179          46 :          IF (explicit) THEN
     180          18 :             DO iparticle_kind = 1, nparticle_kind
     181          12 :                nparticle_local = distribution_1d%n_el(iparticle_kind)
     182         399 :                DO iparticle_local = 1, nparticle_local
     183         381 :                   iparticle = distribution_1d%list(iparticle_kind)%array(iparticle_local)
     184          12 :                   IF (iparticle == distribution_1d%list(iparticle_kind)%array(iparticle_local)) THEN
     185             :                      CALL section_vals_val_get(section_vals=work_section, &
     186             :                                                keyword_name="_DEFAULT_KEYWORD_", &
     187             :                                                i_rep_val=iparticle, &
     188         381 :                                                c_val=rng_record)
     189             :                      distribution_1d%local_particle_set(iparticle_kind)%rng(iparticle_local)% &
     190         381 :                         stream = rng_stream_type_from_record(rng_record)
     191             :                   END IF
     192             :                END DO
     193             :             END DO
     194             :          END IF
     195             :       END IF
     196             : 
     197          46 :    END SUBROUTINE init_local_particle_set
     198             : 
     199             : ! **************************************************************************************************
     200             : !> \brief Create a Wiener process for Langevin dynamics used for
     201             : !>        metadynamics and initialize an
     202             : !>        independent random number generator for each COLVAR.
     203             : !> \param meta_env ...
     204             : !> \date   01.2009
     205             : !> \author Fabio Sterpone
     206             : !>
     207             : ! **************************************************************************************************
     208           6 :    SUBROUTINE create_wiener_process_cv(meta_env)
     209             : 
     210             :       TYPE(meta_env_type), POINTER                       :: meta_env
     211             : 
     212             :       CHARACTER(LEN=40)                                  :: name
     213             :       INTEGER                                            :: i_c
     214           6 :       REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :, :)     :: seed
     215             :       REAL(KIND=dp), DIMENSION(3, 2)                     :: initial_seed
     216             : 
     217           6 :       IF (.NOT. ASSOCIATED(meta_env)) RETURN
     218             : 
     219           6 :       initial_seed = next_rng_seed()
     220             : 
     221             :       ! Each process generates all seeds. The seed generation should be
     222             :       ! quite fast and in this way a broadcast is avoided.
     223             : 
     224          18 :       ALLOCATE (seed(3, 2, meta_env%n_colvar))
     225             : 
     226          54 :       seed(:, :, 1) = initial_seed
     227          10 :       DO i_c = 2, meta_env%n_colvar
     228          10 :          seed(:, :, i_c) = next_rng_seed(seed(:, :, i_c - 1))
     229             :       END DO
     230             : 
     231             :       ! Update initial seed
     232           6 :       initial_seed = next_rng_seed(seed(:, :, meta_env%n_colvar))
     233             : 
     234             :       ! Create a random number stream (Wiener process) for each particle
     235          16 :       DO i_c = 1, meta_env%n_colvar
     236          10 :          WRITE (UNIT=name, FMT="(A,I8)") "Wiener process for COLVAR", i_c
     237          10 :          CALL compress(name)
     238             :          meta_env%rng(i_c) = rng_stream_type(name=name, distribution_type=GAUSSIAN, &
     239          16 :                                              extended_precision=.TRUE., seed=seed(:, :, i_c))
     240             :       END DO
     241           6 :       DEALLOCATE (seed)
     242             : 
     243             :    END SUBROUTINE create_wiener_process_cv
     244             : 
     245             : END MODULE wiener_process

Generated by: LCOV version 1.15