Property Prediction



5       RMMC (RIS Metropolis
Monte Carlo)

 The RMMC module is used to compute the conformational properties of a polymer chain. RMMC performs an RIS Metropolis Monte Carlo (RMMC) simulation of your polymer, calculates the properties that you specify, and can then save the conformations for later viewing or analysis. You can also choose to modify various parameters controlling the simulation to enhance the quality of the results.

Note

Additional documentation for RMMC is available at: HTTP://www.msi.com/doc. Please see Appendix B, "Using MSI Online Documentation" for more information about accessing these documents.  


Using RMMC

This section contains a simple example to get you started with the RMMC module. Steps you must perform for the tutorial to work as designed are in boxes. Explanations are printed in italics. Card deck, card, item, and control panel widget names are given in bold.

The specific example described in this section allows you to estimate the chain statistics of a linear alkane melt at room temperature.

1.   Starting Cerius2

Open a new UNIX shell and type: 

>	cerius2

2.   Accessing POLYMER BUILDER

Within the Visualizer, go to the POLYMER1 card deck and choose the POLYMER BUILDER card.

3.   Build a decane molecule

Select the Homopolymer item on the POLYMER BUILDER card. When the Homopolymer Builder panel appears, click the BUILD button.

You now see a decane molecule in the model window.

4.   Setting the run mode to INTERACTIVE

Go to the RMMC card and select the Job Control item.

This displays the RMMC Job Control panel.

Click the yellow popup labeled Run Mode and change the setting from BACKGROUND to INTERACTIVE.

5.   Accessing the RMMC Run panel

Now select the Run item from the RMMC card.

This displays the RMMC Run panel. The Equilibration and
Simulation check boxes are selected, but the Preparation and
Minimization check boxes are not.

Check the Minimization check box.

In the text entry box labeled File Prefix, type the text: test_decane. Set the number of Equilibration Steps to 1000 and the number of Simulation Steps to 5000. Click the RUN button to launch the RMMC calculation.

You now see messages in the text window indicating that the simulation has started. The logfile for the calculation, test_decane.rislog, is automatically brought up in a separate window and is updated as the simulation proceeds. The calculation should taken only a few minutes to complete. On completion, the displayed logfile contains estimates for several properties calculated in the RMMC simulation, including the average end-to-end distance of the decane molecule.

6.   Display the results of the calculation

Go back to the RMMC card and select the Analyze item.

This brings up a file browser containing the names of files available for analysis from the RMMC simulation run.

Select the file test_decane_dih.tbl. Click the radio button labeled Plot DIHEDRAL Distribution.

Pushing the buttons displays the distribution of dihedral angles averaged over the simulation run. Click on the yellow popups labeled Distribution and Trajectory to select other simulation predictions, and click the appropriate Plot radio button to display the results for the selected property.

Note that this example is an unusually short simulation run. To obtain reliable average chain properties for a typical polymer system, you may need to conduct runs of a million steps or more to sample a large-enough subset of the configuration space of the polymer.


General Methodology

RIS Metropolis Monte Carlo (RMMC) Concepts

Rotatable Bonds

 In an RMMC simulation, you can potentially change the torsion angle of every rotatable bond. By default, all single and partial double bonds in the model are treated as rotatable. If the check box labeled Mark All Bonds Rotatable on the RMMC Monte Carlo Preferences panel is deselected (access this panel by selecting the Preferences->Monte Carlo item), the RMMC program uses backbone atom flags to determine which bonds are rotatable for simulation purposes. A bond is considered rotatable in this context if it satisfies all of the following criteria:

By using the commands on the Polymer Editor panel (accessed by going to the POLYMER BUILDER card and selecting the
Edit-->Polymers item), you can control which atoms are considered backbone atoms, and thus which bonds are considered rotatable. This allows you to treat branched chains, dendrimers, and polymers with flexible side groups using RMMC. In all cases, an unbroken backbone path through the molecule must exist for RMMC to work properly.

Energy Calculation

Unlike traditional RIS calculations, RMMC does not use statistical weights. Instead, it uses the energy as computed from a forcefield in order to calculate chain conformational properties. The only energy terms considered in a RMMC calculation are:

