Monte Carlo simulations of a polymer chain conformation. The effectiveness of local moves algorithms and estimation of entropy.  
Jump to Full Text  
MedLine Citation:

PMID: 23765038 Owner: NLM Status: Publisher 
Abstract/OtherAbstract:

A linear chain on a simple cubic lattice was simulated by the Metropolis Monte Carlo method using a combination of local and nonlocal chain modifications. Kinkjump, crankshaft, reptation and endsegment moves were used for local changes of the chain conformation, while for nonlocal chain rearrangements the "cutandpaste" algorithm was employed. The statistics of local micromodifications was examined. An approximate method for estimating the conformational entropy of a polymer chain, based on the efficiency of the kinkjump motion respecting chain continuity and excluded volume constraints, was proposed. The method was tested by calculating the conformational entropy of the undisturbed chain, the chain under tension and in different solvent conditions (athermal, theta and poor) and also of the chain confined in a slit. The results of these test calculations are qualitatively consistent with expectations. Moreover, the obtained values of the conformational entropy of self avoiding chain with ends fixed over different separations, agree very well with the available literature data. 
Authors:

Agnieszka Mańka; Waldemar Nowicki; Grażyna Nowicka 
Related Documents
:

23689028  On a staggered ifem approach to account for friction in compression testing of soft mat... 7064018  Serum ferritin and mean corpuscular volume as predictors of bone marrow iron stores. 21428508  Blind singlechannel suppression of late reverberation based on perceptual reverberatio... 23505648  A note on the reproduction of geochelone carbonaria and g. denticulata (testudines, tes... 16405468  Vertically integrated flows, discharge potential, and the dupuitforchheimer approximat... 23548298  Development of a genetic evaluation for body condition score for canadian holsteins. 
Publication Detail:

Type: JOURNAL ARTICLE Date: 2013614 
Journal Detail:

Title: Journal of molecular modeling Volume:  ISSN: 09485023 ISO Abbreviation: J Mol Model Publication Date: 2013 Jun 
Date Detail:

Created Date: 2013614 Completed Date:  Revised Date:  
Medline Journal Info:

Nlm Unique ID: 9806569 Medline TA: J Mol Model Country:  
Other Details:

Languages: ENG Pagination:  Citation Subset:  
Affiliation:

Faculty of Chemistry, Adam Mickiewicz University, Umultowska 89b, PL 61714, Poznań, Poland. 
Export Citation:

APA/MLA Format Download EndNote Download BibTex 
MeSH Terms  
Descriptor/Qualifier:

Full Text  
Journal Information Journal ID (nlmta): J Mol Model Journal ID (isoabbrev): J Mol Model ISSN: 16102940 ISSN: 09485023 Publisher: Springer Berlin Heidelberg, Berlin/Heidelberg 
Article Information Download PDF © The Author(s) 2013 OpenAccess: Received Day: 28 Month: 2 Year: 2013 Accepted Day: 30 Month: 4 Year: 2013 Electronic publication date: Day: 14 Month: 6 Year: 2013 pmcrelease publication date: Day: 14 Month: 6 Year: 2013 Print publication date: Year: 2013 Volume: 19First Page: 3659 Last Page: 3670 PubMed Id: 23765038 ID: 3744652 Publisher Id: 1875 DOI: 10.1007/s008940131875z 
Monte Carlo simulations of a polymer chain conformation. The effectiveness of local moves algorithms and estimation of entropy  
Agnieszka MańkaAff1  
Waldemar NowickiAff1 
Address: +48618217376 +4861618291505 gwnow@amu.edu.pl http://www.staff.amu.edu.pl/∼gwnow 
Grażyna NowickaAff1  
Faculty of Chemistry, Adam Mickiewicz University, Umultowska 89b, PL 61714 Poznań, Poland 
The knowledge of the Helmholtz free energy A (or Gibbs free energy G) of a polymer is essential for understanding many complex processes involving macromolecules, as for instance, the adsorption or grafting of polymer chains on interfaces [^{1}–^{6}], the formation of polymer bridges [^{7}, ^{8}] or protective layers in processes of flocculation or steric stabilization of colloidal suspensions [^{8}, ^{9}], the formation of secondary and tertiary structure of proteins [^{10}–^{12}], the release of the genetic material from the viral capsids [^{13}–^{16}] and many others. However, in such complex processes as those mentioned above, it is not possible to analytically compute the changes in the free energy and, particularly, in its component—the conformational entropy. This problem can be overcome by the use of simulation methods which provide the ability to study complex systems in great detail [^{17}–^{19}]. The application of these methods requires a construction of a simulation model which accurately represents the processes being modeled. The polymer model appropriate for investigation of the abovementioned phenomena should include the excluded volume and energetic interactions. It should also account for deformation of a polymer molecule from its equilibrium conformation in the presence of external constraints (e.g., a macromolecule confined in a nanopore or attached to an impenetrable surface) or/and forces (e.g., a macromolecule bridging two particles undergoing Brownian motion). These requirements impose certain restrictions on the choice of the method and algorithm for the estimation of free energy and conformational entropy of a polymer molecule. Aiming at the maximum possible simplicity of the model and taking into account the features it should include, a linear homopolymer chain simulated on a lattice by the Metropolis Monte Carlo (MMC) method [^{20}, ^{21}] was chosen for the study. Constant length of the chain and fixed positions of terminals were assumed. The latter assumption arises from the adopted strategy for the examination of relationship between the Helmholtz free energy of a polymer molecule and the deformation extent, based on simulation of conformation of a chain having arbitrarily imposed fixed distance between its ends.
In this study an attempt was made to develop a method permitting estimation of the Helmholtz free energy of a polymer chain on the basis of probabilities of micromodifications of chain conformation in the MCC simulation. The adopted elementary moves of the chain followed the generalized VerdierStockmayer algorithm [^{21}–^{23}]. In this algorithm a number of local moves can be performed depending on a local chain conformation. The MC simulations using VerdierStockmayer algorithm were shown to reproduce some real kinetic properties of macromolecules [^{24}, ^{25}]. It was proven, however, that all local chain lengthconserving algorithms are nonergodic [^{26}, ^{27}]. Although, it is generally assumed that the systematic errors following from nonergodicity are insignificant when compared with statistical ones [^{28}–^{31}], particularly for relatively short (up to 100 monomers) selfavoiding walks (SAWs) [^{32}]. However, because the model is intended for use in systems in which the nonergodicity error may be more pronounced due to the stretching of polymer chain and/or the presence of space constraints, a hybrid algorithm combining the local and nonlocal moves was employed in simulations. The use of nonlocal algorithm allows elimination of the nonergodicity problem [^{27}, ^{32}, ^{33}]. In this study, the “cutandpaste” method [^{27}] was chosen for nonlocal modifications of the chain. The use of the bond fluctuation algorithms [^{28}] was abandoned because their essential drawback revealed against the model demands. This drawback is ambiguous definition of the monomer neighborhood, which implies problems when the effect of the energy of monomer—monomer or monomer—solvent molecule interactions on the chain conformation has to be considered.
The aim of the study was to examine the statistics of local micromodifications in the simulation of chain conformation by means of the MMC method. The simulations were performed both for free undisturbed chains and for those whose two terminals were fixed over a selected distance. Three types of solvent conditions (i.e., athermal, theta and poor) were considered. Also, the chain in the presence of geometric constraints (the chain in a slit formed by two parallel impenetrable surfaces) was studied. As a probe of a local chain conformation and a local environment the kinkjump algorithm was used. The calculated values of the Helmholtz free energy of a SAW chain were confronted with the results available in literature [^{34}], obtained by the expanded ensemble Monte Carlo method.
The majority of results reported in this paper concern the polymer chain represented by SAW embedded to a simple cubic lattice of the lattice constant b. One monomer is identified with the help of a site on the lattice (bead), so the chain is an ensemble of N consecutive occupied sites. Mutual positions of beads are limited to the vector set being a permutation of (0, 0, ±b), which corresponds to the lattice coordination number ω equal to 6. The simulations were performed for two different bead numbers N, i.e., N = 50 and 100 (the chain of N beads is equivalent to M = N–1 SAW steps, each of the length b) in a box with the edge length considerably greater than the contour length of the longer chain (601b). Periodic boundaries of the simulation box were applied.
The simulations were completed using the Metropolis Monte Carlo (MMC) method [^{27}]. In this study we performed two kinds of simulations: 1) for a free unperturbed chain which was employed as a reference system and 2) for chains of arbitrarily chosen distance L_{H} between terminal beads of fixed coordinates. In the first case the initial chain conformation was generated by the static MC method [^{27}] and next it was randomly modified using algorithms of the following local micromodifications: endbead move (E), kinkjump (K), crankshaft (C) and reptation (R) motions [^{28}] and two nonlocal cutandpaste modifications: the inversion of a chain piece (I) and the reflection of a chain piece through a perpendicular bisector of a straight line joining the first and last bead of this subchain (F) followed by motion (I) (the modification F is necessary for ergodicity [^{27}, ^{33}]). In the second case the simulation started from a regular spirallike conformation which was next equilibrated via four types of modifications: K, C, I and F.
The MC step, defined as the number of modifications needed to give each of the beads the possibility to move once, contained N local micromodifications chosen at random, thus the frequencies of K and C motions were on average the same, whereas those of E or R motions (if applied) were N/2 times lower. The probability of performing the nonlocal cutand paste modification, involving simultaneous changes in coordinates of N/3 beads on average, was 10^{3} times smaller than that of local modifications. Thus, the ratio between frequencies of attempted local “internal” moves (i.e., K or C motions), micromoves of end beads (i.e., E or R motions) and nonlocal cutendpaste moves per a bead was 1:0.02:0.03, respectively. Because the increase in the frequency of cutandpaste modifications up to 10 times involved no changes in the results obtained (i.e., the differences between the results were within the statistical errors of the simulation), it was assumed that the frequency applied was sufficient to ensure the ergodicity of the procedure. This frequency of the nonlocal modifications allowed efficient simulations. On the other hand, the relatively high frequency of local modifications ensured the statistical significance of the probability of K motions, which were used as a probe of chain conformation.
In addition to the simulations performed for pure SAW chains where all molecular interaction energies are assumed as equal to 0, the conformations of chains immersed in the theta solvent (ε_{PS}=ε_{SS} = 0, ε_{PP}/k_{B}T = −0.2693 [^{35}]) and in poor solvent (ε_{PS}=ε_{SS} = 0, ε_{PP}/k_{B}T = −1) were also simulated. Symbols ε_{PS}, ε_{SS} and ε_{PP} denote the energies of bead – solvent molecule, solvent molecule – solvent molecule interactions and the energy assigned to each nonconsecutive pair of beads lying on neighboring lattice sites, respectively. In the course of these simulations, new conformations were accepted or rejected by the Metropolis procedure in the standard way [^{20}, ^{28}].
Some simulations were performed for the chain with terminal beads attached to two parallel surfaces forming a slit of the width equal to L_{H} + 2b. The surfaces were rigid and impenetrable to the polymer. They were purely repulsive and there was no chain adsorption (athermal surface) except for the irreversible attachment of the chain ends. Thus, the influence of surfaces on the polymer conformation was solely of entropic character.
All the results presented in the paper are averaged over 10^{4} conformations. Each conformation is a result of 10^{6} MC steps. The first half of motions (i.e., 0.5⋅10^{6} MC steps) were used for the initial chain thermalization. For the second half the following parameters characterizing the chain conformation were computed: the probability of performing a given move without violation of the continuity of the chain (i.e., the probability that the skeletal constraint is satisfied) and the probability of move execution without violation of the excluded volume constraint.
The probabilities discussed are described by the symbol P_{X}(Y)_{Z}, where X is the movement limitation (skeletal or excluded volume denoted as S and E, respectively), Y is the type of micromodification (E, K, C or R) and Z (if needed) indicates the state of the chain (Z = 0 for the unperturbed chain).
For comparative purposes, the simulations for the conventional random walk were also performed.
In this paragraph we present the results of calculations of probabilities of local micromodifications considered for chains with both free and fixed ends, obtained on the basis of simulation of chain conformations in athermal conditions. Table 1 summarizes the values of probabilities predicted theoretically and obtained from simulations in the case of free chain. Details of considerations on which the theoretical predictions were based are described below. The probabilities of moves satisfying the chain continuity constraint and of those satisfying the excluded volume constraint are considered separately. The results obtained for the chain with fixed distance between ends are only briefly presented here; their detailed discussion is given in further parts of the paper.
The obtained values of P_{S}(R)_{0} and P_{S}(E)_{0} equal to 1 are easy to predict since by no means the move of end bead of the chain can violate the skeleton constraint.
The prediction of P_{S}(K)_{0} and P_{S}(C)_{0} is somewhat more complex. The probability P_{S}(K)_{0} can be evaluated as follows. First, note that only corner located beads can perform kinkjumps, as schematically elucidated in Fig. 1. The probability that three succeeding beads (indexed as i–1, i, i + 1) form the corner is given by the ratio (ω–2)/(ω–1), as can be easily deduced from Fig. 2a. The number of possible corners with vertices on bead i can be additionally reduced by the occupancy of the vertex adjacent sites by other beads. The probability of occupancy of a single such site with a given bead decreases rapidly with its distance from bead i (measured along the chain) and thus, for instance (see Fig. 2b), for beads with indices i–3, i–5…, these probabilities are: 1/(ω–1)^{2}, 5/(ω–1)^{4} and the like. Hence, the value of P_{S}(K)_{0} can be found from the following equation:
1 [
The calculation of the probability P_{S}(C)_{0}, which is equivalent to the probability of occurrence of a local Ushaped conformation, is based on similar considerations as those presented above. At first, let us notice that the formation of such a conformation can be thought as performed in two steps (see Fig. 3a): 1) the formation of a corner by beads i–1, i and i + 1 (the same as in the case of kinkjump movement), 2) the addition of a third edge by suitable location of (i + 2)th bead. Then, the probability of formation of this microconformation is equal to P_{S}(K)_{0}⋅1/(ω–1). However, the second step can be hindered by the occupancy of the “suitable” lattice site by other beads. As it can be deduced from Fig. 3b, the probabilities that this site is blocked by beads with indices i–2, i–4… are equal to 1/(ω–1), 3/(ω–1)^{2}… Thus finally, the probability that the randomly chosen part of the chain can perform the crankshaft move satisfying the skeleton constraint reads:
2 [
The P_{S}(C)_{0} calculated from Eq. 2 is equal to 0.103 and roughly agrees with that obtained from the simulation (see Table 1).
The probabilities of motions not violating the chain continuity constraint, calculated directly from the simulation results for the chain with fixed positions of terminal beads, are presented in Fig. 4. As expected, both probabilities tend to zero when L_{H} reaches the contour length of the chain.
The probabilities of micromodifications not violating the excluded volume constraint are also not difficult to predict in the case of endsegment and reptation moves. The former move, if permitted, changes the orientation of a selected end bead by 90° or 180°, whereas the latter move cutsoff (with the probability equal to 1) one bead from either end of the chain and reattaches it to the opposite chain end. Thus, one can consider these probabilities to be measures of the number of vacant lattice sites adjacent to the selected end bead. One can also notice, that the product P_{E}(R)_{0}(ω–1) is equivalent to the effective coordination number of the lattice
3 [
Prediction of the P_{E}(K)_{0} probability is more complex. As results from Fig. 5, this probability depends on the bead location in the chain and it assumes the highest values (close to P_{E}(R)_{0} and P_{E}(E)_{0}) at the chain ends, rapidly decreasing to a certain almost constant value when the bead distance from the ends (measured along the chain) increases. Such a character of this dependence is consistent with intuitive expectation, and thus allows a conclusion that P_{E}(K)_{0} provides information about the likelihood of vacancy of a site adjacent to ith single bead of the chain.
In order to find a correlation between the probability P_{E}(K)_{0} and the number of vacant sites around the chain let us analyze the occupancy of the site adjacent to a single ith bead. The probability that the site is not empty can be considered as the probability of the union of two independent events: the occupation of the site by a bead with index lower than i and the site occupancy by a bead with index larger than i. Assuming that these two events can occur with the same probability X, the site occupation probability is then equal to 2X–X^{2}. Now, let us consider the probability Y of occupancy of the goal site of a kinkjump move. Since this site has two neighbors (as shown in Fig. 6a), the probability of its occupancy is again equivalent to the probability of the union of two independent events. Thus, on the basis of the above reasoning, we can write:
4 [
Now in turn, let us notice that the number n_{2} of possible chain trajectories passing through a site neighboring two beads (where the kinkjump move is permitted) is lower than the respective number n_{1} for a site neighboring only one bead (where the kinkjump is disallowed). This means that the probability Z of occupancy of a site adjacent to the chain disallowed for the kinkjump move is higher than the probability Y by a factor w:
5 [
Finally, combining Eqs. 4 and 5 and assuming that X = 1–P_{E}(R)_{0}, we get the following relation:
6 [
From Eq. 6 we obtain P_{E}(K)_{0} equal to 0.734 if the value of P_{E}(R)_{0} calculated from simulation results is used or P_{E}(K)_{0} = 0.714 if the value of P_{E}(R)_{0} determined on the basis of Eq. 3 is employed. Taking into account that these results are only rough estimates their agreement with the value of P_{E}(K)_{0} = 0.695 obtained from simulation (Table 1) seems to be satisfactory.
Having evaluated the probability of vacancy of a lattice site neighboring one bead (assumed to be equal to 0.714), one can find the probability P_{E}(C)_{0} that two such adjacent sites are empty (for the crankshaft move to be permitted) as equal to (P_{E}(K)_{0})^{2}, obtaining the value 0.510, which is close to that determined from the simulation.
The results obtained for the chain with fixed ends are presented in Fig. 7. The values of both P_{E}(K) and P_{E}(C) increase with increasing endtoend distance reaching the value of 1 for the fully stretched chain (i.e., L_{H}→bN ⇒ P_{E}(K)→1).
In this section we present a simplified approach to analytical evaluation of the conformational entropy of an SAW chain in a state of conformation distorted by fixing the ends in a selected distance. Because the ideal chain model based on the RW chain can be solved exactly, our approach will thus be to consider the random walk first and then to discuss the impact of the excluded volume effect.
The number of elementary steps M needed to generate the ideal chain is simply the sum of elemental steps performed in different directions:
7 [
8 [
There is no preferential direction in equilibrium. Now, let us assume that the chain is extended into + z direction. This elongation will be compensated by contraction in other directions. Assuming the same extent of contraction in all directions the number Ω of all possible conformations reads:
9 [
Equation 9 can be rewritten in a more general form:
10 [
The number of steps should fulfill the following requirements:
11 [
12 [
Replacing the number of steps, M, by the number of beads N (M = N–1) and generalizing the above reasoning over different lattice types, one can write the following expression for the conformational entropy of a stretched chain
13 [
The entropy of an undeformed chain is maximum (S_{max}) for L = 0 when the chain assumes the relaxed coillike conformation. Applying the Stirling approximation
14 [
15 [
16 [
Table 2 presents the estimates of S_{max} obtained on the basis of Eq. 15 for several models of polymer chains: RW, NRRW and SAW. For the sake of comparison there are also the corresponding values S_{max}^{*} calculated using other approaches (for details see the footnote to Table 2). As seen, in all cases the results show a good agreement (the greater N the better the coincidence between the results). This agreement allows one to assume that Eq. 13 can be also applied for the approximate estimation of conformational entropy of uniaxially deformed (elongated or shrunk) chains, at least those modeled by RW, NRRW and SAWs.
Equation 13 relates the chain entropy to the average coordinate L, or, in other words, the magnitude of the average value of the projection of the endtoend distance L_{H} onto a particular axis (here, onto the zaxis). In the case of ideal chain the L_{H} distance should fulfill two requirements: i) for an unperturbed chain, for which L = 0, the distance L_{H} should be equivalent to the rootmeansquared (rms) endtoend distance 〈R_{H}^{2}〉^{1/2} which is equal to (bM^{1/2}) [^{42}]
20 [
21 [
Assuming that there is a linear dependence between L and L_{H} and generalizing the above considerations to any chain model one obtains:
22 [
23 [
The application of Eq. 22 for determination of the dependence between L and L_{H} requires the knowledge of Flory’s exponent for a given model or polymer/solvent system. For the real chains (i.e., with the excluded volume interactions) this exponent is equal to 3/5 for good solvents, 1/2 for theta solvents and 1/3 for poor solvents [^{27}]. Figure 8 presents the dependencies of the conformational entropy of a chain on the distance L_{H} between its ends, calculated on the basis of Eqs. 13 and 22 for some chosen chain models. As seen, the curves go through maxima at L_{H}=〈R_{H}^{2}〉^{ν}, which is in accordance with expectation.
In the calculation of S = f(L_{H}) for the SAW model, which is presented in Fig. 8, a constant value of ω_{eff} = 4.6839 (valid for N→∞) was assumed. However, for the chains with excluded volume interactions the value of ω_{eff} depends on L_{H}, because the probability of nearestneighbor site occupancy changes with the extent of chain elongation. And thus, one can easily deduce that when the chain is highly stretched the excluded volume effect becomes negligibly small and ω_{eff} tends to 5 if the simple cubic lattice is considered or to ω–1 in a general case.
Unfortunately, except for the limiting case of fully extended chain, the determination of ω_{eff} is not straightforward since there is no analytical way of finding this parameter. In the next section of the paper we propose a method for the estimation of ω_{eff} on the basis of selected observables which are directly accessible from MC simulations. Moreover, to extend the applicability of Eqs. 13 and 22 to the case of chains in constrained environments we propose the employment of such observables in the quantification of the degree of chain deformation instead of L_{H} parameter. The meaning of L_{H} in the case of confined chains is sometimes unclear and, what is more important, not corresponding to the chain deformation degree. For instance, when a linear macromolecule is confined in a L or Ushaped nanochannel, the straight line connecting ends of the molecule may cross the channel walls and its length will not correspond to the extent of the molecule elongation.
As the probability P_{S}(K) is related to a selected separation between chain ends (Fig. 4) one can try to apply the value of this probability for the estimation of L or L_{H} parameters. For all the models studied here (2D RW, 3D RW, 2D SAW and 3D SAW) the dependencies of L versus the ratio P_{S}(K)/P_{S}(K)_{0} are of similar character: at L = 0 (free chain) the value of P_{S}(K)/P_{S}(K)_{0} is approximately 1 and then with increasing L it gradually decreases, approaching zero at L = bM (as seen in Fig. 9a). The shape of the curves suggests an exponential relationship between these two variables. Actually, as can be seen in Fig. 9b, plots of (1–P_{S}(K)/P_{S}(K)_{0}) against L on loglog coordinates give nearly straight lines with slopes dependent on the model considered.
From the detailed analysis of the results collected in Fig. 9 the following relationship was derived:
24 [
25 [
In Fig. 10a the values of L calculated on the basis of probabilities P_{S}(K) obtained from simulations are compared with the corresponding values assumed for different chain models. This comparison indicates a good agreement between the retrieved and assumed values and thus allows a conclusion that Eq. 24 can be used for the calculation of parameter L. Although in the range of small P_{S}(K) values, in which the equation becomes numerically unstable, there is a significant error in the L estimation. The described method of the determination of L is independent of the chain length, as can be concluded from inspection of Fig. 10b.
Equation 3 allows the calculation of the effective coordination number ω_{eff} on the basis of the probability P_{E}(R). However, the applicability of this method of ω_{eff} determination is restricted to a free chain. Searching for an alternative method and having in mind the approximate relation 6 one can try to use P_{E}(K) instead of P_{E}(R) as a local probe of ω_{eff}, and thus omit the requirement of free ends, necessary for P_{E}(R) calculation. Assuming that the validity of Eq. 6 can be extended to deformed chains with fixed ends and by combining it with Eq. 3 one obtains:
26 [
Equation 26 relates ω_{eff} with P_{E}(K). But, as could have been foreseen, its applicability is limited. For instance, it fails in predicting the results for boundary conditions. Namely, the solution of Eq. 26 does not satisfy the relations which can be written for two limiting states of the chain: fully stretched and fully shrunk. The corresponding relations are as follows:
27 [
28 [
In the search for a relationship between ω_{eff} and P_{E}(K), which would be of a more general applicability than Eq. 26, we made use of the results of simulations performed for the models including intrachain interactions, namely, for polymer chains in theta and poor solvents. The results obtained for such models allow direct and independent determination of both ω_{eff} and P_{E}(K). The attractive intrachain interactions, which exist in theta (χ=1/2) and poor (χ>1/2) solvents, modify the value of ω_{eff} as compared to that in athermal environment. Moreover, these interactions change the internal energy U of the model to be nonzero. The average value of U is given by:
29 [
Using the values of U determined from the simulations in which ε_{PS}=ε_{SS} = 0 was assumed, the average value of ω_{eff} can be calculated form the following equation:
30 [
The values of ω_{eff} obtained in this way are summarized in Fig. 11, which shows the relation between ω_{eff} and P_{E}(K) in the form of a function ((ω–1)–ω_{eff})^{1/2} = f(1–P_{E}(K)). As seen, all the data points lie quite well along a straight line whose slope is equal to 2.0, indicating that the searched relationship between ω_{eff} and P_{E}(K) is well approximated by the equation
31 [
It is noteworthy that Eq. 31 fulfills Eqs. 27 and 28.
It follows from Fig. 12a that, as expected, when L_{H} tends to the contour length of the chain, the value of ω_{eff} calculated on the basis of Eq. 30 tends to ω–1. The most distinct decrease in ω_{eff} accompanying the decrease in L_{H} is seen in the case of a chain in poor solvent conditions, which is also as expected. The results collected in Fig. 12b indicate that the chain length has only a minor influence on the obtained ω_{eff} values, at least in the studied N range (N = 50 and N = 100).
In this section, the results of calculations of the conformational entropy and energy of a polymer, obtained by means of the method based on the use of probabilities P_{S}(K) and P_{E}(K), are presented. The calculations were performed for different selected endtoend distances and solvent conditions, as well as in the presence of geometrical constraints.
Figure 13a presents the entropies S of a SAW chain of N = 50 obtained at different L_{H} values. These results are compared with the corresponding results available in literature and obtained with the expanded ensemble Monte Carlo method [^{34}]. Practically, in the whole range of L_{H} values a very good agreement is observed. It is also noteworthy that as expected, the maximum on the curve S(L_{H}) occurs at the L_{H} value corresponding to the rms distance between ends of a free chain: R_{H}/bM ≈0.20 for SAW of N = 50 and v = 3/5 (Eq. 23). For the sake of comparison with the results reported in ref. [^{34}], the stretching forces F required to extend the chain to assumed endtoend separations were also calculated. The corresponding values of F were found from the following relation:
32 [
As illustrated in Fig. 13b, the obtained results are in fairly good agreement with those of ref. [^{34}].
Using the method described in this paper, the values of Helmholtz free energies A of the chain of N = 100 in different solvent conditions (i.e., of different values of Flory’s interaction parameter: χ=0, χ=1/2 and 1/2<χ≈2) were calculated from:
33 [
The obtained results, together with corresponding results for U and S, are shown in Fig. 14. As seen, for large extensions (large L_{H}) the conformational entropy of the chain is practically independent of the solvent condition (Fig. 14a). It can be explained by the fact that in the case of strongly stretched chain the probability of beadbead contacts (which mimic intermonomer interactions) is very small and thus its influence on the chain entropy is also small. For intermediate and low extensions the entropy of the chain in solvents of χ≥1/2 is smaller than in the athermal condition, because of smaller values of ω_{eff} (compare Fig. 12a). As seen in Fig. 14b, in nonathermal conditions the decrease in L_{H} is accompanied by a monotonic decrease in internal energy of the chain. The superposition of energetic and entropic contributions produces the A(L_{H}) dependencies exhibiting minima (Fig. 14c). In spite of the low precision of estimation of positions of these minima (which is because of relatively small length of the examined chain, N = 100), it can be found out that they occur at L_{H} values roughly equal to rms endtoend distances of the free chain in corresponding solvent conditions (the values of R_{H}/bM calculated from Eq. 23 are 0.16, 0.10 and 0.05 for athermal, theta and poor solvent conditions respectively).
The described method of estimation of the conformational entropy of macromolecules was also tested for a chain in the space with geometrical constraints. The entropies of SAW chain with its ends attached to two opposite parallel plane surfaces which were also assumed to be impenetrable, infinite and noninteracting with the chain except its terminals were calculated at different separations between the surfaces. The results obtained are collected in Fig. 15 and compared with the entropies found for the chain in the open space when selected distance L_{H} between ends was equal to that in the presence of the surfaces. For the confined chain in the range of small plane separations (L_{H} + 2b ≤ 〈R_{H}^{2}〉^{v}) the calculated entropies are clearly smaller as compared to those of the chain in the open space. This can be elucidated by the excluded volume effect of the surfaces, which reduces the value of ω_{eff}. From Fig. 15, it can also be deduced that for a chain with free ends the average endtoend distance in the equilibrium conformation is longer in the presence of the surfaces.
The probabilities of different types of micromodifications of chain conformation which were found from the performed MMC simulations agree reasonably well with the results obtained on the basis of approximate theoretical considerations regarding the nearest environment of a local chain conformation. The proposed method for the estimation of conformational entropy of a polymer chain using the probabilities of kinkjump moves, P_{S}(K) and P_{E}(K), provides results consistent with literature [^{34}] when applied to a deformed chain (i.e., whose two ends are fixed over a distance L_{H} different from the most probable endtoend separation of an undisturbed chain). In the range of small deformations the entropy estimation accuracy was low because of the numerical instability of Eq. 31. The dependencies of the conformational entropy S and free energy A of the chain on the imposed endtoend distance L_{H}, calculated on the basis of probabilities P_{S}(K) and P_{E}(K) for different solvent conditions as well as in the presence of the geometric constraints in the simulation space, qualitatively capture major features of the expected behavior of a polymer chain under such “test” conditions.
The work was partly supported from the resources for science in 2010–2013 by the Ministry of Science and Higher Education, Poland, grant number: N N204 015238.
References
1..  Fleer G,Cohen Stuart M,Leermakers F. Lyklema JEffect of polymers on the interaction between colloidal particlesFundamentals of interface and colloid science, Vol. V: soft colloids.Year: 2005New YorkElsevier 
2..  Milner ST,Witten TA,Cates ME. Theory of the grafted polymer brushMacromoleculesYear: 1988212610261910.1021/ma00186a051 
3..  Lodge TP,Muthukumar M. Physical chemistry of polymers: entropy, interactions, and dynamicsJ Phys ChemYear: 1996100132751329210.1021/jp960244z 
4..  Kutner I,Srebnik S. Conformational behavior of semiflexible polymers confined to a cylindrical surfaceChem Phys LettYear: 2006430848810.1016/j.cplett.2006.08.085 
5..  PoncinEpaillard F,Vrlinic T,Debarnot D,Mozetic M,Coudreuse A,Legeay G,El Moualij B,Zorzi W. Surface treatment of polymeric materials controlling theadhesion of biomoleculesJ Funct BiomaterYear: 2012352854310.3390/jfb3030528 
6..  Nowicki W,Nowicka G,Dokowicz M,Mańka A. Conformational entropy of a polymer chain grafted to rough surfacesJ Mol ModelYear: 20131933734822918701 
7..  Arya G. Energetic and entropic forces governing the attraction between polyelectrolytegrafted colloidsJ Phys Chem BYear: 2009113157601577010.1021/jp908007z19842639 
8..  Charlaganov M,Košovan P,Leermakers FMA. New ends to the tale of tails: adsorption of comb polymers and the effect on colloidal stabilitySoft MatterYear: 200951448145910.1039/b816832f 
9..  Tadros T. Interaction forces between particles containing grafted or adsorbed polymer layersAdv Colloid and Interface SciYear: 200310419122610.1016/S00018686(03)00042312818496 
10..  Petrey D,Honig B. Free energy determinants of tertiary structure and the evaluation of protein modelsProtein SciYear: 200092181219110.1110/ps.9.11.218111152128 
11..  Liu SQ,Ji XL,Tao Y,Tan DY,Zhang KQ,Fu YX. Kaumaya PProtein folding, binding and energy landscape: a synthesisProtein engineeringYear: 2012RijekaInTech 
12..  Bereau T,Bachmann M,Deserno M. Interplay between secondary and tertiary structure formation in protein folding cooperativityJ Am Chem SocYear: 2010132131291313110.1021/ja105206w20822175 
13..  Tzlil S,Kindt JT,Gelbart WM,BenShaul A. Forces and pressures in DNA packaging and release from viral capsidsBiophys JYear: 20038431616162710.1016/S00063495(03)74971612609865 
14..  Zappa E, Indelicato G, Albano A, Cermelli P (2013) A GinzburgLandau model for the expansion of dodecahedral viral capsid. Int J NonLinear Mech. doi:10.1016/j.ijnonlinmec.2013.03.003 
15..  Šiber A,Zandi R,Podgornik R. Thermodynamics of nanospheres encapsulated in virus capsidsPhys Rev EYear: 20108105191905193010.1103/PhysRevE.81.051919 
16..  Petrov AS,Boz MB,Harvey SC. The conformation of doublestranded dna inside bacteriophages depends on capsid size and shapeJ Struct BiolYear: 2007160224124810.1016/j.jsb.2007.08.01217919923 
17..  Binder K,Paul W. Recent developments in Monte Carlo simulations of lattice models for polymer systemsMacromoleculesYear: 20082008414537455010.1021/ma702843z 
18..  Meirovitch H. Recent developments in methodologies for calculating the entropy and free energy of biological systems by computer simulationCurr Opin Struct BiolYear: 20071718118610.1016/j.sbi.2007.03.01617395451 
19..  Muthukumar M (2011) Polymer translocation. CRC, New York 
20..  Metropolis N,Ulam S. The Monte Carlo methodJ Am Stat AssocYear: 19494433534110.1080/01621459.1949.1048331018139350 
21..  Landau DP,Binder K. A guide to Monte Carlo simulations in statistical physicsYear: 2000New YorkCambridge University Press 
22..  Verdier PH,Stockmayer WH. Monte Carlo calculations on the dynamics of polymers in dilute solutionJ Chem PhysYear: 19623622723510.1063/1.1732301 
23..  Verdier PH. A simulation model for the study of the motion of randomcoil polymer chainsJ Comput PhysYear: 1969420421010.1016/00219991(69)900679 
24..  Binder K,Paul W. Monte Carlo simulations of polymer dynamics: recent advancesJ Polym Sci Part B: Polym PhysYear: 19973513110.1002/(SICI)10990488(19970115)35:1<1::AIDPOLB1>3.0.CO;2# 
25..  Kloczkowski A,Koliński A. Mark JETheoretical models and simulations of polymer chainPhysical properties of polymers handbookYear: 20072New YorkSpringer 
26..  Madras M,Sokal AD. Nonergodicity of local, lengthconserving Monte Carlo algorithms for the selfavoiding walkJ Stat PhysYear: 19874757359510.1007/BF01007527 
27..  Sokal AD. Binder KMonte Carlo methods for the selfavoiding walksMonte Carlo and molecular dynamics simulations in polymer scienceYear: 1995New YorkOxford University Press 
28..  Carmesin I,Kremer K. The bond fluctuation method: a new effective algorithm for the dynamics of polymers in all spatial dimensionsMacromoleculesYear: 1988212819282310.1021/ma00187a030 
29..  SalesPardo M,Guimerà R,Moreira AA,Widom J,Amaral LA. Mesoscopic modeling for nucleic acid chain dynamicsPhys Rev EYear: 20057105190205191510.1103/PhysRevE.71.051902 
30..  Binder K. Applications of Monte Carlo methods to statistical physicsRep Prog PhysYear: 19976048755910.1088/00344885/60/5/001 
31..  Quake SE. Topological effects of knots in polymersPhys Rev LettYear: 1994733317332010.1103/PhysRevLett.73.331710057346 
32..  Baschnagel J, Wittmer JP, Meyer H (2004) Monte Carlo simulation of polymers: coarsegrained models. In: Attig N, Binder K, Grubmüller H, Kremer K (eds) Computational soft matter: from synthetic polymers to proteins, lecture notes, John von Neumann institute for computing, NIC Series, Vol. 23, Jülich 
33..  van Rensburg EJJ. Monte Carlo methods for the selfavoiding walkJ Phys A: Math TheorYear: 20094232300132309810.1088/17518113/42/32/323001 
34..  VorontsovVelyaminov PN,Ivanov DA,Ivanov SD,Broukhno AV. Expaded ensemble Monte Carlo calculations of free energy for closed, stretched and confined lattice polymersColloid Surface Physicochem Eng AspectYear: 199914817117710.1016/S09277757(98)005482 
35..  Panagiotopoulos AZ,Wong V,Floriano MA. Phase Equilibria of lattice polymers from histogram reweighting Monte Carlo simulationsMacromoleculesYear: 19983191291810.1021/ma971108a 
36..  Zhao D,Huang Y,He Z,Qian R. Monte Carlo simulation of the conformational entropy of polymer chainsJ Chem PhysYear: 19961041672167410.1063/1.470753 
37..  Nowicki W. Sructure and entropy of a long polymer chain in the presence of nanoparticlesMacromoleculesYear: 2002351424143610.1021/ma010058d 
38..  Watts MG. Application of the method of Pade approximants to the excluded volume problemJ. Phys. AYear: 19758616610.1088/03054470/8/1/012 
39..  Sykes MF,Guttman J,Watts MG,Robberts PD. The asymptotic behaviour of selfavoiding walks and returns on a latticeJ Phys AYear: 1972565366010.1088/03054470/5/5/006 
40..  Barber MN,Guttman AJ,Middlemiss KM,Torie GM,Whittington SG. Some tests of scaling theory for a selfavoiding walkJ Phys A Math GenYear: 1978111833184210.1088/03054470/11/9/017 
41..  Őttinger HC. Computer simulation of threedimensional multiplechain systems: scaling laws and virial coefficientsMacromoleculesYear: 198518939810.1021/ma00143a015 
42..  Teraoka I (2002) Polymer solutions. an introduction to physical properties. Wiley, New York 
Figures
[Figure ID: Fig1] 
Fig. 1
The kinkjump motion trials for two different locations of a bead in the chain a the forbidden move, b the permitted move 
[Figure ID: Fig2] 
Fig. 2
Illustration of estimation of P_{S}(K)_{0}: a possible conformations of a chain fragment composed of three beads with indices i–1, i and i + 1, b a possible location of (i–3)th bead which reduces the number of corner forming conformations (in red – the chain fragment whose formation is disallowed) 
[Figure ID: Fig3] 
Fig. 3
Illustration of the estimation of P_{S}(C)_{0}: a a possible conformation of a chain fragment composed of four beads with indices i–1, i, i + 1 and i + 2, b an exemplary possible (i–4)th bead location disallowing the completion of Ushaped conformation by suitable orientation of (i + 2)th bead (in red – the chain fragment whose formation is disallowed) 
[Figure ID: Fig4] 
Fig. 4
The influence of L_{H} value on probabilities P_{S}(K) and P_{S}(C). The horizontal lines indicate corresponding probabilities determined for the free chain (SAW, N = 100) 
[Figure ID: Fig5] 
Fig. 5
The dependence of the relative probability that a bead adjacent site is empty on the position of the bead along the chain (SAW model, N = 100, results obtained from a single simulation) 
[Figure ID: Fig6] 
Fig. 6
Visualization of the neighborhood of two local microconformations containing a bead (indicated by the arrow) which a can and b cannot be moved by the kinkjump. Note that for the move to be permitted the target site has to neighbor two other beads, while the other “disallowed” sites neighbor one bead only (see blue lines in the figure). The probability of the site occupancy by other, more distant beads is different in both cases. In red there are marked the possible trajectories of chain fragment composed of three beads passing through the site a allowed and b disallowed for the kinkjump move, whose numbers are equal to 20 and 25, respectively 
[Figure ID: Fig7] 
Fig. 7
The dependencies of P_{E}(K) and P_{E}(C) on L_{H}. The horizontal lines indicate the corresponding probabilities determined for the free chain—P_{E}(K)_{0} and P_{E}(C)_{0}, respectively (SAW model, N = 100) 
[Figure ID: Fig8] 
Fig. 8
Dependencies of conformational entropies vs. fixed endtoend distance calculated for 2D RW, 3D RW and 3D SAW (N = 100). Perpendicular lines indicate rms endtoend distances of RW and SAW chains equal to 10b and 15.8b, respectively. The rms endtoend distance for RW model is independent on the dimensionality of the lattice 
[Figure ID: Fig9] 
Fig. 9
Plots of: aP_{S}(K)/P_{S}(K)_{0} against L (inset: enlarged part of the dependence around L = 0, coordinates are not marked) and b (1–P_{S}(K)/P_{S}(K)_{0}) against L in loglog scale for positive values of variables; N = 100 
[Figure ID: Fig10] 
Fig. 10
The plots of the values of L calculated on the basis of probabilities P_{S}(K) against the assumed L values: a for different chain models, b for two different chain lengths (3D SAW). The slope of the straight line is equal to 1 
[Figure ID: Fig11] 
Fig. 11
The values of the expression ((ω–1)–ω_{eff})^{1/2} plotted against 1–P_{E}(K), obtained for a chain of N = 100 in theta and poor solvents. The regression coefficient of the straight line indicated in the figure is equal to 2 
[Figure ID: Fig12] 
Fig. 12
The plots of the values of ω_{eff} calculated from Eq. 30 against the assumed values of L_{H} for a different solvent conditions, N = 100 and b two different lengths of SAWs 
[Figure ID: Fig13] 
Fig. 13
The relationships of aS vs L_{H} and bF vs L_{H}; SAW chain of N = 50. For comparison, the same dependencies calculated by means of expanded ensemble MC method [^{34}] (squares) are shown 
[Figure ID: Fig14] 
Fig. 14
The comparison of relationships of aS vs L_{H}, bU vs L_{H} and c) A vs L_{H}, which were determined for three different solvent conditions: χ=0 (athermal), χ=1/2 (theta) and 1/2<χ≈2 (poor), N = 50. The L_{H} values corresponding to minima in the free energy, predicted by Eq. 23, are marked by vertical lines 
[Figure ID: Fig15] 
Fig. 15
The entropy of SAW chain (N = 100) with ends attached to the opposite parallel surfaces versus the surface separation (circles). For comparison, the corresponding results for the chain with fixed ends but in the open space are indicated (squares) 
Tables
The probabilities of different types of motions estimated theoretically (theor) and calculated on basis of the simulation (sim) performed for unperturbed chains of N = 100 and N = 50 in the athermal conditions
Quantity  Micromodification  

(R)_{0,theor}  (R)_{0,sim}  (E)_{0,theor}  (E)_{0,sim}  (K)_{0,theor}  (K)_{0,sim}  (C)_{0,theor}  (C)_{0,sim}  
P_{S}N = 100  1.000  1.000  1.000  1.000  0.763  0.774  0.103  0.101 
N = 50  1.000  1.000  0.776  0.100  
P_{E}N = 100  0.937  0.942  0.937  0.938  0.714  0.695  0.510  0.506 
N = 50  0.943  0.940  0.696  0.511 
Comparison of polymer chain entropies S_{max} calculated from Eq. 13 for various chain models and two different chain lengths with the corresponding values [
System  ω  ω_{eff}  N = 1000  N = 10,000  

S_{max}/k_{B} 
[ 
100σ  S_{max}/k_{B} 
[ 
100σ  
2D RW  4  4  1377  1385  0.58  13850  13863  0.09 
3D RW  6  6  1776  1790  0.78  17896  17918  0.12 
2D NRRW  4  3  1092  1098  0.54  10977  10986  0.08 
3D NRRW  6  5  1597  1608  0.68  16077  16094  0.11 
3D SAW  6  4.6839 [^{38}]  1533  1544  0.71  15426  15441  0.10 
*) In calculation of [
a) for the RW model: [
b) for the NRRW model: [
c) for the SAW model: [
where C’ = 1.17, γ=7/6 and ω_{eff,∞}=4.6838 for the 3D linear, infinitely long macromolecule in the athermal solution [^{38}–^{41}]
Article Categories:
Keywords: Keywords Lattice polymers, Self avoiding walk. 
Previous Document: A Novel Technique for Sculpting Costal Cartilage in Microtia Repair and Rhinoplasty.
Next Document: A coupled twodimensional main chain torsional potential for protein dynamics: generation and implem...