J. Phys. Chem. 1994, 98, 10303-10306
10303
Metropolis Monte Carlo Brownian Dynamics Simulation of the Ion Atmosphere Polarization around a Rodlike Polyion M. Yoshida Toppan Printing Company, Sugito-machi, Kitakatsushika-gun, Saitama 345, Japan
K. Kikuchi* Department of Chemistry, College of Arts and Sciences, The University of Tokyo, Komaba, Meguro-ku, Tokyo 153, Japan Received: March 8, 1994; In Final Form: August 1, 1994@
The Metropolis Monte Carlo (MC) method when applied to part of the degrees of freedom of a system is shown to turn into a Brownian dynamics simulation by solving the master equation which describes the sampling process in the Metropolis MC by the S2-expansion approximation. Using the method, Brownian dynamics simulations are carried out on model DNA solutions to determine relaxation times of the ion atmosphere polarization. The MC approach to the problem eliminates the distinction between “free” and “bound” counterions often assumed in analytical studies. Molecular interpretation of the relaxation of the counterion cloud based on a simple model is attempted. In the absence of a complete theory, computer simulation is a promising method to assist in the interpretation of the experimental data.
Introduction There have been presented many calculations of the dipole moment induced in the ion atmosphere of a rodlike polyion by external electric The theoretical description of the response of polyelectrolyte molecules in solution to applied electric fields is, however, an extremely complex and difficult problem even if simple models are assumed for the geometry of the polyions. This is because there are involved many fields: concentrations of small ions, the electrical potential, and the solvent velocity, to be determined as functions of space and time which are coupled with each other through essentially nonlinear equations.23 Consequently in many cases recourse has been had to somewhat ad hoc approximations such that counterions are classified into “free” and “bound” ions, only the latter contributing to the induced moment. Recently we have used computer simulation to approach the problem in a practical way.24’25We have carried out Metropolis Monte Carlo (MC) simulations on counterion distributions around rodlike polyions (actually oligoions) of various length and their deformation in electric fields. The calculations have provided some new insights into the behavior of small ions in external electric fields: the density of counterions in the immediate vicinity of the oligoion is depleted near the ends of the oligoion rod already in the absence of external electric fields, the induced dipole moment as calculated from the shift of the center of the distribution of all the counterions in a simulation volume reached a maximum and then decreased as the applied field strength increased. Olmsted et a1.26have also found the end effect on the local concentration of counterions at various axial positions along the polyion surface by the grand canonical MC method. The high-field behavior of the induced moment is due to progressive stripping of the ion atmosphere from the electrostatic influence of the polyion leaving fewer charges left to polarize in the field as predicted by Rau and Charney16who, however, treated the electric polarization of an ion atmosphere in the Debye-Huckel approximation. @
Abstract published in Advance ACS Abstracts, September 1, 1994.
0022-365419412098-10303$04.5OlO
The Metropolis MC method2’ has been conventionally regarded as applicable only to the study of equilibrium properties and unable generally to evaluate dynamic properties such as transport coefficients and time correlation functions. We have recently s h o ~ nthat ~ ~if we , ~restrict ~ ourselves to a part of the degrees of freedom of a system, e.g., Brownian particles embedded in a solvent continuum, the MC procedure amounts to simulating Brownian movement of the particles, thus justifying our application of the method to the study of such a nonequilibrium process as the response of polyelectrolytes to external electric field^.^^,^^ We further extended the method to include the effects of hydrodynamic interaction^.^^ Thus while in the previous paper^^^^^^ we have studied only the steady-state electric polarization properties of polyelectrolyte solutions, the Metropolis MC method allows us to simulate the relaxation process of the counterion polarization i t ~ e l f . ~ We *,~~ first prepare an equilibrium distribution of the counterions around the polyion and then follow the time course of the change of the distribution upon application of an electric field pulse to the solution. Finite polarization relaxation times of the ion atmosphere have often been postulated in theoretical formulations of transient electric birefringence or dichroism of axially symmetric macromolecules as a slow induced d i p ~ l e . ~ l -Their ~ ~ reliable determination has, however, been limited not only theoretically but also technically by nonideal electric field pulse shape. Recently P O r s ~ h k has e ~ ~succeeded to determine relaxation times of the counterion dynamics for rodlike DNA restriction fragments by analyzing electric dichroism rise curves measured with high time resolution by efficient numerical procedures. Without invoking any of the theoretical treatments of the dynamics of counterion motion,1,2,4,6,7,10,13,14,17,18 he explained the polarization process by an ion dissociation reaction model and the mechanism described the dependences of the observed relaxation times on the salt concentration and the electric field strength semiquantitatively. In this paper we first formulate the theoretical basis of our simulation method by, following van K a m ~ e n solving , ~ ~ the master equation which describes the sampling process in the 0 1994 American Chemical Society
10304 J. Phys. Chem., Vol. 98, No. 40, 1994
Yoshida and Kikuchi
Metropolis MC by the Q-expansion approximation. We then apply the method to the determination of relaxation times of the ion atmosphere polarization around some model DNA molecules. Molecular interpretation of the relaxation of the counterion cloud based on a simple model is attempted. It is discussed that in the absence of a complete theory, computer simulation can be most useful for studying dynamics of the ion atmosphere polarization. Method
The Metropolis MC27was originally developed as a method suited to electronic computers for calculating statistical mechanical configurational integrals. Since the MC sampling is a Markovian process, if we introduce a “time” scale t which actually labels the order of subsequent configurations X , the “dynamic” evolution of the probability distribution function P(X,t) is governed by the master equation a p w at
-= J {W(XlX’) P(X’,t) - W(X’1X) P(X,t)} dx‘
(1)
where W(X1X’) is the transition probability per unit “time” for a transition from configuration X‘ to X . The master equation may be solved approximately by an expansion method in powers of a parameter S2, the size of the system.35 Given that W(X(X’)has the canonical form
w(x~x’) = Qo(x’;r) + Q-’Q,(x’;r)
+ Q-2Q2(x’;r) + ...
(2)
with X‘ = X’/Q and r = X - X’, the jump moments are defined by
(3)
a,,(x) = Jr”QA(x;r) d r
When al,o = 0, the Q expansion yields as the lowest approximation a nonlinear Fokker-Planck equation:
with z = Q-2t. In the Metropolis MC, the reciprocal of the maximum displacement allowed for an MC move during a “time” interval At may be taken as the parameter 8. Then the transition probability has the canonical form with Q&r)
Irl I1
= 1/(2At)
for
=0
Irl > 1
for
and, supposing, for example, the derivative of the potential energy of the system U ( x ) with respect to x is negative, Le., U’(x) < 0, with Q1(x;r)
= -[1/(2kBTAt)]U’(x)r =O
for
for
-1 Ir < 0
r