Roger Edberg and Drew Whitehouse Supercomputer Facility Australian National University
Visualization techniques are used to follow the progress of two types of polymer calculations that are performed on the ANU's VP2200 supercomputer. First, we look at a simulated annealing calculation that unravels a random walk on a tetrahedral lattice and produces a realistic polymer conformation with excluded volume. Second, we take the resulting polymer and simulate its motion by molecular dynamics.
Random walk conformations of any length are quite easy to generate on the computer. On the other hand, conformations with excluded volume can be difficult to generate for even a small molecule of say 20 atoms (except for a few trivial conformations such as all-trans.)
Here, we start with a random walk
conformation on a tetrahedral lattice, and proceed to optimize it,
i.e. to remove all of the excluded volume interactions, by a simulated
annealing [1] procedure. The function to be optimized is the 'cost
function':
where n0 is the number of self-intersections, n1 is the number of nearest- neighbor contacts, and n2 is the number of next-nearest-neighbor contacts. The optimization algorithm is:
Steps 2-5 are iterated until a conformation with C=0 is obtained. Note
that the quantity T plays the role that temperature does in the Metropolis
Monte Carlo algorithm; here it controls the acceptance of moves that result
in increased energy. For our definition of the cost function, the algorithm
should converge for T<1.0.
The annealing calculation shown on film was performed for a 200 atom
molecule at T=0.8. Colors indicate the localization of the cost function
within the molecule, i.e. excluded volume contacts. The unravelling process
is not a physically realistic one because portions of the molecule may pass
through one another as bond vectors are changed. Some additional work is
necessary to improve the convergence of the algorithm; the cost function
decreases quickly in the initial stages of the calculation, while in latter
stages it takes much time to unravel the last few 'knots' in the molecule.
Once we have the unravelled conformation we can use it as a starting configuration for a molecular dynamics simulation. We choose a model with harmonic atom-atom bonds and Lennard-Jones (12-6) non-bonded interactions for atoms more than four apart on the chain. The velocity version of the Verlet algorithm is used to integrate the equations of motion [2]. Atomic velocities were randomized with an initial kinetic energy per particle of 3.0 (reduced units; initial temperature=2.0), and the system was allowed to evolve isoenergetically.
The film shows the molecule relaxing once the lattice geometry restriction is removed. The Lennard-Jones non-bonded potential causes the molecule to relax into a globular conformation. The mean squared radius of gyration decreases from its initial value of 26 to about 5.
In future work we plan to apply similar visualization techniques to Monte
Carlo and molecular dynamics simulations of polymeric liquids and mixtures.
[2] M.P. Allen and D. Tildesley, "Computer Simulation of Liquids", Clarendon Press, Oxford (1987)