Flexible polymer chain dynamics

Download Report

Transcript Flexible polymer chain dynamics

Flexible polymer chain
dynamics
Yong-Gu Lee
Bead-spring chain model
Pg 29
A double-stranded DNA molecule is described by a bead-spring chain model composed of
Nb beads of hydrodynamic radius a = 77 nm connected by Ns = Nb − 1 entropic springs.
Each bead represents a DNA segment of 4850 base pairs, i.e., Nb = 11 corresponds to a
stained -DNA, which has a contour length of 22 μm and radius of gyration of 730 nm.
The springs connecting the beads obey a worm-like chain force law. Here, bk is the
Kuhn length for DNA and Nk,s = 20 is the number of Kuhn length per spring.
Increased from natural 16.3um by TOTO-1 dye
Pg 5
Yu Zhang, “Brownian dynamics simulation of DNA in complex geometries,”
PhD thesis, University of Wisconsin – Madison, 2011
Diffusion coefficient
• The projection of the Brownian motion in the x-y
plane was tracked using video fluorescent
microscopy and the increase in the mean square
displacments with time t were recorded for an
ensemble of molecule paths.
x
2
 y
2
 2Dt
Douglas E. Smith, Thomas T. Perkins, and Steven Chu, “Dynamical scaling of DNA diffusion coefficents,” Macromolecules, 1996, pp. 1372-1373
Zimm model for the radius of
Gyration
>> 22^(3/5)
ans =
6.3893
Persistent length
• Informally, for pieces of the polymer that are shorter than the
persistence length, the molecule behaves rather like a flexible
elastic rod, while for pieces of the polymer that are much longer
than the persistence length, the properties can only be described
statistically, like a three-dimensional random walk.
• Formally, the persistence length, P, is defined as the length over
which correlations in the direction of the tangent are lost. In a
more chemical based manner it can also be defined as the
average sum of the projections of all bonds j ≥ i on bond i in an
indefinitely long chain.
• Let us define the angle θ between a vector that is tangent to the
polymer at position 0 (zero) and a tangent vector at a distance L
away from position 0, along the contour of the chain. It can be
shown that the expectation value of the cosine of the angle falls
off exponentially with distance
cos   e
http://en.wikipedia.org/wiki/Persistence_length
pL
Persistent length
• For a chain possessing some degree of stiffness, the
tangents to two segments of the chain contour will
tend to be pointed in the same direction, provided
that the segments are sufficiently close to each
other, In other words, the local contour tends to
persist in a given direction (segment directions are
correlated). If the chain contour is represented by a
virtual chain comprising a string of segment vectors,
each of length b, the stiffness of the chain can be
represented in a quantitative fashion as the sum ,
𝑃𝑛
of the average projections of each segment vectors
on the first segment vector.
n
Pn  b  cos 
i
i 1
Paul J. Haggerman, “Flexibility of DNA,” Annual Review of Biophysics, Vol 17. pp. 265-286, 1988
Radius of gyration
• It is defined as the root-mean-square distance of
the collection of atoms, or groups, from their
common center of gravity.
s   n  1
2
1
n
s
2
i
i0
s i is the distance of atom i from the center of gravity of the chain in a specifie d configuration.
Paul J. Flory, “Statistical mechanics of chain molecules,” pp. 4, Hanser Publishers, 1988
Kuhn length
• The Kuhn length is a theoretical treatment, developed
by Werner Kuhn, in which a real polymer chain is
considered as a collection of N Kuhn segments each
with a Kuhn length b. Each Kuhn segment can be
thought of as if they are freely jointed with each other.
Each segment in a freely jointed chain can randomly
orient in any direction without the influence of any
forces, independent of the directions taken by other
segments. Instead of considering a real chain consisting
of n bonds and with fixed bond angles, torsion angles,
and bond lengths, Kuhn considered an equivalent ideal
chain with N connected segments, now called Kuhn
segments, that can orient in any random direction.
Kramers freely jointed bead-rod
chain model
P. M. Saville and E. M. Sevick, “collision of a field-driven polymer with a finite-sized obstacle: A
Brownian dynamics simulation,” Macromolecules, Vol. 32, pp. 892-899, 1999
Length based sorting
P. M. Saville and E. M. Sevick, “collision of a field-driven polymer with a finite-sized obstacle: A
Brownian dynamics simulation,” Macromolecules, Vol. 32, pp. 892-899, 1999
Experimental parameters for Kramers
freely joined bead-rod chain model
• As in previous work, fluorescently tagged Lambda-DNA is modeled using the
Kramers freely jointed bead-rod chain model. This model has demonstrated
quantitative agreement with experiments in single molecule flow studies of DNA
even though the model does not reproduce the correct bending forces or
effective entropic spring force of the wormlike chain or Kratky-Porod model. To
simulate the Kratky-Porod model, however, in the present application, 10 subKuhn step level discretization points are needed in order to produce the correct
persistence length in a discrete model appropriate for Brownian dynamics. Thus,
the number of dynamics variables becomes prohibitively large for the long time,
large ensemble simulations considered in this work (see Ensemble Size section
below). Therefore, we choose to focus on the qualitative comparisons to the
experimental data with our simulations. As we shall discuss below, errors
associated with ensemble size in these comparisons are at least as important as
those associated with the shortcomings of the model. In the Kramers model, the
polyelectrolyte is divided into N beads, at which the mass and hydrodynamic
drag are concentrated, connected by N - 1 massless rods. The entire molecule
has a contour length of 21 um and is modeled using 150 beads, yielding a Kuhn
step size, bk, of 0.142 um. The drag is assumed to be isotropic, and
>> 150* 0.142
hydrodynamic interactions (HI) between segments of the chain are ignored.
ans =
21.3000
Nerayo P. Teclemariam, Victor A. Beck, Eric S. G. Shaqfeh,,§ and Susan J. Muller,
“Dynamics of DNA Polymers in Post Arrays: Comparison of Single Molecule Experiments
and Simulations,” Macromolecules, Vol. 40, pp. 3848-3859, 2007
Polymer chain dynamics
T he B row nian dynam ics algorithm is im ple m ented by integrating the stochastic dif ferential equations to first order in  t
F or a tim e step  t the displacem ent of be ad  consists of tw o parts:
r  t   t   r  t   t    r .
T he first term r  t   t  is the position w hich w ould have been reached in the absence of constratints. T he second term
 r is the correction to r  t   t  due to the constraint forces.
