3. Monte Carlo Simulations Understanding Molecular Simulation - PowerPoint PPT Presentation

3. Monte Carlo Simulations Understanding Molecular Simulation

Molecular Simulations ➡ Molecular dynamics : solve equations of motion r 1 r 2 r n ➡ Monte Carlo : importance sampling r 1 r 2 r n Understanding Molecular Simulation

Monte Carlo Simulations 3. Monte Carlo 3.1.Introduction 3.2.Statistical Thermodynamics (recall) 3.3.Importance sampling 3.4.Details of the algorithm 3.5.Non-Boltzmann sampling 3.6.Parallel Monte Carlo Understanding Molecular Simulation

3. Monte Carlo Simulations 3.2 Statistical Thermodynamics Understanding Molecular Simulation

Canonical ensemble: statistical mechanics Consider a small system that can exchange energy E 1 with a big reservoir =1/k B T ⎛ ⎞ ( ) − ∂ ln Ω ( ) = ln Ω E ln Ω E 1 , E − E 1 ⎟ E 1 + ! ⎜ ∂ E ⎝ ⎠ If the reservoir is very big we can ignore the higher order terms: ( ) ln Ω E 1 , E − E 1 = − E 1 ( ) ln Ω E k B T Hence, the probability to find E 1 : ( ) ( ) Ω E ( ) ( ) Ω E 1 , E − E 1 Ω E 1 , E − E 1 Ω E 1 , E − E 1 ( ) = = = C P E 1 ( ) ( ) Ω E ( ) ( ) ∑ ∑ Ω E i , E − E i Ω E i , E − E i Ω E i i β =1/k B T ⎡ ⎤ ) ∝ exp − E 1 ( ⎡ ⎤ ⎥ ∝ exp − β E 1 P E 1 ⎢ ⎣ ⎦ k B T ⎣ ⎦ Understanding Molecular Simulation

Summary: Canonical ensemble (N,V,T) ( ) Partition function: − U r 1 k B T dr 3 N ∫ Q N , V , T = Λ 3 N N ! e ( ) Probability to find a particular U r 3 N ( ) ∝ e − k B T P r 3 N configuration: Ensemble average: 1 ( ) e ( ) dr 3 N ( ) e ∫ ( ) dr 3 N − β U r A r ∫ − β U r A r Λ 3 N N ! N , V , T = = A ( ) dr 3 N − β U r Q N , V , T ∫ e Free energy: β F = − ln Q NVT Understanding Molecular Simulation

3.Monte Carlo Simulations 3.3 Importance Sampling Understanding Molecular Simulation

Numerical Integration Understanding Molecular Simulation

Monte Carlo simulations Generate M configurations using { } 3 N , r 2 3 N , r 3 3 N , r 4 3 N , ! , r M 3 N r Monte Carlo moves: 1 We can compute the average: ( ) ∑ M A = 3 N A r i i = 1 The probability to generate a ( ) P MC r 3 N ( ) dr 3 N ∫ A r 3 N configuration in our MC A = ( ) dr 3 N P MC r 3 N ∫ scheme: P MC Question: how to chose P MC such that we sample the canonical ensemble? Understanding Molecular Simulation

Ensemble Average ( ) dr 3 N ( ) e 1 1 − β U r 3 N ∫ NVT = A r 3 N A N ! Λ 3 N Q NVT We can rewrite this using the ( ) P r 3 N ( ) dr 3 N probability to find a particular ∫ NVT = A r 3 N A configuration ( ) − β U r 3 N ( ) = e with P r 3 N Λ 3 N N ! Q NVT ( ) dr 3 N ( ) e − β U r 3 N ∫ A r 3 N ( ) P r 3 N ( ) dr 3 N ∫ NVT = = A r 3 N A ( ) dr 3 N − β U r 3 N ∫ e Understanding Molecular Simulation

Monte Carlo - canonical ensemble ( ) dr 3 N ( ) e Canonical ensemble: − β U r 3 N ∫ A r 3 N ( ) P r 3 N ( ) dr 3 N ∫ NVT = = A r 3 N A ( ) dr 3 N − β U r 3 N ∫ e ( ) − β U r 3 N with ( ) = e 2. No need to know P r 3 N Λ 3 N N ! Q NVT the partition function ( ) P MC r 3 N ( ) dr 3 N ∫ A r 3 N Monte Carlo: ( ) ∑ M A = A = 3 N A r i ( ) dr 3 N P MC r 3 N ∫ i = 1 Hence, we need 1. No need ( ) ( ) = Ce P MC r 3 N − β U r 3 N to know C to sample: ( ) dr 3 N ( ) dr 3 N ( ) e ( ) e − β U r 3 N − β U r 3 N ∫ ∫ C A r 3 N A r 3 N A = = = A ( ) dr 3 N ( ) dr 3 N − β U r 3 N − β U r 3 N NVT ∫ ∫ C e e Understanding Molecular Simulation

Importance Sampling: what got lost? Understanding Molecular Simulation

3.Monte Carlo Simulation 3.4 Details of the algorithm Understanding Molecular Simulation

Algorithm 1 (Basic Metropolis Algorithm) basic Metropolis algorithm PROGRAM mc perform ncycl MC cycles do icycl=1,ncycl displace a particle call mcmove if (mod(icycl,nsamp).eq.0) sample averages + call sample enddo end Comments to this algorithm: 1. Subroutine mcmove attempts to displace a randomly selected particle (see Algorithm 2). 2. Subroutine sample samples quantities every nsamp th cycle. Understanding Molecular Simulation

Algorithm 2 (Attempt to Displace a Particle) attempts to displace a particle SUBROUTINE mcmove select a particle at random o=int(ranf()*npart)+1 energy old configuration call ener(x(o),eno) give particle random displacement xn=x(o)+(ranf()-0.5)*delx energy new configuration call ener(xn,enn) acceptance rule (2.2.1) if (ranf().lt.exp(-beta accepted: replace x(o) by xn + *(enn-eno)) x(o)=xn return end Comments to this algorithm: 1. Subroutine ener calculates the energy of a particle at the given position. 2. Note that, if a configuration is rejected, the old configuration is retained. 3. The ranf() is a random number uniform in [ 0, 1 ] . Understanding Molecular Simulation

Questions • How can we prove that this scheme generates the desired distribution of configurations? • Why make a random selection of the particle to be displaced? • Why do we need to take the old configuration again? • How large should we take: delx? Understanding Molecular Simulation

3.Monte Carlo Simulations 3.4.1 Detailed balance Understanding Molecular Simulation

Questions • How can we prove that this scheme generates the desired distribution of configurations? • Why make a random selection of the particle to be displaced? • Why do we need to take the old configuration again? • How large should we take: delx? canonical ensembles Understanding Molecular Simulation

Markov Processes Markov Process Next step only depends on the current state • Ergodic: all possible states can be reached by a set of • single steps Detailed balance • Process will approach a limiting distribution Understanding Molecular Simulation

Ensembles versus probability • P(o) : probability to find the state o • Ensemble: take a very large number ( M ) of identical systems: N(o) = M x P(o) ; the total number of systems in the state o Understanding Molecular Simulation

Markov Processes - Detailed Balance K(o → n) o n K(o → n) : total number of systems in our ensemble that move o → n ( ) = N o ( ) × α o → n ( ) × acc o → n ( ) K o → n • N(o) : total number of systems in our ensemble in state o • α (o → n) : a priori probability to generate a move o → n • acc(o → n) : probability to accept the move o → n Understanding Molecular Simulation

Markov Processes - Detailed Balance K(o → n) n o K(n → o) Condition of detailed balance: ( ) = K n → o ( ) K o → n ( ) = N o ( ) × α o → n ( ) × acc o → n ( ) K o → n ( ) = N n ( ) × α n → o ( ) × acc n → o ( ) K n → o ( ) ( ) × α n → o ( ) acc o → n N n = ( ) ( ) × α o → n ( ) acc n → o N o Understanding Molecular Simulation

NVT-ensemble In the canonical ensemble the number ( ) ∝ e ( ) − β U n N n of configurations in state n is given by: We assume that in our Monte Carlo moves the a priori probability ( ) = α n → o ( ) = α α o → n to perform a move is independent of the configuration: ( ) ( ) × α n → o ( ) ( ) acc o → n N n N n = = ( ) ( ) × α o → n ( ) ( ) acc n → o N o N o Which gives as condition for ( ) ( ) acc o → n − β U n = e the acceptance rule: ( ) ( ) acc n → o − β U o e Understanding Molecular Simulation

Understanding Molecular Simulation

Metropolis et al. ( ) ( ) acc o → n − β U n = e Many acceptance ( ) ( ) acc n → o rules that satisfy: − β U o e Metropolis et al. introduced: ( ) = min 0, e − β Δ U ( ) ( ) = min 0, e ( ) − U o ( ) ⎡ ⎤ − β U n acc o → n ⎣ ⎦ If: Δ U < 0 acc( o → n ) = 1 accept the move If: acc( o → n ) = e − β Δ U Δ U > 0 draw a uniform random number [0;1] and accept the new configuration if: Δ U o → n ranf < e − β Δ U Understanding Molecular Simulation

3.Monte Carlo Simulation 3.4.2 Particle selection Understanding Molecular Simulation

Questions • How can we prove that this scheme generates the desired distribution of configurations? • Why make a random selection of the particle to be displaced? • Why do we need to take the old configuration again? • How large should we take: delx ? Understanding Molecular Simulation

Detailed Balance Understanding Molecular Simulation

3.Monte Carlo Simulation 3.4.3 Selecting the old configuration Understanding Molecular Simulation

Questions • How can we prove that this scheme generates the desired distribution of configurations? • Why make a random selection of the particle to be displaced? • Why do we need to take the old configuration again? • How large should we take: delx ? Understanding Molecular Simulation

Understanding Molecular Simulation

Mathematical Transition probability ( ) = α o → n ( ) × acc o → n ( ) π o → n from o → n: As by definition we ( ) ∑ make a transition: π o → n = 1 n The probability we do not make a move: ( ) = 1 − ( ) ∑ π o → o π o → n n ≠ 0 This term ≠ 0 Understanding Molecular Simulation

Model Let us take a spin system: (with energy U ↑ = +1 and U ↓ = -1) ( ) = e Probability to find ↑ : ( ) − β U ↑ P ↑ A possible configuration: If we do not keep the old configuration: (independent of the temperature) Understanding Molecular Simulation