There are two parameters--Minimum Bonds and Maximum Bonds--that determine the interacting pairs of atoms for the purpose of calculating nonbond (van der Waals and Coulomb) energies. (Access these parameters on the RMMC Interactions panel, by going to the RMMC card deck and selecting the
Preferences->Interactions item.) Nonbond energies are not computed for atoms closer than Minimum BONDS bonds away from each other. (See Figure 4.) The usual value for Minimum BONDS is 3. Nonbond energies are also neglected for atoms further than Maximum BONDS bonds away from each other. Reasonable values for Maximum BONDS range from 4 to about 6 for polymer chains in theta conditions. (Larger values of Maximum BONDS may greatly increase CPU and memory requirements, although some excluded volume effects might, in principle, be treated in this way.

Figure 4 . Nonbond interactions in RMMC energy calculation

Black atoms interact with the central atom (indicated by the arrow); white atoms do not. In this example, Minimum BONDS is 3 and Maximum BONDS is 5.  

Parameters Controlling the Simulation

A number of parameters can affect the outcome of a RMMC simulation. These include the following:

Minimum BONDS and Maximum BONDS have already been mentioned. The temperature enters the Boltzmann factor that is used to determine whether to accept or reject a conformation. The forcefield provides the parameters for the torsional and nonbond energy calculation. The inclusion or exclusion of atomic partial charges and the dielectric constant determine the calculated Coulomb energy.

Whether articulated side groups are treated as flexible in a RMMC simulation is determined by whether the side group atoms have their backbone flags set. (The commands on the Monomer Editor and Polymer Editor panels can be used to set these flags. (Access these panels by going to the POLYMER BUILDER card and selecting the Edit-->Polymers and Edit--> Monomers items) It is best to set these flags on the repeat units before building the polymer.) In principle, it is more realistic to treat side groups as flexible rather than rigid. However, doing so increases the computation time and might not be necessary if the groups are small (as, for example, in polypropylene).

The "energy scaling factor" is a number that multiplies all computed energies. If experimental data are available, this factor can be used to refine the properties calculated by RMMC so that they match experimental properties as closely as possible. The same factor can then be used for calculations on related chains for which no experimental data exist.

The reason for the relatively large number of parameters controlling the energy calculation is as follows. A RMMC simulation is performed on an isolated polymer chain. Yet, the properties desired are those for a chain in solution or in the melt. If solvent molecules or other chain molecules were present in the simulation, then a high quality forcefield such as COMPASS should be capable of predicting the correct properties. Because these other molecules are not explicitly present, their effects must be mimicked by altering the way the energy is computed. (In the language of statistical mechanics, it is not "bare" interaction energies but potentials of mean force (PMF) that determine the polymer conformations in a solvent or in the presence of other chains. See McQuarrie (1976) for more on the PMF concept. The goal in a RMMC simulation is to have the computed energy come as close to the potential of mean force as possible.)

Using Maximum BONDS to limit the range of nonbond interactions is a way of mimicking theta conditions. However, there is no known a priori way to determine the ideal value of Maximum BONDS for a particular polymer. Thus, calculations with differing values of this parameter should be performed and the results evaluated according to their reasonableness or agreement with established data.

The dielectric constant provides another way of mimicking the presence of solvent. But it must be kept in mind that, at a molecular level, a solvent is not a dielectric continuum, so that this too is an approximation. Consequently, the best value for the dielectric constant in a RMMC calculation might not be the same as the system's bulk dielectric constant. In the event that variation of Maximum BONDS and the dielectric constant within reasonable limits does not give adequate results, the energy scaling parameter is available as a last resort.

"RIS" Metropolis Monte Carlo (RMMC) Simulations

If your chain is branched, or if its statistical weights are unknown, RMMC simulation is an alternative to traditional RIS approaches. To calculate the properties of a polymer chain using RMMC, perform the following steps.

1.   Build the chain using the POLYMER BUILDER card.

2.   Go to the RMMC card and select the Properties item; choose the properties you wish to calculate on the control panel that opens.

3.   Select various preferences for the simulation using the panels triggered by the Preferences item on the RMMC card.

4.   On the RMMC Run panel, define the length of each stage to conduct in the simulation and select a temperature. Choose a force field for use in computing the nonbond and coulombic interactions in the simulation.

5.   Use the Run button (the RUN_RMMC command) to execute the RMMC simulation.

6.   When the simulation finishes, use the Analyze item on the RMMC card to plot the computed chain properties.

7.   After looking over the results, perform the calculation with a different value of Maximum BONDS and/or the dielectric constant in order to test the sensitivity of the results to the simulation parameters.

Computing Dihedral Distribution Functions

RMMC calculates distribution functions for the dihedral angle that are averaged over all rotatable bonds in the model. There is currently no provision for the estimation of distribution functions for specific torsional types, or for the computation of bond pair distribution functions.

Output Files

An RMMC calculation may produce a number of output files. Listed below are the extensions and contents of these files.


Theory

RIS Metropolis Monte Carlo Simulation

 Conventional RIS methods require that statistical weights exist for the polymer of interest. The RIS Metropolis Monte Carlo (RMMC) simulation does not require statistical weights. Instead, RMMC uses a forcefield (such as CVFF, PCFF or COMPASS) to determine the conformational properties of a chain. Because it does not use statistical weights (or even discrete torsion angles), RMMC is not, strictly speaking, a rotational isomeric state method. However, it shares several assumptions with RIS theory, and is used to compute the same types of properties as RIS. In this sense, it is closely related to RIS theory.

Properties Calculated with RMMC

The RMMC program can directly calculate the following properties:

You can also save the chain conformations in a file for analysis by other programs.

The RMMC Algorithm

The RMMC method was developed at Biosym (now Molecular Simulations). Although most of the ideas underlying the method are not new, the particular way these ideas have been put together is. The closest method to RMMC that has been published was developed by Dodd and Theodorou (1994). As with all RIS-related methods, RMMC is appropriate for theta chains--whether in solution or in the melt. Long-range excluded volume interactions are neglected.

As in traditional RIS methods, only torsional degrees of freedom are considered in determining a chain's conformation; bond lengths and angles are fixed. Unlike these methods, RMMC allows torsion angles to vary continuously; it does not impose the assumption of discrete rotational states.

As its name indicates, RMMC is a Metropolis (Metropolis 1953) Monte Carlo method. (This can be contrasted to the Markovian approach used in a RIS Monte Carlo calculation to build independent chains.) In a Metropolis simulation of a polymer, one begins with a chain in an arbitrary conformation. A Monte Carlo step consists of making a small change to that conformation--e.g., by rotating a bond--and then deciding whether or not to retain that change, based on the temperature and the energy of the new conformation relative to the old one. This process is repeated many times in order to yield a set of conformations characteristic of that chain at the specified temperature.

In outline, a RMMC simulation proceeds as follows:

1.   Perform an energy minimization on the molecule (so that bond lengths and angles adopt reasonable values).

2.   Randomly select a rotatable backbone bond.

3.   Select a random torsion value for this bond between -180 and +180 degrees.

4.   Rotate the bond to its new torsion value and compute the new energy of the chain.

5.   Generate a random number, R, between 0 and 1. If exp[-(Enew -Eold)/kT] > R, keep the new torsion value. Otherwise, restore the old value.

6.   After enough steps have been done that the chain is equilibrated: compute the properties of the chain conformation (e.g., r2, s2), and update the running averages of these properties.

7.   Repeat from step (2) until the desired number of iterations has been performed.

Treatment of Constraints

 In a RMMC simulation, bond lengths and bond angles are constrained. (For this reason, "pre-minimization" is recommended in order that the bonds lengths and angles adopt reasonable values.) The next three paragraphs deal with some of the technical subtleties of simulations with constraints. Knowledge of these subtleties is not required for a basic understanding of the RMMC simulation approach, but is included here for the sake of completeness.

In a Monte Carlo simulation that incorporates constraints (such as fixed bond lengths and angles) and uses internal coordinates (rather than Cartesian coordinates) in making its moves, care must be taken that the conformations are sampled correctly. Simulations may be performed with constraints that are regarded as "rigid" or as "flexible." In the rigid case, the momenta conjugate to the constrained coordinates are neglected. In the flexible case, it is imagined that the force constants for the constrained coordinates are so large that the coordinates do not move significantly from their equilibrium values, but that the conjugate momenta are activated. It is hard to say a priori which is a more reasonable approximation for a polymer. In a real polymer, bond lengths and angles are, of course, not rigidly constrained. On the other hand, according to quantum mechanics, very stiff degrees of freedom are not activated at ordinary temperatures. Go and Scheraga (1976) argue that, in practice, flexible constraints should be more quantitatively accurate.

It should be emphasized that the distinction between "rigid" and "flexible" constraints is in the assumptions behind the physical model. Thus, even in the "flexible" case, the constrained coordinates do not move during the simulation.

A flexible constraint requires no special treatment in a Monte Carlo simulation. This is because the configuration integral, as expressed in terms of torsional degrees of freedom, contains a Jacobian-like term that is independent of these degrees of freedom. Thus, this term can be moved outside the integral, and is divided out of any conformational averages. In the rigid case, this is not true, and in principle the Boltzmann factor used in the Metropolis acceptance criterion needs to be multiplied by an extra term (Fixman 1974). In practice, this extra term has been found to change the conformational statistics only slightly (Almarza et al. 1990). For simplicity, the flexible constraint assumption, with no extra terms in the Boltzmann factor, is made for the RMMC algorithm. Any errors introduced by this approach are likely to be much smaller than those due to other approximations in the RMMC algorithm--such as approximating the potential of mean force by a forcefield with a truncated interaction range. For more on the issue of constraints in simulations, see Fixman (1974), Allen and Tildesley (1987), and Almarza et al. (1990).


References

Allen, M.P.; Tildesley, D.J. Computer Simulation of Liquids, Oxford University Press (1987).

Almarza, N.G.; Enciso, E.; Alonso, J.; Bermejo, F.J.; Alvarez, M. "Monte Carlo simulations of liquid n-butane," Mol. Phys., 70, 485 (1990).

Dodd, L.R.; Theodorou, D.N. "Atomistic Monte-Carlo Simulation and Continuum Mean-Field Theory of the Structure and Equation of State Properties of Alkane and Polymer Melts," Advances in Polymer Science, 116, 249 (1994).

Fixman, M., "Classical Statistical Mechanics of Constraints: A Theorem and Application to Polymers," Proc. Nat. Acad. Sci., 71, 3050 (1974).

Go, N.; Scheraga, H.A., "On the Use of Classical Statistical Mechanics in the Treatment of Polymer Chain Conformation," Macromolecules, 9, 535 (1976).

Metropolis, N.; Rosenbluth, A.W.; Rosenbluth, M.N.; Teller, A.H.; Teller, E. J. Chem. Phys., 21, 1087 (1953).

McQuarrie, D.A. Statistical Mechanics, Harper and Rowe (1976).



Last updated December 08, 1998 at 07:26PM Pacific Standard Time.
Copyright © 1998, Molecular Simulations, Inc. All rights reserved.