Monte Carlo simulations of a semi-flexible polymer chain: a first glance

Monte Carlo simulations of a semi-flexible polymer chain: a first glance

D. T. Seaton; S. J. Mitchell; D. P. Landau

Center for Simulational Physics, University of Georgia, Athens, GA 30602, USA

ABSTRACT

We present preliminary results for Monte Carlo simulations of a three dimensional semi-flexible polymer chain with continuous monomer positions. In these simulations, standard Metropolis Monte Carlo methods are used to examine the basic properties of the model, such as equilibration configurations, overall size, and transition temperatures.

Keywords: Monte Carlo simulations; Polymer chain

I. INTRODUCTION

There are many examples of polymer systems, ranging in complexity from strands of DNA and protein molecules to examples of simple plastics. Regardless of the level of complexity, the three dimensional structure of a polymer chain will determine many of its most important properties. Monte Carlo methods have shown great promise in studying polymer conformation problems [1]. In particular, a wide variety of Monte Carlo methods have been applied to single chain polymer systems that include flexibility[2-4]. Many of these types of studies consider the bond fluctuation model [5], in which the monomer positions are confined to a lattice. Through the use of a polymer model with continuous monomer positions, we hope to add more realism to the study of general polymer conformation problems. In this preliminary work, we concentrate on understanding the basic features of our polymer model, such as equilibration time, overall size of equilibrated configurations, and location of transition temperatures. In particular, we examine the transition between folded and unfolded states.

II. MODEL

The polymer model consists of a chain of N identical monomers, where each monomer position varies continuously in three dimensions. There are three types of interactions: a harmonic bond length interaction between nearest neighbor bonded monomers, a bond angle interaction between adjacent bonds, and a non-bonded interaction between all non-bonded monomers. The Hamiltonian for this model is

where N is the number of monomers, J_L is the bond length interaction parameter, l_i is the i^th bond length, where the preferred bond length is l_i = 1.0, and ^-1 is the sum over all N-1 bonds. In the second term, J_A is the bond angle interaction parameter, ^-2 is the sum over all N-2 bond angles, and q_i is the angle between adjacent bonds, which has a preferred bond angle of q_i = 180. In the third term, J_NB is the non-bonded interaction parameter, å_i,j is the sum over all non-bonded pairs of monomers, and f(r_ij) is the non-bonded interaction potential, where r_ij is the distance between two non-bonded monomers. The function f(r) represents a quasi-Lennard-Jones potential [6] and is intended to model variable solvent quality [7,8]. This function is given by

where the minimum of this function is at r = 2. Fig. 1 compares f(r) to a standard Lennard-Jones potential, showing that f(r) has a softer repulsion as monomers come closer together.

We work in reduced units, in which J_NB = 1, and J_L and J_A can be varied. For the current study we examine the behavior of the model for

Initial simulations indicated that when J_L is too small ( ~ 1), the monomers distribute in a bimodal fashion, where bonds are either very short or very long. The current choice of J_L = 30 appears to produce relatively stiff bonds. In regulating the bond angle stiffness, we set J_A = 30 to allow for a semi-flexible polymer chain with enough freedom to find globular states.

III. METHOD

In simulating the polymer model, standard Metropolis Monte Carlo methods [9] have been used. Changes are made to the polymer configurations using the following types of moves: monomer diffusion, reptation, crank-shaft, and random pivot. One Monte Carlo sweep (MCS) consists of N attempted monomer diffusion moves and a single attempt of each of the other types of moves. A monomer diffusion move randomly selects a monomer, and then displaces it a random distance within a small cubic volume. A reptation move randomly selects one monomer from either end of the chain and then reattachs it to the opposite end of the chain. A crank-shaft move randomly selects two non-bonded monomers and then rotates the portion of the chain between those monomers by a small random angle. A random pivot move randomly selects a single monomer and the portion of the connected chain on either side of the monomer. A random direction is then chosen and the selected portion of the chain is then rotated around this random axis by a small random angle.

To characterize the behavior of the polymer, we observe quantities such as the radius of gyration and the average energy per particle. These two quantities offer much insight into the general properties of the polymer model, such as configuration types and equilibration time. Energies, E/N, are calculated using Eq. 1, and the radius of gyration is calculated using

where

_i is the i^th monomer position, and _cm is the position of the center of mass. The radius of gyration represents the overall size of the polymer chain. For a given N, low values of R_g represent more compact folded configurations, and higher values of R_g represent unfolded configurations.

IV. SIMULATION RESULTS