T he bead positi on after an unconstrained m ove, r  t   t  ,
is calculated in a straightforw ard m anne r:
r  t   t   r  t    κ  r  t    t 


 W 
F

t
κ : transpose of the velocity gradient ten sor
 = 2 k BT 
 2t 
 a   i B i u i  t 
   i
 r   
 W is chosen from a G aussian distribtuion cha racterized by the first and second m om en ts:
 W  t   0,
 W  t   W   t     δ  t.
a is the rod like link length.
Tony W. Liu. “Flexible polymer chain dynamics and rheological properties in steady flows,” Journal of Chemical Physics, Vol 90. No. 5826, 1989
•
•
•
•
•
•
•
•
•
•
•
/* Rouse matrix *** Aij*/
MatrixNminus1d* getRouseMatrix()
/* Modified Rouse matrix *** AHatij*/
MatrixNminus1d* getModifiedRouseMatrix(Matrix3Nminus1d& umatrix)
/* Kramers matrix *** Cij*/
MatrixNminus1d* getKramersMatrix()
/* Kramers tensors matrix is second rank tensor *** Kij*/
Matrix3d* getKramersTensorMatrix(int j, int k, Matrix3Nminus1d& umatrix)
/* B_bar_ij*/
MatrixNminus1_Nd* getB_BarMatrix()
const double physical_coefficient = 1.0 / a * sqrt(2.0 * boltzman_constant *
temperature/friction_coefficient);
• Matrix3d kappa;//Transpose of the velocity gadient tensor del(v)
int _tmain(int argc, _TCHAR* argv[])
{
Matrix3d kappa;//Transpose of the velocity gadient tensor del(v)
kappa << 0,0,0, 0,0,0, 0,0,0;
const double maxt = 100;
const double deltat = double(1)/double(maxt)*1000.0e-7;//1000e-12;
////The whole simulation duration is 1000 pico second
const double boltzman_constant=1.3806503e-23;
//1.3806503 × 10-23 m2 kg s-2 K-1
const double temperature = 298.0;
const double friction_coefficient = 6.0 * M_PI * 8.9e-4 * 10.0e-9;
// 6 PI Mew R, 10 nm radius bead
// The dynamic viscosity of water is 8.90 × 10 4 Pa·s at about 25 °C.
const double a = 100.0e-9; // rod length
const double alpha = sqrt(2.0 * boltzman_constant * temperature *
friction_coefficient);
const double physical_coefficient = 1.0 / a * sqrt(2.0 * boltzman_constant
* temperature/friction_coefficient);
VectorNminus1d ux, uy, uz;
for( int i = 0; i < N-1; ++ i )
{
ux(i) = 0;
uy(i) = -1;
uz(i) = 0;
}
Matrix3Nminus1d umatrix;
Matrix3Nminus1d u_primematrix;
Matrix3Nd rmatrix;
Matrix3Nd r_primematrix;
Matrix3Nd compensating_matrix;
MatrixNminus1d* RouseMatrix;
RouseMatrix = getRouseMatrix();
MatrixNminus1_Nd* B_Bar_Matrix;
B_Bar_Matrix = getB_BarMatrix();
MatrixN_Nminus1d* B_Matrix;
B_Matrix = getB_Matrix();
for( int m = 0; m < N-1; ++m )
{
umatrix(0,m) = ux(m);
umatrix(1,m) = uy(m);
umatrix(2,m) = uz(m);
}
rmatrix = (a*(*B_Matrix) * umatrix.transpose()).transpose();
// 3 x N
//std::cout << rmatrix << std::endl;
//print the initial coordinates of the beads
ofstream file;
file.open("E:/temp/beads.xyz");
//open a file
double mu = 0.0;
double sigma = 1.0;
std::tr1::normal_distribution<double> normal(mu, sigma);
for( int t = 0; t < maxt; ++t)
{
for( int m = 0; m < N; ++m)
{
_Vector3d *delta_w;
delta_w = get_delta_w(normal, deltat);
r_primematrix.col(m) = rmatrix.col(m) + kappa *
rmatrix.col(m) * deltat + alpha / friction_coefficient * (*delta_w);
delete delta_w;
}//m
//std::cout << r_primematrix << std::endl; //print the
unconstrained movement of the beads (Eq. 65)
for( int m = 0; m < N-1; ++m)
{
u_primematrix.col(m) = (r_primematrix.col(m+1)r_primematrix.col(m))/a;
umatrix.col(m) = (rmatrix.col(m+1)-rmatrix.col(m))/a;
}
VectorNminus1d gamma;
{
u_primematrix.col(m) = (r_primematrix.col(m+1)r_primematrix.col(m))/a;
umatrix.col(m) = (rmatrix.col(m+1)-rmatrix.col(m))/a;
}
VectorNminus1d gamma;
get_gamma(gamma, umatrix, u_primematrix, deltat); //Eq. 70
for( int n = 0; n < N; ++n )
{
rmatrix.col(n) = r_primematrix.col(n);
for( int m = 0; m < N-1; ++m)
{
rmatrix.col(n) += -2.0 * deltat /
friction_coefficient * a *
gamma(m)*(*B_Bar_Matrix).coeff(m,n)*umatrix.col(m);
}
}
//std::cout << "t = " << t << "-------------" << std::endl;
_Vector3d babo;
file << N << "\n";//No. of files
file << "beads.xyz time [ps] " << float(t) << "\n";//Time
stamp
for( int m = 0; m < N; ++m)
{
}
//std::cout << "t = " << t << "-------------" << std::endl;
_Vector3d babo;
file << N << "\n";//No. of files
file << "beads.xyz time [ps] " << float(t) << "\n";//Time
stamp
for( int m = 0; m < N; ++m)
{
//std::cout << "\tm = " << m << std::endl;
//std::cout << "\t" << rmatrix.col(m).coeff(0) << ' ' <<
rmatrix.col(m).coeff(1) << ' ' << rmatrix.col(m).coeff(2) << ' ' <<
std::endl;
file << "beads " << rmatrix.col(m).coeff(0) << ' ' <<
rmatrix.col(m).coeff(1) << ' ' << rmatrix.col(m).coeff(2) << ' ' <<
"\n";
}
}//t
file.close();
return 0;
}