|
CP2K 2.4 (Revision 12889)
|
00001 !-----------------------------------------------------------------------------! 00002 ! CP2K: A general program to perform molecular dynamics simulations ! 00003 ! Copyright (C) 2000 - 2013 CP2K developers group ! 00004 !-----------------------------------------------------------------------------! 00005 00006 ! ***************************************************************************** 00013 MODULE wiener_process 00014 00015 USE atomic_kind_list_types, ONLY: atomic_kind_list_type 00016 USE cp_para_types, ONLY: cp_para_env_type 00017 USE cp_subsys_types, ONLY: cp_subsys_get,& 00018 cp_subsys_type 00019 USE distribution_1d_types, ONLY: distribution_1d_type,& 00020 init_local_particle_set 00021 USE f77_blas 00022 USE force_env_types, ONLY: force_env_get,& 00023 force_env_type 00024 USE input_section_types, ONLY: section_vals_get_subs_vals,& 00025 section_vals_type 00026 USE kinds, ONLY: dp 00027 USE md_environment_types, ONLY: get_md_env,& 00028 md_environment_type,& 00029 need_per_atom_wiener_process 00030 USE metadynamics_types, ONLY: meta_env_type 00031 USE parallel_rng_types, ONLY: GAUSSIAN,& 00032 create_rng_stream,& 00033 next_rng_seed 00034 USE particle_list_types, ONLY: particle_list_type 00035 USE simpar_types, ONLY: simpar_type 00036 USE string_utilities, ONLY: compress 00037 #include "cp_common_uses.h" 00038 00039 IMPLICIT NONE 00040 00041 PRIVATE 00042 00043 ! Global parameters in this module 00044 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'wiener_process' 00045 00046 ! Public subroutines 00047 PUBLIC :: create_wiener_process, create_wiener_process_cv 00048 00049 CONTAINS 00050 00051 ! ***************************************************************************** 00058 SUBROUTINE create_wiener_process(md_env,error) 00059 00060 TYPE(md_environment_type), POINTER :: md_env 00061 TYPE(cp_error_type), INTENT(INOUT) :: error 00062 00063 CHARACTER(LEN=*), PARAMETER :: routineN = 'create_wiener_process', 00064 routineP = moduleN//':'//routineN 00065 00066 CHARACTER(LEN=40) :: name 00067 INTEGER :: iparticle, iparticle_kind, iparticle_local, nparticle, 00068 nparticle_kind, nparticle_local, stat 00069 LOGICAL :: failure 00070 REAL(KIND=dp), ALLOCATABLE, 00071 DIMENSION(:, :, :) :: seed 00072 REAL(KIND=dp), DIMENSION(3, 2) :: initial_seed 00073 TYPE(atomic_kind_list_type), POINTER :: atomic_kinds 00074 TYPE(cp_para_env_type), POINTER :: para_env 00075 TYPE(cp_subsys_type), POINTER :: subsys 00076 TYPE(distribution_1d_type), POINTER :: local_particles 00077 TYPE(force_env_type), POINTER :: force_env 00078 TYPE(particle_list_type), POINTER :: particles 00079 TYPE(section_vals_type), POINTER :: force_env_section, 00080 subsys_section, work_section 00081 TYPE(simpar_type), POINTER :: simpar 00082 00083 failure = .FALSE. 00084 NULLIFY(work_section,force_env) 00085 CPPrecondition (ASSOCIATED(md_env),cp_failure_level,routineP,error,failure) 00086 00087 CALL get_md_env(md_env=md_env, force_env=force_env, para_env=para_env,& 00088 simpar=simpar, error=error) 00089 00090 ![NB] shouldn't the calling process know if it's needed 00091 IF (need_per_atom_wiener_process(md_env, error=error)) THEN 00092 00093 ! Load initial seed (not needed for a restart) 00094 00095 initial_seed = next_rng_seed(error=error) 00096 00097 CALL force_env_get(force_env,force_env_section=force_env_section,& 00098 subsys=subsys,error=error) 00099 00100 subsys_section => section_vals_get_subs_vals(force_env_section,"SUBSYS",error=error) 00101 00102 CALL cp_subsys_get(subsys=subsys,atomic_kinds=atomic_kinds,local_particles=local_particles,& 00103 particles=particles,error=error) 00104 00105 nparticle_kind = atomic_kinds%n_els 00106 nparticle = particles%n_els 00107 00108 ! Allocate the (local) data structures for the Wiener process 00109 00110 ALLOCATE(local_particles%local_particle_set(nparticle_kind),STAT=stat) 00111 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00112 00113 DO iparticle_kind=1,nparticle_kind 00114 nparticle_local = local_particles%n_el(iparticle_kind) 00115 ALLOCATE (local_particles%local_particle_set(iparticle_kind)%rng(nparticle_local),STAT=stat) 00116 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00117 DO iparticle_local=1,nparticle_local 00118 NULLIFY (local_particles%local_particle_set(iparticle_kind)%rng(iparticle_local)%stream) 00119 END DO 00120 END DO 00121 00122 ! Each process generates all seeds. The seed generation should be 00123 ! quite fast and in this way a broadcast is avoided. 00124 00125 ALLOCATE (seed(3,2,nparticle),STAT=stat) 00126 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00127 00128 seed(:,:,1) = initial_seed 00129 DO iparticle=2,nparticle 00130 seed(:,:,iparticle) = next_rng_seed(seed(:,:,iparticle-1),error=error) 00131 END DO 00132 00133 ! Update initial seed 00134 00135 initial_seed = next_rng_seed(seed(:,:,nparticle),error=error) 00136 00137 ! Create a random number stream (Wiener process) for each particle 00138 00139 DO iparticle_kind=1,nparticle_kind 00140 nparticle_local = local_particles%n_el(iparticle_kind) 00141 DO iparticle_local=1,nparticle_local 00142 iparticle = local_particles%list(iparticle_kind)%array(iparticle_local) 00143 WRITE (UNIT=name,FMT="(A,I8)") "Wiener process for particle",iparticle 00144 CALL compress(name) 00145 CALL create_rng_stream(rng_stream=local_particles%local_particle_set(iparticle_kind)%& 00146 &rng(iparticle_local)%stream,name=name,distribution_type=GAUSSIAN,& 00147 extended_precision=.TRUE., seed=seed(:,:,iparticle),error=error) 00148 END DO 00149 END DO 00150 00151 DEALLOCATE (seed,STAT=stat) 00152 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00153 00154 ! Possibly restart Wiener process 00155 NULLIFY (work_section) 00156 work_section => section_vals_get_subs_vals(section_vals=subsys_section,& 00157 subsection_name="RNG_INIT", error=error) 00158 CALL init_local_particle_set(distribution_1d=local_particles,& 00159 nparticle_kind=nparticle_kind, & 00160 work_section=work_section,error=error) 00161 END IF 00162 00163 END SUBROUTINE create_wiener_process 00164 00165 ! ***************************************************************************** 00173 SUBROUTINE create_wiener_process_cv(meta_env,error) 00174 00175 TYPE(meta_env_type), POINTER :: meta_env 00176 TYPE(cp_error_type), INTENT(INOUT) :: error 00177 00178 CHARACTER(LEN=*), PARAMETER :: routineN = 'create_wiener_process_cv', 00179 routineP = moduleN//':'//routineN 00180 00181 CHARACTER(LEN=40) :: name 00182 INTEGER :: i_c, stat 00183 LOGICAL :: failure 00184 REAL(KIND=dp), ALLOCATABLE, 00185 DIMENSION(:, :, :) :: seed 00186 REAL(KIND=dp), DIMENSION(3, 2) :: initial_seed 00187 00188 failure=.FALSE. 00189 IF (.NOT.ASSOCIATED(meta_env)) RETURN 00190 00191 initial_seed = next_rng_seed(error=error) 00192 00193 DO i_c=1,meta_env%n_colvar 00194 NULLIFY (meta_env%rng(i_c)%stream) 00195 END DO 00196 00197 ! Each process generates all seeds. The seed generation should be 00198 ! quite fast and in this way a broadcast is avoided. 00199 00200 ALLOCATE (seed(3,2,meta_env%n_colvar),STAT=stat) 00201 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00202 00203 seed(:,:,1) = initial_seed 00204 DO i_c=2,meta_env%n_colvar 00205 seed(:,:,i_c) = next_rng_seed(seed(:,:,i_c-1),error=error) 00206 END DO 00207 00208 ! Update initial seed 00209 initial_seed = next_rng_seed(seed(:,:,meta_env%n_colvar),error=error) 00210 00211 ! Create a random number stream (Wiener process) for each particle 00212 DO i_c=1,meta_env%n_colvar 00213 WRITE (UNIT=name,FMT="(A,I8)") "Wiener process for COLVAR",i_c 00214 CALL compress(name) 00215 CALL create_rng_stream(rng_stream=meta_env%rng(i_c)%stream,name=name,distribution_type=GAUSSIAN,& 00216 extended_precision=.TRUE., seed=seed(:,:,i_c),error=error) 00217 END DO 00218 DEALLOCATE (seed,STAT=stat) 00219 CPPostcondition(stat==0,cp_failure_level,routineP,error,failure) 00220 00221 END SUBROUTINE create_wiener_process_cv 00222 00223 END MODULE wiener_process
1.7.3