[MMTK] Defining a molecule using cartesian coordinates

[Harry Zuzan]

I have tried to define a few simple molecules and minimize them using
the AmberForceField.  Below is an example of my attempt to define a
methane molecule with the initial positions of atoms given in cartesian
coordinates.   The command lines are followed by the error messages
and my definition of the methane molecule is appended.

How should the definition of the cartesian coordinates be revised in
order to let the molecule be minimized?

Harry Zuzan


>>> from mmtk import *
>>> U = InfiniteUniverse(AmberForceField())
>>> U.molecule = Molecule('methane')
>>> m = SteepestDescentMinimizer(U)
>>> m(steps = 100)
Traceback (innermost last):
  File "<stdin>", line 1, in ?
  File "/usr/local/src/mmtk/mmtk/Minimization.py", line 42, in __call__
    evaluator = self.universe.energyEvaluator().evaluator.CEvaluator()
  File "/usr/local/src/mmtk/mmtk/Universe.py", line 237, in
energyEvaluator
    eval = ForceField.EnergyEvaluator(self, self._forcefield,
  File "/usr/local/src/mmtk/mmtk/ForceField.py", line 204, in __init__
    self.evaluator = self.ff.evaluator(self.universe,
  File "/usr/local/src/mmtk/mmtk/AmberForceField.py", line 313, in
evaluator
    return CompoundForceField.evaluator(self, universe, subset1,
subset2,
  File "/usr/local/src/mmtk/mmtk/ForceField.py", line 96, in evaluator
    eval_objects.append(ff.evaluator(system, subset1, subset2,
  File "/usr/local/src/mmtk/mmtk/AmberForceField.py", line 70, in
evaluator
    return BondedForceField.evaluator(self, universe, subset1, subset2,
  File "/usr/local/src/mmtk/mmtk/BondedInteractions.py", line 52, in
evaluator
    self.addBondAngleTerm(data, angle, bu, global_data)
  File "/usr/local/src/mmtk/mmtk/AmberForceField.py", line 102, in
addBondAngleTerm
    raise KeyError, 'No parameters for angle ' + `a1` +  \
KeyError: No parameters for angle Atom methane.H1--Atom methane.C--Atom
methane.H2


-------------------------------------
name = "methane"

structure = "CH4"

C = Atom('C')
H1 = Atom('H')
H2 = Atom('H')
H3 = Atom('H')
H4 = Atom('H')


bonds = [Bond(C, H1), Bond(C, H2), Bond(C, H3), Bond(C, H4)]

amber_atom_type = {C: 'CT', H1: 'H3', H2: 'H3', H3: 'H3', H4: 'H3'}

amber_charge = {C: 0.0, H1: 0.1, H2: 0.1, H3: 0.1, H4: 0.1}

configurations = {
'default': Cartesian({
C: (1.006, 0.581, 0.411),
H1: (1.006, 0.581, 0.257),
H2: (1.006, 0.436, 0.462),
H3: (0.880, 0.653, 0.462),
H4: (1.132, 0.653, 0.462)})
}