Typical simulations ran between 10⁵ and 10⁶ MCS, when starting from a straight configuration. In each simulation, R_g and E/N are calculated after each MCS. For the following results, radius of gyration and energy per particle are reported every 10³ MCS. The chain size in these simulations was usually in the range of N = 25 to N = 80, however, simulations have been run for sizes up to N = 256. Typical configurations for a N = 40 chain can be seen in Fig. 2. The first chain represents a folded configuration at a low temperature (T = 1.0), and the second chain represents an unfolded configuration at a higher temperature (T = 5.0).

In Fig. 3 the radius of gyration and energy per particle are both plotted as functions of time for three different chain sizes, all at a low temperature (T = 1.0). Symbols are added every 10⁴ MCS, to help distinguish the three curves. In this figure, each chain relaxes into a folded state, and at this temperature, the relaxation time scale is much shorter than 2.5 ×10⁵ MCS.

Fig. 4 is similar to Fig. 3 but shows data at a slightly higher temperature (T = 2.0). The N = 40 curve has been removed for clarity. The principal feature in this figure is that fluctuations occur on time scales much larger than those at T = 1.0. We see that the fluctuations in Fig. 3 are quite rapid on this time scale, but for T = 2.0 in Fig. 4 fluctuations occur on much larger time scales. Another feature seen in this figure is that larger fluctuations occur for shorter chains.

Observing how the radius of gyration changes with respect to temperature allows us to gain insight into the "transition" between folded and ufolded states. In Fig. 5, three R_g time series are given for a polymer of size N = 25, each at a different temperature. For the low temperature, T = 1.0, the polymer reaches a folded state with no large fluctuations. For the highest temperature, T = 4.0, the polymer rapidly fluctuates between slightly folded states and unfolded states. For the temperature T = 2.5, which is near the apparent transition temperature for this chain size, the polymer slowly fluctuates between folded and unfolded states. We can also see the relationship between radius of gyration and temperature in Fig. 6, where the chain sizes are N = 25, 40, and 80. As chain size increases, the apparent transition temperature increases. The shifting and broadening of the transition region is thought to be related to the stiffnes of the chain. By exploring different values of J_A, we expect to achieve a better understanding of this apparent transition.

V. CONCLUSIONS

Using standard Monte Carlo methods, we have studied the basic properties of a semi-flexible polymer chain. Through the observation of such quantities as the radius of gyration and the energy per particle, we have been able to understand the basic physical behaviors of our polymer model. However, these standard Monte Carlo methods become inefficient for larger systems, creating a need for a faster sampling method. In the near future, we would like to increase the efficiency of our simulations through the application of Wang-Landau sampling [10-12] techniques. This will not only "speed up" our simulations, but it should also allow us to gain a better understanding of the model.

Acknowledgments

We thank W. Paul and F. Rampf for helpful discussions. This work is supported by NSF grants DMR-0307082 and DMR-0341874, and LLNL grant B551576.

Received on 23 December, 2005

• [1] See e.g., K. Binder and A. Milchev, J. Comput. Aided Mater. Des. 9, 33 (2002).
• [2] M. Stukan, V. Ivanov, A. Y. Grosber, W. Paul, and K. Binder, J. Chem. Phys. 118, 3392 (2003).
• [3] F. Rampf, W. Paul, and K. Binder, Europhys. Lett. 70, 628 (2005).
• [4] J. Martemyanova, M. Stukan, V. Ivanov, M. M¨uller, W. Paul, and K. Binder, J. Chem. Phys. 122, 174907 (2005).
• [5] I. Carmesin and K. Kremer, Macromolecules 21, 2819 (1988).
• [6] A. Milchev, W. Paul, and K. Binder, J. Chem. Phys. 99, 4786 (1993).
• [7] V. Ivanov, W. Paul, and K. Binder, J. Chem. Phys. 109, 5659 (1998).
• [8] V. Ivanov, M. Stukan, V. Vasilevskaya, and K. Binder, Macromol. Theory Simul. 9, 488 (2000).
• [9] D. P. Landau and K. Binder, A Guide to Monte Carlo Simulations in Statistical Physics (Cambridge University Press, New York, NY, 2000), ISBN 0-521-65366-5.
• [10] F. Wang and D. Landau, Phys. Rev. Lett. 86, 2050 (2000).
• [11] F. Wang and D. P. Landau, Phys. Rev. E, 64, 056101 (2001).
• [12] F. Wang and D. P. Landau, Braz. J. Phys., 34, 354 (2004).

Publication Dates

• Publication in this collection
  16 Oct 2006
• Date of issue
  Sept 2006

History

• Received
  23 Dec 2005