微信客服：xiaoxionga100
微信客服：ITCS521
MATLAB代写-BD1
时间：2021-10-04
Instructions BD1 Exercise 1: Random numbers are usually generated via a sequence of large positive integers Xi, each generated from the previous one by some simple operation (as multiplication) conducted in finite modulus arithmetic, 1 mod( )i iX aX M+ = (1), mod( ) ( / )X M X M X M= − ⋅ int (2), a and M are both large positive integers. The generated by (2) integer thus lies in (0,M-1), the desired random numbers ξi = iX /M. The period of a random number generator is the number of random numbers generated before the sequence begins to repeat itself. What is the (maximal) period of this random number generator? Exercise 2: Let ξ be a number of points on a side of a dice. What is an average number which can be produced after throw the dice once? Exercise 3: Prove that the normal distribution with mean a and variance σ2, 2 2 1 ( )( , ) exp 22 x aN a σ σπσ − = − , is properly normalized, i.e., its integral over the whole sample space is equal to unity. Exercise 4: Calculate the characteristic function of the random variable distributed according to the normal distribution N(a,σ) with mean a and dispersion σ2. Modelling Exercise 1: In this exercise we will check the random number generator which we will use later to simulate Brownian motion. You have to use MatLab for this exercise. a) Using Matlab random-numbers generator, generate (some amount of) random numbers x uniformly in the interval (0,1) and plot the corresponding distribution function (histogram). Show that the generated numbers are (indeed) uniformly distributed. Check the quality (the deviations from the uniformity, the fluctuations, for different amount of just generated random numbers) of the distribution upon changing the amount of generated numbers. Calculate the standard deviation for different amount of realizations; b) Using these just obtained numbers, generate two sets of normally distributed random numbers, one set with fixed mean
the Central Limit Theorem for this purpose. Plot the histograms, fit them with the corresponding Gaussian
functions – on the same plot - and check the quality of this fit. The simplest check is with the eye inspection.
The goodness of fit can also be calculated by summarizing the (absolute) discrepancies between the
numerical values and the values expected from the Gaussian fit.
c) In 2D, generate L (L can vary from 1 to ~10000, and can be 1, 10, 100, 1000, 10000) different random
walks ON A SQUARE LATTICE. “Random” means not self-avoiding, i.e. the lattice walker walks with equal
probability in one of the four different directions. Each random walk consists of N (N=100) unit (length =1)
steps. Study the statistical properties of these random walks, i.e. calculate the mean-squared end-to-end
distance
theoretical prediction, see Lecture 1 notes. How does this depend on L? Calculate error bars for each L.
Instructions BD2
Exercise 5: The X-axis projection x(t) of the Brownian particle radius-vector r(t) is time-dependent
random variable (the same is valid for Y and Z projections). Let G(t,t’)=
function of the stochastic variable x(t) at two different time moments, t and t’.
a) Why G(t,t’)=G(t-t’)?
b) Show that G(t) has a maximum at t=0.
c) If ( )G ω is the Fourier transform of G(t), for the inverse Fourier transform we have
1( ) ( )
2
i tG t G e dωω ω
π
+∞
−
−∞
= ∫
Prove that for the (double) Fourier transform ( , )G ω ω′ the following is valid
( , ) 2 ( ) ( )G Gω ω πδ ω ω ω′ ′= +
Use the Fourier representation of δ-function here.
d) prove the Wiener-Khintchin theorem
1( ) ( )cos
2
G t D t dω ω ω
π
+∞
−∞
= ∫
where ( ) 2 ( )D Gω ω= .
Exercise 6: Find the characteristic time of momentum relaxation for a Brownian particle of mass m
and friction coefficient ζ. Estimate the times when the friction term in the equation of motion is much
higher than the inertial term. Hint: use the equation of motion )()( trtrm
ζ−= and try to find the solution
in the form of )exp( tω− .
Modelling Exercise 2: Brownian-dynamics simulation of a polymer chain.
Your task in this exercise is to start to build a MatLab program which simulates a simple (so-called freely-
joint) model of a short polymer chain: N = 17 spheres, all with friction coefficients ς , connected by the
Nb=16=N-1 of the so-called Hookean springs (chemical bonds li= 1i i+−r r , i=1,...,Nb) and having spring
constant K (which is a number, in some dimensionless units), see below.
ζ
ζ
In this case the potential energy of the chain is very simple, and is given by the sum of potential energies of
each individual spring,
1 1
2 2
1
1 1
1 / 2 ( ) 1 / 2
N N
i i i
i i
U K K l
− −
+
= =
= − =∑ ∑r r
where li is the instant (i.e., fluctuating) length of the i-th spring, and 2il is obviously its square. We do not
take into account the more complicated interactions between the spheres (the so-called excluded-volume
interactions, which can be modelled by the Lennard-Jones (LJ) potential, or Coulombic interactions if
spheres are charged), so, no other potential forces. We study the Brownian motion of this polymer chain in
the thermostat with temperature T. The polymer chain diffuses in infinite space, no boundary conditions are
used here.
a) Please refresh your knowledge about ensembles in statistical physics (for example, from D. Chandler or A.
Carter books). Assume the NVT (canonical) ensemble conditions. The distribution function 1( , ,..., )Nρ 2l l l for
the bond vectors li to have some length, | |i l=l , and
2
il =l2 is proportional to the Boltzmann factor,
1 2exp( ( , ,..., ) / )N BU k T− l l l . Show that this distribution function can be factrorized, 1
1
( , ,..., ) ( )
b
b
N
N i
i
ρ ρ
=
=∏2l l l l .
b) Now, take into account that ( )iρ l does not depend on the orientation of the ith-bond vector
, , ,( , , )i x i y i zi l l l=l , but on its length instead. Using the spherical coordinates, we have in this case
2
, , ,( ) ( ) sini i x i y i z i i idl dl dl l l dl d dρ ρ θ θ ϕ=l ,
and natural normalization condition,
2
2 2
0 0 0 0
( ) sin 4 ( ) 1i i i i i il l dl d d l l dl
π π
ρ θ θ ϕ π ρ
∞ ∞
= =∫ ∫ ∫ ∫
If so, any average property ( )iA l which depends on the bond length il only, can be calculated using these
spherical coordinates and (properly normalized!) distribution function 2( ) exp( / 2 )i i Bl C Kl k Tρ = − which
depends on the bond length only,
2
0
( ) ( ) ( )i i i i iA l l l A l dlρ
∞
= ∫
The constant C should be determined from the normalization. Write down this properly normalized
expression for ( )ilρ .
c) Having this distribution function, calculate the squared average bond length 2il for the spring in your
chain. What is the value of this length at T=0?
d) Calculate the characteristic diffusion time τ in this problem. What you have as physical quantities here
are:
- friction coefficient ς ;
- average bond length 2 0il l= , which can be served as a unit of length in your simulations;
- thermostat temperature T and the corresponding thermal energy kBT.
Using analysis of dimensions create the unit of time. What is the physical meaning of this time? ALL THE SO
CHOSEN UNITS (FRICTION, LENGTH, TEMPERATURE, ENERGY AND TIME) CAN BE PUT EQUAL TO UNITY.
e) Having the potential energy, calculate the total potential force acting on a particle number i. Show the
forces for the terminal particles (number 1 and number 17), and for the particle inside the chain.
f) Write down the corresponding Langevin equation for the velocity of particle i at time t, in a diffusive (non-
inertial) limit for this case, as function of the remaining forces (i.e., potential force and random force). Write
down the corresponding finite-difference scheme for each Cartesian component of the particle’s position.
Note the difference between the equation for a particle inside the chain, and for two terminal particles.
g) Choose the proper value of the integration time step Δt in the problem, having the value of the
characteristic physical time of the problem (part c)). Program the non-inertial Langevin equation for the
time-dependent positions of all the particles in the finite-difference approximation, starting from some
initial time, and performing M (large number, ~ 104 - 105) integration steps.
h) All the fundamental units in this problem:
- friction coefficient ς
- thermal energy Bk T
- length 2 0il l=
can be taken equal to unity in your simulations. All other physical properties can be measured in these
units. What is the value of the dimensionless spring constant in your problem? What is the value of the
dimensionless time step?
i) Perform simulations without random forces. What will be the final configuration of the chain? Explain the
physical meaning of the result.
Instructions BD3
Modelling Exercise 3: Brownian-dynamics simulation of a polymer chain, continued.
a) We continue with the simulations of the chain. It is time now to have random forces in the code!
What is the value of the random force dispersion (in dimensionless units) which should be used? For each
particle, and for each Cartesian coordinate, generate random forces from the normal distribution with zero
mean and proper dispersion (you have calculated this dispersion earlier).
b) Read once more paragraphs 2.6 and 2.7 of the manual. Place the particles initially along the line,
i.e., create the absolutely stretched polymer chain. Start simulations, check stability.
c) Perform simulations with different time steps, and check the stability.
d) In simulations, calculate for each bond its instant length, plot (for few selected bonds) its value as
function of the simulated time, see fluctuations. Calculate the time-averaged squared bond length, also
averaged for all the bonds in the chain. Compare the simulated values with the theoretical predictions (see
Exercise 2).
e) Save the positions of the terminal (number 1 and number 17) polymer beads for different time
moments, and use them to calculate the mean squared end-to-end chain length as a function of time, and
averaged over time. Calculate the time average and standard deviation, for different simulation times. For
this model of a polymer chain of N bonds (N random walks of some length l which is fluctuating, but is
constant on average!) the mean-squared end-to-end distance is easily calculated (see Lecture 1 and Exercise
1). Compare your simulation results with these theoretical predictions.
Instructions BD4
Exercise 7: prove that Maxwell distribution
3/2 2
( ) exp , ( ) 1
2 2B B
m mf f d
k T k Tπ
= − =
∫
vv v v
is the solution of the Fokker-Planck equation
2 2
2 2( ) ( )B
f f f fZ f k T f
t m
γξ γ∂ ∂ ∂ ∂ ∂ ∂+ = + = +
∂ ∂ ∂ ∂ ∂ ∂
v v v
r v v v v
in equilibrium, i.e., when ( , , ) ( )f t f v=r v , does not depend on time and coordinate.
Exercise 8: Using the Green function
[ ]
2
0
0 0 3/2
00
( )( , | , ) exp
4 ( )4 ( )
Nn t t
D t tD t tπ
−
= − −−
r rr r
for Einstein diffusion equation, calculate the first and second moments
2
0 0 00, ( ) 6 ( )D t t< − >= < − >= −r r r r of the particle’s displacement r-r0 from the formal definition
of the averaging <>.
Exercise 9: The Einstein equation
= 〈Δ2()〉 = 6
relates macroscopic transport diffusion coefficient D with microscopic information on the mean-square
displacement of a particle during its Brownian diffusion.
a)What is the physical dimension of this coefficient? Which other examples of such transport coefficients
can you mention? Name two at least.
b) How this equation is changed in the N-dimensional space?
c) Consider ballistic (inertial) motion of a particle with a constant velocity v. What will be the expression for
MSD in this case? Compare the power-law dependences (so-called scaling laws) of MSD on time for ballistic
and diffusive motion and explain the difference.
Exercise 10. The normal modes normalization
As discussed in the lecture, the transition from the projections of the bond vectors to normal modes can be
done using the following linear transformation between two Nb-dimensional vectors; vector
}{ 1 2,, ..., bNu u u u= and vector }{ 1 2, ,... bNq q q q= are connected as
1
q u
u q−
=
=
B
B
or, using the i and p indexes,
1 1
2sin sin
1 ( 1) 1
b bN N
j p p
p pb b b
jp jpu A q q
N N N
π π
= =
= =
+ + +∑ ∑
2 sin , , 1,...
1 ( 1)jp bb b
jpB j p N
N N
π
= =
+ +
Show that with this choice of the normalization constant A, 2( ) 1jj =B . Without proof,
2( ) 0,jk j k= ≠B .
These two facts simply mean that B2 is the unit matrix, B2 = I, and B=B-1.
Modelling Exercise 4: Brownian-dynamics simulation of a polymer chain, normal modes analysis.
Use the polymer chain which you have already created in Exercises 2 and 3. Start your Brownian dynamics
simulations. In the course of the simulation create and save to the disk Nb bond vectors (you probably did it
already),
1 , 1,...i i i bb r r i N+= − =
, and their X projections, 1 , 1,...i i i bu x x i N+= − = . For each time moment, i.e.
every integration step, calculate the normal modes
1
2( ) sin ( ), 1,... , 1,...
( 1) 1
bN
p j b b
jb b
jpq t u t j N p N
N N
π
=
= = =
+ +∑ (1)
a) Show that the linear transformation (1) is orthonormal, by direct simulation and comparison of
2 2 2 2
1 2| | ... bNq q q q= + + +
and 2 2 2 21 2| | ... bNu u u u= + + +
b) From the lecture material we know that each mode relaxes independently,
/
2
( ) (0)
,
2 (1 cos ) 14 sin
2
pt
p p
p p
pp b
q t q e
pkkK k NK
τ
ς ς πτ
−=
= = =
− +
For each mode with a number p, calculate the following autocorrelation function:
/
2
(0) ( )
pp p t
p
q q t
e
q
τ−=
directly in your simulations. The angular brackets, as usual, denote the time averaging. Show that, indeed:
- the modes relax exponentially,
- the relaxation time pτ is equal to what was predicted analytically.
- Plot these autocorrelation functions. Explain how the relaxation times were extracted from the fits,
calculate the relaxation times and compare with theory.
学霸联盟
