Molecular Mechanics in Molmodel

Once you have built a system of Compounds (molecules) in Molmodel, and modeled them with the degrees of freedom you want, you can apply molecular forces using DuMM, the Molmodel molecular mechanics force field. More...

Classes

class  CustomBondStretch
 Abstract base class for custom bond stretch terms, that is, functional forms that apply a distance-based force along a 1-2 bond between a pair of atoms. More...
class  DuMMForceFieldSubsystem
 This is a concrete subsystem that provides basic molecular mechanics functionality for coarse-grained molecules built in the SimTK framework. More...
class  Amber99ForceSubsystem
 This class is just a DuMMForceFieldSubsystem for which the constructor pre-loads the definitions of the Amber99 force field. More...

Namespaces

namespace  SimTK::DuMM
 

This namespace is used for symbols which are useful in conjunction with Molmodel's DuMMForceFieldSubsystem.


Modules

 User defined molecular mechanics force terms
 

These classes are used with DuMMForceFieldSubsystem to permit user-defined functional forms for bonded terms.


 Handy conversion constants
 

These are a set of multiplicative factors for use in converting among the various commonly used molecular mechanics units.


Functions

 SimTK_DEFINE_UNIQUE_INDEX_TYPE (AtomClassIndex)
 SimTK_DEFINE_UNIQUE_INDEX_TYPE (ChargedAtomTypeIndex)
 SimTK_DEFINE_UNIQUE_INDEX_TYPE (BondIndex)
 SimTK_DEFINE_UNIQUE_INDEX_TYPE (ClusterIndex)
virtual Real calcEnergy (Real distance) const =0
virtual Real calcForce (Real distance) const =0
virtual Real calcEnergy (Real bendAngle) const =0
virtual Real calcTorque (Real bendAngle) const =0
virtual Real calcEnergy (Real dihedralAngle) const =0
virtual Real calcTorque (Real dihedralAngle) const =0
 DuMMForceFieldSubsystem ()
 DuMMForceFieldSubsystem (MolecularMechanicsSystem &)
void setTraceOpenMM (bool shouldTrace)
 OBSOLETE NAME -- just use setTracing().
 SimTK_PIMPL_DOWNCAST (DuMMForceFieldSubsystem, ForceSubsystem)
 SimTK_PIMPL_DOWNCAST (DuMMForceFieldSubsystem, Subsystem)
void defineIncompleteAtomClass (DuMM::AtomClassIndex classIx, const char *name, int elementNumber, int valence)
void defineIncompleteAtomClass_KA (DuMM::AtomClassIndex classIx, const char *name, int elementNumber, int valence)
void setAtomClassVdwParameters (DuMM::AtomClassIndex atomClassIx, Real vdwRadiusInNm, Real vdwWellDepthInKJPerMol)
void setAtomClassVdwParameters_KA (DuMM::AtomClassIndex atomClassIx, Real radiusInAng, Real wellDepthInKcal)
bool isValidAtomClass (DuMM::AtomClassIndex) const
void defineIncompleteChargedAtomType (DuMM::ChargedAtomTypeIndex typeIx, const char *name, DuMM::AtomClassIndex classIx)
void defineIncompleteChargedAtomType_KA (DuMM::ChargedAtomTypeIndex typeIx, const char *name, DuMM::AtomClassIndex classIx)
void setChargedAtomTypeCharge (DuMM::ChargedAtomTypeIndex, Real charge)
void setChargedAtomTypeCharge_KA (DuMM::ChargedAtomTypeIndex chargedAtomTypeIx, Real charge)

Friends

class MolecularMechanicsSystem

Runtime methods

Methods in this group are used to obtain state-dependent information about this molecular system at run time.



const Vec3 & getAtomLocationInGround (const State &state, DuMM::AtomIndex atomIx)
 Obtain the current location of a particular atom in the Ground frame.
const Vec3 & getAtomVelocityInGround (const State &state, DuMM::AtomIndex atomIx)
 Obtain the current velocity of a particular atom in the Ground frame.
const Vec3 & getAtomBodyStationInGround (const State &state, DuMM::AtomIndex atomIx)
 Obtain the station vector of an atom measured from its body's origin, but re-expressed in the Ground frame.

Define particular molecules

Methods in this group are used to define the atoms and bonds in the particular set of molecules being simulated.



DuMM::AtomIndex addAtom (DuMM::ChargedAtomTypeIndex chargedAtomTypeIx)
 Add a new atom to the model.
DuMM::BondIndex addBond (DuMM::AtomIndex atom1Ix, DuMM::AtomIndex atom2Ix)
 Declare that there is a covalent bond between two atoms.
DuMM::AtomIndex getBondAtom (DuMM::BondIndex bond, int which) const
 For a given 1-2 bond, return the atoms which are connected by that bond.
int getNumAtoms () const
 How many atoms are currently in the model?
int getNumBonds () const
 How many 1-2 bonds are currently in the model?
Real getAtomMass (DuMM::AtomIndex atomIx) const
 Obtain the mass in Daltons (g/mol) of the atom indicated by the given AtomIndex.
int getAtomElement (DuMM::AtomIndex atomIx) const
 Obtain the element (by atomic number) of the atom indicated by the given AtomIndex.
Real getAtomRadius (DuMM::AtomIndex atomIx) const
 Obtain the van der Waals radius of the atom indicated by the given AtomIndex.
MobilizedBodyIndex getAtomBody (DuMM::AtomIndex atomIx) const
 Obtain the Simbody MobilizedBodyIndex of the rigid body on which a particular atom has been fixed.
Vec3 getAtomStationOnBody (DuMM::AtomIndex atomIx) const
 Obtain the station at which a particular atom is fixed on its body.
Vec3 getAtomStationInCluster (DuMM::AtomIndex atomIx, DuMM::ClusterIndex clusterIx) const
 Obtain the station at which a particular atom is fixed within a particular Cluster (an atom can be in more than one Cluster).
Vec3 getElementDefaultColor (int atomicNumber) const
 For display purposes, return the RGB value of a suggested color for an element given by atomic number.
Vec3 getAtomDefaultColor (DuMM::AtomIndex atomIx) const
 For display purposes, return the RGB value of a suggested color with which to display a particular atom.

Define clusters and bodies

Methods in this group control the grouping of atoms into rigid clusters and the placement of such clusters onto the rigid bodies of the underlying multibody system.

Note: we use the term "station" to refer to a fixed location with respect to a cluster or body frame, that is, a point "stationary" on that cluster or body.



DuMM::ClusterIndex createCluster (const char *clusterName)
 Create an empty Cluster (rigid group of atoms).
void placeAtomInCluster (DuMM::AtomIndex atomIx, DuMM::ClusterIndex clusterIx, const Vec3 &station)
 Place an existing atom at a particular station in the local frame of a Cluster.
void placeClusterInCluster (DuMM::ClusterIndex childClusterIndex, DuMM::ClusterIndex parentClusterIndex, const Transform &X_PC)
 Place a Cluster (the child) in another Cluster (the parent).
MassProperties calcClusterMassProperties (DuMM::ClusterIndex clusterIx, const Transform &X_BC=Transform()) const
 Calcuate the composite mass properties of a Cluster, either in its own reference frame C or in reference frame B with the Cluster placed relative to B using the indicated Transform X_BC.
void attachClusterToBody (DuMM::ClusterIndex clusterIx, MobilizedBodyIndex body, const Transform &X_BC=Transform())
 Place a Cluster's local frame C at a particular location and orientation with respect to a MobilizedBody's frame B.
void attachAtomToBody (DuMM::AtomIndex atomIx, MobilizedBodyIndex body, const Vec3 &station=Vec3(0))
 Place an individual atom at a particular station on a body without an intervening cluster.
MobilizedBodyIndex getClusterBody (DuMM::ClusterIndex clusterIx) const
 Find the MobilizedBody on which a Cluster has been fixed in place.
Transform getClusterPlacementOnBody (DuMM::ClusterIndex clusterIx) const
 Find where on its body a Cluster has been placed by returning the transform X_BC giving the orientation and position of Cluster frame C in body frame B.
Transform getClusterPlacementInCluster (DuMM::ClusterIndex childClusterIndex, DuMM::ClusterIndex parentClusterIndex) const
 Find where on parent cluster P a child cluster C has been placed, by returning the transform X_PC.

Define atom categories

An AtomClass is used to collect together a set of properties which are expected to be shared by many individual atoms.

The properties are: the element (as atomic number), expected valence, and van der Waals parameters. Charge is not included in AtomClass but in a second more detailed classification level called ChargedAtomType.

This fails if the atom class already exists. Van der waals radius is given as Rmin, the radius at which the energy well minimum is seen (actually it is 1/2 the distance between atom centers for a pair of atoms of this class). This is *not* Sigma, which we define as the radius (half distance) at which the energy crosses zero, that is, a little closer together than when the energy well is at maximum depth. CAUTION: it is also common in force fields that these terms are given by diameter rather than radius -- be sure you know what convention was use in your source so that you can convert to radius here. To convert for LJ: Rmin = 2^(1/6) * Sigma. The radius is in nm, the well depth in kJ/mol.



void defineAtomClass (DuMM::AtomClassIndex atomClassIx, const char *atomClassName, int elementNumber, int expectedValence, Real vdwRadiusInNm, Real vdwWellDepthInKJ)
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
void defineAtomClass_KA (DuMM::AtomClassIndex atomClassIx, const char *atomClassName, int element, int valence, Real vdwRadiusInAng, Real vdwWellDepthInKcal)
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
void defineAtomClass_KA (int atomClassIx, const char *atomClassName, int element, int valence, Real vdwRadiusInAng, Real vdwWellDepthInKcal)
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
bool hasAtomClass (DuMM::AtomClassIndex) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
bool hasAtomClass (const String &atomClassName) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
DuMM::AtomClassIndex getAtomClassIndex (const String &atomClassName) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
DuMM::AtomClassIndex getNextUnusedAtomClassIndex () const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
DuMM::AtomClassIndex getAtomClassIndex (DuMM::AtomIndex atomIx) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
Real getVdwRadius (DuMM::AtomClassIndex atomClassIx) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
Real getVdwWellDepth (DuMM::AtomClassIndex atomClassIx) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
void defineChargedAtomType (DuMM::ChargedAtomTypeIndex atomTypeIx, const char *atomTypeName, DuMM::AtomClassIndex atomClassIx, Real partialChargeInE)
 PartialCharge in units of e (charge on a proton); same in MD & KA.
void defineChargedAtomType_KA (DuMM::ChargedAtomTypeIndex atomTypeIx, const char *atomTypeName, DuMM::AtomClassIndex atomClassIx, Real partialChargeInE)
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
void defineChargedAtomType_KA (int atomTypeIx, const char *atomTypeName, int atomClassIx, Real partialChargeInE)
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
bool hasChargedAtomType (DuMM::ChargedAtomTypeIndex) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
bool hasChargedAtomType (const String &chargedTypeName) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
DuMM::ChargedAtomTypeIndex getChargedAtomTypeIndex (const String &chargedTypeName) const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.
DuMM::ChargedAtomTypeIndex getNextUnusedChargedAtomTypeIndex () const
 Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

Control force field nonbonded behavior in special circumstances

These methods permit setting overall force field behavior in special circumstances, including van der Waals mixing behavior for dissimilar atom pairs and scaling of non-bonded terms for closely-bonded atoms.



const char * getVdwMixingRuleName (VdwMixingRule) const
 Obtain a human-readable name for one of our van der Waals mixing rules.
void setVdwMixingRule (VdwMixingRule)
 Set the van der Waals mixing rule -- our default is Waldman-Hagler.
VdwMixingRule getVdwMixingRule () const
 Get the van der Waals mixing rule currently in effect.
void setVdw12ScaleFactor (Real)
 default 0
void setVdw13ScaleFactor (Real)
 default 0
void setVdw14ScaleFactor (Real)
 default 1
void setVdw15ScaleFactor (Real)
 default 1
void setCoulomb12ScaleFactor (Real)
 default 0
void setCoulomb13ScaleFactor (Real)
 default 0
void setCoulomb14ScaleFactor (Real)
 default 1
void setCoulomb15ScaleFactor (Real)
 default 1

Tinker biotypes and pre-defined force field parameter sets

DuMM understands Tinker-format parameter files that can be used to load in a whole force field description.

This requires assigning Tinker Biotypes to the atoms in the molecule.



void loadAmber99Parameters ()
 Use Amber99 force field parameters.
void populateFromTinkerParameterFile (std::istream &)
 Load force field parameters from a TINKER format force field parameter file.
void setBiotypeChargedAtomType (DuMM::ChargedAtomTypeIndex chargedAtomTypeIndex, BiotypeIndex biotypeIx)
 Associate a Tinker Biotype with a ChargedAtomType in this subsystem.
DuMM::ChargedAtomTypeIndex getBiotypeChargedAtomType (BiotypeIndex biotypeIx) const
 Get charged atom type index in this force field associated with a particular Biotype.
std::ostream & generateBiotypeChargedAtomTypeSelfCode (std::ostream &os) const
 Generate C++ code from the current contents of this DuMM force field object.

Bond stretch terms

Bond stretch parameters (between 2 atom classes).

You can use the standard, built-in functional form (a harmonic) or define your own.



void defineBondStretch (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, Real stiffnessInKJperNmSq, Real nominalLengthInNm)
 Define a harmonic bond stretch term between two atom classes using the built-in functional form.
void defineBondStretch_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, Real stiffnessInKcalPerAngSq, Real nominalLengthInAng)
 Same as defineBondStretch() except that for convenience this takes stiffness in (kcal/mol)/A^2, and nominal length is in A (angstroms).
void defineBondStretch_KA (int class1, int class2, Real stiffnessInKcalPerAngSq, Real nominalLengthInAng)
 Same as defineBondStretch_KA() but takes integer class arguments for backwards compatibility.
void defineCustomBondStretch (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::CustomBondStretch *bondStretchTerm)
 Define a custom bond stretch term to be applied to the indicated pair of atom classes whenever they are found in a 1-2 bond.

Bond bending terms

Bond bending parameters (for 3 atom classes bonded 1-2-3).

You can use the standard, built-in functional form (a harmonic based on the 1-2-3 angle) or define your own.



void defineBondBend (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, Real stiffnessInKJPerRadSq, Real nominalAngleInDeg)
 Define a harmonic bond bending term applying to a triple of atom classes using the built-in functional form.
void defineBondBend_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, Real stiffnessInKcalPerRadSq, Real nominalAngleInDeg)
 Same as defineBondBend() except that for convenience this takes stiffness in (kcal/mol)/rad^2 (nominal angle is still in degrees).
void defineBondBend_KA (int class1, int class2, int class3, Real stiffnessInKcalPerRadSq, Real nominalAngleInDeg)
 Same as defineBondBend_KA() but takes integer class arguments for backwards compatibility.
void defineCustomBondBend (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::CustomBondBend *bondBendTerm)
 Define a custom bond bend term to be applied to the indicated triple of atom classes whenever they are found in a 1-2-3 bonded sequence.

Bond torsion terms

Bond torsion (dihedral) parameters (for 4 atom classes bonded 1-2-3-4).

You can use the standard, built-in functional form (a combination of sinusoids) or define your own. Bond torsion terms produce energy and a scalar torque that are dependent only on the rotation angle about the 2-3 bond in the specified atom class sequence.

Theory

Label the four bonded atoms r-x-y-s. Rotation occurs about the axis v=y-x, that is, a vector from x to y. We define a torsion angle theta using the "polymer convention" rather than the IUPAC one which is 180 degrees different. Ours is like this:

             r                         r      s
   theta=0    \             theta=180   \    / 
               x--y                      x--y
                   \
                    s
 

The sign convention is the same for IUPAC and polymer: A positive angle is defined by considering r-x fixed in space. Then using the right hand rule around v (that is, thumb points from x to y) a positive rotation rotates y->s in the direction of your fingers.

We use a periodic energy function like this:

       E(theta) = sum E_n(1 + cos(n*theta - theta0_n))
 

where n is the periodicity, E_n is the amplitude (kJ/mol) for term n, and theta0_n is the phase offset for term n. The torque term (applied about the v axis) is then

       T(theta) = -[sum -n*E_n*sin(n*theta - theta0_n)]
 

Note that in the functional forms above the amplitude paramter E_n could be considered a "half amplitude" since the total excursion of the energy and torque functions is 2*E_n. When entering this parameter, be sure you understand whether the data you have represents the (half) amplitude as above (most common) or the full excursion range of the function, in which case you should provide only 1/2 the excursion as the value for DuMM's amplitude parameter.

For a custom torsion bond, DuMM will provide the angle theta to the user-written energy and torque routines, expecting a scalar value E(theta) and T(theta) back as above. However, the functional form is specified by the user in that case, and it is the user's responsibility to ensure that the returned torque is the negative gradient of the energy with respect to the angle parameter theta.

In either case, DuMM takes care of translating the returned scalar torque into an appropriate mobility torque when possible, or an appropriate set of forces on each of the atoms r-x-y-s such that the desired pure torque is realized as the resultant of the forces.



void defineBondTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity, Real ampInKJ, Real phaseInDegrees)
 Define bond torsion terms applying to a quadruple of atom classes using the built-in functional form, a combination of sinusoids with different periods (see full discussion in the group header for bond torsion terms).
void defineBondTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKJ, Real phase1InDegrees, int periodicity2, Real amp2InKJ, Real phase2InDegrees)
 Same as defineBondTorsion() but permits two torsion terms (with different periods) to be specified simultaneously.
void defineBondTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKJ, Real phase1InDegrees, int periodicity2, Real amp2InKJ, Real phase2InDegrees, int periodicity3, Real amp3InKJ, Real phase3InDegrees)
 Same as defineBondTorsion() but permits three torsion terms (with different periods) to be specified simultaneously.
void defineBondTorsion_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees)
 Same as defineBondTorsion() but permits takes the amplitude in kcal/mol (but note that this is converted immediately to our MD unit system of kJ/mol).
void defineBondTorsion_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees, int periodicity2, Real amp2InKcal, Real phase2InDegrees)
 Same as defineBondTorsion_KA() but permits two torsion terms (with different periods) to be specified simultaneously.
void defineBondTorsion_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees, int periodicity2, Real amp2InKcal, Real phase2InDegrees, int periodicity3, Real amp3InKcal, Real phase3InDegrees)
 Same as defineBondTorsion_KA() but permits three torsion terms (with different periods) to be specified simultaneously.
void defineBondTorsion_KA (int class1, int class2, int class3, int class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees)
 Same as defineBondTorsion_KA() but takes integer class arguments for backwards compatibility.
void defineBondTorsion_KA (int class1, int class2, int class3, int class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees, int periodicity2, Real amp2InKcal, Real phase2InDegrees)
 Same as defineBondTorsion_KA() but takes integer class arguments for backwards compatibility.
void defineBondTorsion_KA (int class1, int class2, int class3, int class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees, int periodicity2, Real amp2InKcal, Real phase2InDegrees, int periodicity3, Real amp3InKcal, Real phase3InDegrees)
 Same as defineBondTorsion_KA() but takes integer class arguments for backwards compatibility.
void defineCustomBondTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, DuMM::CustomBondTorsion *bondTorsionTerm)
 Define a custom bond torsion term to be applied to the indicated quadruple of atom classes whenever they are found in a 1-2-3-4 bonded sequence.

Amber-style improper torsions

As with normal torsions, (see defineBondTorsion()), only one term may have a given periodicity.

The amplitudes are in kJ/mol. The third atom is the central one to which the other three are bonded; this is not the same in reverse order.



void defineAmberImproperTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity, Real ampInKJ, Real phaseInDegrees)
 Provide one torsion term in MD units, using kilojoules for amplitude.
void defineAmberImproperTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKJ, Real phase1InDegrees, int periodicity2, Real amp2InKJ, Real phase2InDegrees)
 Provide two torsion terms in MD units, using kilojoules/mole for amplitude.
void defineAmberImproperTorsion (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKJ, Real phase1InDegrees, int periodicity2, Real amp2InKJ, Real phase2InDegrees, int periodicity3, Real amp3InKJ, Real phase3InDegrees)
 Provide three torsion terms in MD units, using kilojoules/mole for amplitude.
void defineAmberImproperTorsion_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees)
 Provide one torsion term in KA units, using kilocalories/mole for amplitude.
void defineAmberImproperTorsion_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees, int periodicity2, Real amp2InKcal, Real phase2InDegrees)
 Provide two torsion terms in KA units, using kilocalories/mole for amplitude.
void defineAmberImproperTorsion_KA (DuMM::AtomClassIndex class1, DuMM::AtomClassIndex class2, DuMM::AtomClassIndex class3, DuMM::AtomClassIndex class4, int periodicity1, Real amp1InKcal, Real phase1InDegrees, int periodicity2, Real amp2InKcal, Real phase2InDegrees, int periodicity3, Real amp3InKcal, Real phase3InDegrees)
 Provide three torsion terms in KA units, using kilocalories/mole for amplitude.

GBSA implicit solvation

These methods are used to set GBSA terms.



void setSolventDielectric (Real)
void setSoluteDielectric (Real)
Real getSolventDielectric () const
Real getSoluteDielectric () const
void setGbsaIncludeAceApproximation (bool)
void setGbsaIncludeAceApproximationOn ()
void setGbsaIncludeAceApproximationOff ()

Global scale factors

These non-physical parameters can be used to weaken or disable (or magnify) individual force field terms.

These are always 1 for correct implementation of any force field; other values are primarily useful for testing the effects of individual terms on results or performance. Set to 0 to disable the corresponding term altogether.



void setVdwGlobalScaleFactor (Real)
 scale all van der Waals terms
void setCoulombGlobalScaleFactor (Real)
 scale all Coulomb terms
void setGbsaGlobalScaleFactor (Real)
 scale all GBSA terms
void setBondStretchGlobalScaleFactor (Real)
 scale all built-in bond stretch terms
void setBondBendGlobalScaleFactor (Real)
 scale all built-in bond bending terms
void setBondTorsionGlobalScaleFactor (Real)
 scale all built-in bond torsion terms
void setAmberImproperTorsionGlobalScaleFactor (Real)
 scale all improper torsion terms
void setCustomBondStretchGlobalScaleFactor (Real)
 scale all custom bond stretch terms
void setCustomBondBendGlobalScaleFactor (Real)
 scale all custom bond bending terms
void setCustomBondTorsionGlobalScaleFactor (Real)
 scale all custom bond torsion terms
Real getVdwGlobalScaleFactor () const
 get current scale factor for van der Waals terms
Real getCoulombGlobalScaleFactor () const
 get current scale factor for Coulomb terms
Real getGbsaGlobalScaleFactor () const
 get current scale factor for GBSA terms
Real getBondStretchGlobalScaleFactor () const
 get current scale factor for built-in bond stretch terms
Real getBondBendGlobalScaleFactor () const
 get current scale factor for built-in bond bending terms
Real getBondTorsionGlobalScaleFactor () const
 get current scale factor for built-in bond torsion terms
Real getAmberImproperTorsionGlobalScaleFactor () const
 get current scale factor for improper torsion terms
Real getCustomBondStretchGlobalScaleFactor () const
 get current scale factor for custom bond stretch terms
Real getCustomBondBendGlobalScaleFactor () const
 get current scale factor for custom bond bending terms
Real getCustomBondTorsionGlobalScaleFactor () const
 get current scale factor for custom bond torsion terms
void setAllGlobalScaleFactors (Real s)
 Set all the global scale factors to the same value.

Computational options

These methods control how DuMM performs its computations.



void setTracing (bool)
 For debugging, you can ask DuMM to dump some information to std::clog about its attempts to figure out the best way to use the available hardware.
void setUseMultithreadedComputation (bool)
 Enable or disable multithreaded computation (enabled by default).
bool getUseMultithreadedComputation () const
 Is multithreaded computation enabled?
void setNumThreadsRequested (int)
 Request how many threads we want to use if multithreading is enabled; zero (the default) means let DuMM choose, which will likely be one thread per processor.
int getNumThreadsRequested () const
 What was the last value passed to setNumThreadsRequested()? The default is zero meaning DuMM chooses the number of threads.
bool isUsingMultithreadedComputation () const
 Is DuMM using the multithreaded code? This could return true even if there is just one thread, if you forced it with setNumThreadsToUse().
int getNumThreadsInUse () const
 Find out how many threads DuMM is actually using; will be zero until after realizeTopology().
void setUseOpenMMAcceleration (bool)
 This determines whether we use OpenMM GPU acceleration if it is available.
bool getUseOpenMMAcceleration () const
 Return the current setting of the flag set by setUseOpenMMAcceleration().
void setAllowOpenMMReference (bool)
 This allows us to use OpenMM even if only the Reference platform is available.
bool getAllowOpenMMReference () const
 Return the current setting of the flag set by setAllowOpenMMReference().
bool isUsingOpenMM () const
 Return true if DuMM is currently using OpenMM for its computations.
std::string getOpenMMPlatformInUse () const
 Return the OpenMM Platform currently in use, or the empty string if we're not using OpenMM.

Bookkeeping, debugging, and internal-use-only methods

Hopefully you won't need these.



long long getForceEvaluationCount () const
 How many times has the forcefield been evaluated?
void dump () const
 Produce an ugly but comprehensive dump of the contents of DuMM's internal data structures, sent to std::cout (stdout).
void dumpCForceFieldParameters (std::ostream &os, const String &methodName="loadParameters") const
 Generate C++ code to reproduce forceField parameters presently in memory.
void loadTestMoleculeParameters ()
 Load test parameters.

Detailed Description

Once you have built a system of Compounds (molecules) in Molmodel, and modeled them with the degrees of freedom you want, you can apply molecular forces using DuMM, the Molmodel molecular mechanics force field.


Function Documentation

DuMM::AtomIndex addAtom ( DuMM::ChargedAtomTypeIndex  chargedAtomTypeIx  )  [inherited]

Add a new atom to the model.

The AtomIndex number is returned; you don't get to pick your own. Use the AtomIndex to identify this particular atom subsequently.

DuMM::BondIndex addBond ( DuMM::AtomIndex  atom1Ix,
DuMM::AtomIndex  atom2Ix 
) [inherited]

Declare that there is a covalent bond between two atoms.

Note that these are AtomIndex numbers, not AtomClasses or ChargedAtomTypes.

void attachAtomToBody ( DuMM::AtomIndex  atomIx,
MobilizedBodyIndex  body,
const Vec3 station = Vec3(0) 
) [inherited]

Place an individual atom at a particular station on a body without an intervening cluster.

void attachClusterToBody ( DuMM::ClusterIndex  clusterIx,
MobilizedBodyIndex  body,
const Transform X_BC = Transform() 
) [inherited]

Place a Cluster's local frame C at a particular location and orientation with respect to a MobilizedBody's frame B.

All the atoms within the Cluster will become fixed to the body while maintaining the same relative positions as they had in the Cluster.

MassProperties calcClusterMassProperties ( DuMM::ClusterIndex  clusterIx,
const Transform X_BC = Transform() 
) const [inherited]

Calcuate the composite mass properties of a Cluster, either in its own reference frame C or in reference frame B with the Cluster placed relative to B using the indicated Transform X_BC.

virtual Real calcEnergy ( Real  dihedralAngle  )  const [pure virtual, inherited]
virtual Real calcEnergy ( Real  bendAngle  )  const [pure virtual, inherited]
virtual Real calcEnergy ( Real  distance  )  const [pure virtual, inherited]
virtual Real calcForce ( Real  distance  )  const [pure virtual, inherited]
virtual Real calcTorque ( Real  dihedralAngle  )  const [pure virtual, inherited]
virtual Real calcTorque ( Real  bendAngle  )  const [pure virtual, inherited]
DuMM::ClusterIndex createCluster ( const char *  clusterName  )  [inherited]

Create an empty Cluster (rigid group of atoms).

The Cluster index number is returned; you don't get to pick your own. The name is just for display; you must use the index to reference the Cluster. Every Cluster has its own reference frame C.

void defineAmberImproperTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKJ,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKJ,
Real  phase2InDegrees,
int  periodicity3,
Real  amp3InKJ,
Real  phase3InDegrees 
) [inherited]

Provide three torsion terms in MD units, using kilojoules/mole for amplitude.

void defineAmberImproperTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKJ,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKJ,
Real  phase2InDegrees 
) [inherited]

Provide two torsion terms in MD units, using kilojoules/mole for amplitude.

void defineAmberImproperTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity,
Real  ampInKJ,
Real  phaseInDegrees 
) [inherited]

Provide one torsion term in MD units, using kilojoules for amplitude.

void defineAmberImproperTorsion_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKcal,
Real  phase2InDegrees,
int  periodicity3,
Real  amp3InKcal,
Real  phase3InDegrees 
) [inline, inherited]

Provide three torsion terms in KA units, using kilocalories/mole for amplitude.

void defineAmberImproperTorsion_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKcal,
Real  phase2InDegrees 
) [inline, inherited]

Provide two torsion terms in KA units, using kilocalories/mole for amplitude.

void defineAmberImproperTorsion_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees 
) [inline, inherited]

Provide one torsion term in KA units, using kilocalories/mole for amplitude.

void defineAtomClass ( DuMM::AtomClassIndex  atomClassIx,
const char *  atomClassName,
int  elementNumber,
int  expectedValence,
Real  vdwRadiusInNm,
Real  vdwWellDepthInKJ 
) [inline, inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

void defineAtomClass_KA ( int  atomClassIx,
const char *  atomClassName,
int  element,
int  valence,
Real  vdwRadiusInAng,
Real  vdwWellDepthInKcal 
) [inline, inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

void defineAtomClass_KA ( DuMM::AtomClassIndex  atomClassIx,
const char *  atomClassName,
int  element,
int  valence,
Real  vdwRadiusInAng,
Real  vdwWellDepthInKcal 
) [inline, inherited]
void defineBondBend ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
Real  stiffnessInKJPerRadSq,
Real  nominalAngleInDeg 
) [inherited]

Define a harmonic bond bending term applying to a triple of atom classes using the built-in functional form.

You are only allowed to specify a built-in term once per class triple, regardless of the order they are listed. It is permitted to redefine the same term as long as the parameters are the same both times, otherwise this method will throw an exception if (class1,class2,class3) or (class3,class2,class1) has already been assigned. Stiffness (energy per angle^2) must be given in (kJ/mol)/rad^2, and nominal bond angle is in degrees(!). Note that the nominal angle is in degrees while the stiffness is in radians. Odd, I know, but that seems to be how it's done! CAUTION: Energy is ka^2 using this definition, while torque is -2ka; note factor of 2 in torque -- conventions vary.

Parameters:
[in] class1,class2,class3 The atom class triple to which this term applies; 1-2-3 or 3-2-1 order doesn't matter.
[in] stiffnessInKJPerRadSq The bond stiffness in kilojoules per radian^2.
[in] nominalAngleInDeg Bond angle at which energy and force are zero, in degrees(!).
See also:
defineBondBend_KA() for kcal-angstrom units
defineCustomBondBend() to define your own functional form
void defineBondBend_KA ( int  class1,
int  class2,
int  class3,
Real  stiffnessInKcalPerRadSq,
Real  nominalAngleInDeg 
) [inline, inherited]

Same as defineBondBend_KA() but takes integer class arguments for backwards compatibility.

References DuMMForceFieldSubsystem::defineBondBend_KA().

Referenced by DuMMForceFieldSubsystem::defineBondBend_KA().

void defineBondBend_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
Real  stiffnessInKcalPerRadSq,
Real  nominalAngleInDeg 
) [inline, inherited]

Same as defineBondBend() except that for convenience this takes stiffness in (kcal/mol)/rad^2 (nominal angle is still in degrees).

Note that the stiffness is immediately converted to our standard MD units internally.

See also:
defineBondBend() for a complete description
Parameters:
[in] class1,class2,class3 The atom class triple to which this term applies; 1-2-3 or 3-2-1 order doesn't matter.
[in] stiffnessInKcalPerRadSq The bond stiffness in kilocalories per radian^2.
[in] nominalAngleInDeg Bond angle at which energy and force are zero, in degrees(!).
void defineBondStretch ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
Real  stiffnessInKJperNmSq,
Real  nominalLengthInNm 
) [inherited]

Define a harmonic bond stretch term between two atom classes using the built-in functional form.

You are only allowed to specify a built-in term once per pair of classes, regardless of the order they are listed. It is permitted to redefine the same term as long as the parameters are the same both times, otherwise this method will throw an exception if (class1,class2) or (class2,class1) has already been assigned. Stiffness (energy per length^2) must be given in (kJ/mol)/nm^2, and nominal bond length is in nanometers. CAUTION: Energy is kx^2 using this definition, while force is -2kx; note factor of 2 in force -- conventions vary.

Parameters:
[in] class1,class2 The atom class pair to which this term applies; order doesn't matter.
[in] stiffnessInKJperNmSq The bond stiffness in kilojoules per nm^2.
[in] nominalLengthInNm Bond length at which energy and force are zero, in nm.
See also:
defineBondStretch_KA() for kcal-angstrom units
defineCustomBondStretch() to define your own functional form
void defineBondStretch_KA ( int  class1,
int  class2,
Real  stiffnessInKcalPerAngSq,
Real  nominalLengthInAng 
) [inline, inherited]

Same as defineBondStretch_KA() but takes integer class arguments for backwards compatibility.

References DuMMForceFieldSubsystem::defineBondStretch_KA().

Referenced by DuMMForceFieldSubsystem::defineBondStretch_KA().

void defineBondStretch_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
Real  stiffnessInKcalPerAngSq,
Real  nominalLengthInAng 
) [inline, inherited]

Same as defineBondStretch() except that for convenience this takes stiffness in (kcal/mol)/A^2, and nominal length is in A (angstroms).

Note that these are immediately converted to our standard MD units internally.

See also:
defineBondStretch() for a complete description
Parameters:
[in] class1,class2 The atom class pair to which this term applies; order doesn't matter.
[in] stiffnessInKcalPerAngSq The bond stiffness in kilocalories per Angstrom^2.
[in] nominalLengthInAng Bond length at which energy and force are zero, in Angstrom.

References SimTK::square().

void defineBondTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKJ,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKJ,
Real  phase2InDegrees,
int  periodicity3,
Real  amp3InKJ,
Real  phase3InDegrees 
) [inherited]

Same as defineBondTorsion() but permits three torsion terms (with different periods) to be specified simultaneously.

void defineBondTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKJ,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKJ,
Real  phase2InDegrees 
) [inherited]

Same as defineBondTorsion() but permits two torsion terms (with different periods) to be specified simultaneously.

void defineBondTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity,
Real  ampInKJ,
Real  phaseInDegrees 
) [inherited]

Define bond torsion terms applying to a quadruple of atom classes using the built-in functional form, a combination of sinusoids with different periods (see full discussion in the group header for bond torsion terms).

You are only allowed to specify a built-in term of a given period once per class quadruple, whether specified 1-2-3-4 or 4-3-2-1. It is permitted to redefine the same term as long as the parameters are the same both times, otherwise this method will throw an exception if (class1,class2,class3,class4) or (class4,class3,class2,class1) has already been assigned. The amplitude of the energy functions must be given in (kJ/mol), the period is an integer which will multiply the dihedral angle, and the phase (offset) angle is given in degrees (not radians). CAUTION: various conventions exist for specifying the parameters for bond torsion terms. Pay particular attention to the amplitude -- in our convention it is really a half-amplitude since the sinusoids range from -1 to 1 making the full energy and torque "excursion" twice the amplitude.

Parameters:
[in] class1,class2,class3,class4 The atom class quadruple to which this term applies; 1-2-3-4 or 4-3-2-1 order doesn't matter.
[in] periodicity The periodicity (dependence on input angle theta) of the sinusoidal term.
[in] ampInKJ The amplitude (half-range) of the energy function, in kilojoules.
[in] phaseInDegrees Angle at which periodicity*theta yields maximum energy (2*ampInKJ).
See also:
defineBondTorsion_KA() for kcal-angstrom units
defineCustomBondTorsion() to define your own functional form
void defineBondTorsion_KA ( int  class1,
int  class2,
int  class3,
int  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKcal,
Real  phase2InDegrees,
int  periodicity3,
Real  amp3InKcal,
Real  phase3InDegrees 
) [inline, inherited]

Same as defineBondTorsion_KA() but takes integer class arguments for backwards compatibility.

void defineBondTorsion_KA ( int  class1,
int  class2,
int  class3,
int  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKcal,
Real  phase2InDegrees 
) [inline, inherited]

Same as defineBondTorsion_KA() but takes integer class arguments for backwards compatibility.

void defineBondTorsion_KA ( int  class1,
int  class2,
int  class3,
int  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees 
) [inline, inherited]

Same as defineBondTorsion_KA() but takes integer class arguments for backwards compatibility.

void defineBondTorsion_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKcal,
Real  phase2InDegrees,
int  periodicity3,
Real  amp3InKcal,
Real  phase3InDegrees 
) [inline, inherited]

Same as defineBondTorsion_KA() but permits three torsion terms (with different periods) to be specified simultaneously.

void defineBondTorsion_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees,
int  periodicity2,
Real  amp2InKcal,
Real  phase2InDegrees 
) [inline, inherited]

Same as defineBondTorsion_KA() but permits two torsion terms (with different periods) to be specified simultaneously.

void defineBondTorsion_KA ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
int  periodicity1,
Real  amp1InKcal,
Real  phase1InDegrees 
) [inline, inherited]

Same as defineBondTorsion() but permits takes the amplitude in kcal/mol (but note that this is converted immediately to our MD unit system of kJ/mol).

void defineChargedAtomType ( DuMM::ChargedAtomTypeIndex  atomTypeIx,
const char *  atomTypeName,
DuMM::AtomClassIndex  atomClassIx,
Real  partialChargeInE 
) [inline, inherited]
void defineChargedAtomType_KA ( int  atomTypeIx,
const char *  atomTypeName,
int  atomClassIx,
Real  partialChargeInE 
) [inline, inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

void defineChargedAtomType_KA ( DuMM::ChargedAtomTypeIndex  atomTypeIx,
const char *  atomTypeName,
DuMM::AtomClassIndex  atomClassIx,
Real  partialChargeInE 
) [inline, inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

void defineCustomBondBend ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::CustomBondBend bondBendTerm 
) [inherited]

Define a custom bond bend term to be applied to the indicated triple of atom classes whenever they are found in a 1-2-3 bonded sequence.

Pass in a pointer to a newly- allocated CustomBondBend object; the subsystem takes over ownership of the object so it should NOT be deleted by the caller.

Parameters:
[in] class1,class2,class3 The atom class triple to which this term applies; 1-2-3 or 3-2-1 order doesn't matter.
[in] bondBendTerm Pointer to a heap-allocated object of a concrete type derived from DuMM::CustomBondBend.
See also:
DuMM::CustomBondBend
defineBondBend() to use the built-in harmonic functional form
void defineCustomBondStretch ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::CustomBondStretch bondStretchTerm 
) [inherited]

Define a custom bond stretch term to be applied to the indicated pair of atom classes whenever they are found in a 1-2 bond.

Pass in a pointer to a newly- allocated CustomBondStretch object; the subsystem takes over ownership of the object so it should NOT be deleted by the caller.

Parameters:
[in] class1,class2 The atom class pair to which this term applies; order doesn't matter.
[in] bondStretchTerm Pointer to a heap-allocated object of a concrete type derived from DuMM::CustomBondStretch.
See also:
DuMM::CustomBondStretch
defineBondStretch() to use the built-in harmonic functional form
void defineCustomBondTorsion ( DuMM::AtomClassIndex  class1,
DuMM::AtomClassIndex  class2,
DuMM::AtomClassIndex  class3,
DuMM::AtomClassIndex  class4,
DuMM::CustomBondTorsion bondTorsionTerm 
) [inherited]

Define a custom bond torsion term to be applied to the indicated quadruple of atom classes whenever they are found in a 1-2-3-4 bonded sequence.

Pass in a pointer to a newly- allocated CustomBondTorsion object; the subsystem takes over ownership of the object so it should NOT be deleted by the caller.

Parameters:
[in] class1,class2,class3,class4 The atom class triple to which this term applies; 1-2-3-4 or 4-3-2-1 order doesn't matter.
[in] bondTorsionTerm Pointer to a heap-allocated object of a concrete type derived from DuMM::CustomBondTorsion.
See also:
DuMM::CustomBondTorsion
defineBondTorsion() to use the built-in sinusoidal functional form
void defineIncompleteAtomClass ( DuMM::AtomClassIndex  classIx,
const char *  name,
int  elementNumber,
int  valence 
) [protected, inherited]
void defineIncompleteAtomClass_KA ( DuMM::AtomClassIndex  classIx,
const char *  name,
int  elementNumber,
int  valence 
) [inline, protected, inherited]
void defineIncompleteChargedAtomType ( DuMM::ChargedAtomTypeIndex  typeIx,
const char *  name,
DuMM::AtomClassIndex  classIx 
) [protected, inherited]
void defineIncompleteChargedAtomType_KA ( DuMM::ChargedAtomTypeIndex  typeIx,
const char *  name,
DuMM::AtomClassIndex  classIx 
) [inline, protected, inherited]
DuMMForceFieldSubsystem ( MolecularMechanicsSystem  )  [explicit, inherited]
DuMMForceFieldSubsystem (  )  [inherited]
void dump (  )  const [inherited]

Produce an ugly but comprehensive dump of the contents of DuMM's internal data structures, sent to std::cout (stdout).

void dumpCForceFieldParameters ( std::ostream &  os,
const String methodName = "loadParameters" 
) const [inherited]

Generate C++ code to reproduce forceField parameters presently in memory.

std::ostream& generateBiotypeChargedAtomTypeSelfCode ( std::ostream &  os  )  const [inherited]

Generate C++ code from the current contents of this DuMM force field object.

bool getAllowOpenMMReference (  )  const [inherited]

Return the current setting of the flag set by setAllowOpenMMReference().

Real getAmberImproperTorsionGlobalScaleFactor (  )  const [inherited]

get current scale factor for improper torsion terms

MobilizedBodyIndex getAtomBody ( DuMM::AtomIndex  atomIx  )  const [inherited]

Obtain the Simbody MobilizedBodyIndex of the rigid body on which a particular atom has been fixed.

An exception will be thrown if this atom is not fixed to any body.

const Vec3& getAtomBodyStationInGround ( const State state,
DuMM::AtomIndex  atomIx 
) [inherited]

Obtain the station vector of an atom measured from its body's origin, but re-expressed in the Ground frame.

That is, this is the vector which, if added to the spatial location of the atom body's origin, would yield the atom's spatial location (as returned by getAtomLocationInGround()). This requires that state has been realized to Stage::Position.

DuMM::AtomClassIndex getAtomClassIndex ( DuMM::AtomIndex  atomIx  )  const [inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

DuMM::AtomClassIndex getAtomClassIndex ( const String atomClassName  )  const [inherited]
Vec3 getAtomDefaultColor ( DuMM::AtomIndex  atomIx  )  const [inherited]

For display purposes, return the RGB value of a suggested color with which to display a particular atom.

int getAtomElement ( DuMM::AtomIndex  atomIx  )  const [inherited]

Obtain the element (by atomic number) of the atom indicated by the given AtomIndex.

const Vec3& getAtomLocationInGround ( const State state,
DuMM::AtomIndex  atomIx 
) [inherited]

Obtain the current location of a particular atom in the Ground frame.

This requires that state has been realized to Stage::Position.

Real getAtomMass ( DuMM::AtomIndex  atomIx  )  const [inherited]

Obtain the mass in Daltons (g/mol) of the atom indicated by the given AtomIndex.

Real getAtomRadius ( DuMM::AtomIndex  atomIx  )  const [inherited]

Obtain the van der Waals radius of the atom indicated by the given AtomIndex.

Vec3 getAtomStationInCluster ( DuMM::AtomIndex  atomIx,
DuMM::ClusterIndex  clusterIx 
) const [inherited]

Obtain the station at which a particular atom is fixed within a particular Cluster (an atom can be in more than one Cluster).

An exception will be thrown if this atom is not fixed to the cluster.

Vec3 getAtomStationOnBody ( DuMM::AtomIndex  atomIx  )  const [inherited]

Obtain the station at which a particular atom is fixed on its body.

An exception will be thrown if this atom is not fixed to any body.

const Vec3& getAtomVelocityInGround ( const State state,
DuMM::AtomIndex  atomIx 
) [inherited]

Obtain the current velocity of a particular atom in the Ground frame.

This requires that state has been realized to Stage::Velocity.

DuMM::ChargedAtomTypeIndex getBiotypeChargedAtomType ( BiotypeIndex  biotypeIx  )  const [inherited]

Get charged atom type index in this force field associated with a particular Biotype.

DuMM::AtomIndex getBondAtom ( DuMM::BondIndex  bond,
int  which 
) const [inherited]

For a given 1-2 bond, return the atoms which are connected by that bond.

You select one atom at a time by setting parameter which to 0 or 1. 0 will return the lower-numbered AtomIndex, regardless of the order in which the atoms were specified to addBond().

Real getBondBendGlobalScaleFactor (  )  const [inherited]

get current scale factor for built-in bond bending terms

Real getBondStretchGlobalScaleFactor (  )  const [inherited]

get current scale factor for built-in bond stretch terms

Real getBondTorsionGlobalScaleFactor (  )  const [inherited]

get current scale factor for built-in bond torsion terms

DuMM::ChargedAtomTypeIndex getChargedAtomTypeIndex ( const String chargedTypeName  )  const [inherited]
MobilizedBodyIndex getClusterBody ( DuMM::ClusterIndex  clusterIx  )  const [inherited]

Find the MobilizedBody on which a Cluster has been fixed in place.

Transform getClusterPlacementInCluster ( DuMM::ClusterIndex  childClusterIndex,
DuMM::ClusterIndex  parentClusterIndex 
) const [inherited]

Find where on parent cluster P a child cluster C has been placed, by returning the transform X_PC.

Transform getClusterPlacementOnBody ( DuMM::ClusterIndex  clusterIx  )  const [inherited]

Find where on its body a Cluster has been placed by returning the transform X_BC giving the orientation and position of Cluster frame C in body frame B.

Real getCoulombGlobalScaleFactor (  )  const [inherited]

get current scale factor for Coulomb terms

Real getCustomBondBendGlobalScaleFactor (  )  const [inherited]

get current scale factor for custom bond bending terms

Real getCustomBondStretchGlobalScaleFactor (  )  const [inherited]

get current scale factor for custom bond stretch terms

Real getCustomBondTorsionGlobalScaleFactor (  )  const [inherited]

get current scale factor for custom bond torsion terms

Vec3 getElementDefaultColor ( int  atomicNumber  )  const [inherited]

For display purposes, return the RGB value of a suggested color for an element given by atomic number.

For example, if the atomicNumber is 8 (Oxygen) the suggested color will be Red (1,0,0).

long long getForceEvaluationCount (  )  const [inherited]

How many times has the forcefield been evaluated?

Real getGbsaGlobalScaleFactor (  )  const [inherited]

get current scale factor for GBSA terms

DuMM::AtomClassIndex getNextUnusedAtomClassIndex (  )  const [inherited]
DuMM::ChargedAtomTypeIndex getNextUnusedChargedAtomTypeIndex (  )  const [inherited]
int getNumAtoms (  )  const [inherited]

How many atoms are currently in the model?

int getNumBonds (  )  const [inherited]

How many 1-2 bonds are currently in the model?

int getNumThreadsInUse (  )  const [inherited]

Find out how many threads DuMM is actually using; will be zero until after realizeTopology().

int getNumThreadsRequested (  )  const [inherited]

What was the last value passed to setNumThreadsRequested()? The default is zero meaning DuMM chooses the number of threads.

This doesn't tell you how many threads are actually in use -- see getNumThreadsInUse() for that.

std::string getOpenMMPlatformInUse (  )  const [inherited]

Return the OpenMM Platform currently in use, or the empty string if we're not using OpenMM.

Real getSoluteDielectric (  )  const [inherited]
Real getSolventDielectric (  )  const [inherited]
bool getUseMultithreadedComputation (  )  const [inherited]

Is multithreaded computation enabled?

bool getUseOpenMMAcceleration (  )  const [inherited]

Return the current setting of the flag set by setUseOpenMMAcceleration().

Real getVdwGlobalScaleFactor (  )  const [inherited]

get current scale factor for van der Waals terms

VdwMixingRule getVdwMixingRule (  )  const [inherited]

Get the van der Waals mixing rule currently in effect.

See also:
getVdwMixingRuleName() for a human-readable version
const char* getVdwMixingRuleName ( VdwMixingRule   )  const [inherited]

Obtain a human-readable name for one of our van der Waals mixing rules.

Real getVdwRadius ( DuMM::AtomClassIndex  atomClassIx  )  const [inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

Real getVdwWellDepth ( DuMM::AtomClassIndex  atomClassIx  )  const [inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

bool hasAtomClass ( const String atomClassName  )  const [inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

bool hasAtomClass ( DuMM::AtomClassIndex   )  const [inherited]
bool hasChargedAtomType ( const String chargedTypeName  )  const [inherited]

Same routine as defineAtomClass() but in Kcal/Angstrom (KA) unit system, i.e., radius (still not sigma) is in Angstroms, and well depth in kcal/mol.

bool hasChargedAtomType ( DuMM::ChargedAtomTypeIndex   )  const [inherited]
bool isUsingMultithreadedComputation (  )  const [inherited]

Is DuMM using the multithreaded code? This could return true even if there is just one thread, if you forced it with setNumThreadsToUse().

bool isUsingOpenMM (  )  const [inherited]

Return true if DuMM is currently using OpenMM for its computations.

bool isValidAtomClass ( DuMM::AtomClassIndex   )  const [protected, inherited]
void loadAmber99Parameters (  )  [inherited]

Use Amber99 force field parameters.

This duplicates the Tinker Amber99 parameter set in pre-built code, so you don't need to load the parameters from a file.

Referenced by Amber99ForceSubsystem::Amber99ForceSubsystem().

void loadTestMoleculeParameters (  )  [inherited]

Load test parameters.

void placeAtomInCluster ( DuMM::AtomIndex  atomIx,
DuMM::ClusterIndex  clusterIx,
const Vec3 station 
) [inherited]

Place an existing atom at a particular station in the local frame of a Cluster.

It is fine for an atom to be in more than one Cluster as long as only one of them ends up attached to a body.

void placeClusterInCluster ( DuMM::ClusterIndex  childClusterIndex,
DuMM::ClusterIndex  parentClusterIndex,
const Transform X_PC 
) [inherited]

Place a Cluster (the child) in another Cluster (the parent).

The child's local frame C is placed at a given Transform with respect to the parent's frame P. All the atoms in the child Cluster maintain their relative positioning.

void populateFromTinkerParameterFile ( std::istream &   )  [inherited]

Load force field parameters from a TINKER format force field parameter file.

(Only the Amber99 force field has been tested.)

void setAllGlobalScaleFactors ( Real  s  )  [inline, inherited]

Set all the global scale factors to the same value.

This is commonly used to turn everything on or off, followed by selectively disabling or enabling individual terms.

void setAllowOpenMMReference ( bool   )  [inherited]

This allows us to use OpenMM even if only the Reference platform is available.

This is for testing/debugging; one should never use the Reference platform in production since it will likely be slower (and less accurate) than the CPU implementation.

void setAmberImproperTorsionGlobalScaleFactor ( Real   )  [inherited]

scale all improper torsion terms

void setAtomClassVdwParameters ( DuMM::AtomClassIndex  atomClassIx,
Real  vdwRadiusInNm,
Real  vdwWellDepthInKJPerMol 
) [protected, inherited]
void setAtomClassVdwParameters_KA ( DuMM::AtomClassIndex  atomClassIx,
Real  radiusInAng,
Real  wellDepthInKcal 
) [inline, protected, inherited]
void setBiotypeChargedAtomType ( DuMM::ChargedAtomTypeIndex  chargedAtomTypeIndex,
BiotypeIndex  biotypeIx 
) [inherited]

Associate a Tinker Biotype with a ChargedAtomType in this subsystem.

A Biotype association is required for every Biotype found on atoms in the current System.

Parameters:
chargedAtomTypeIndex Prexisting charged atom type index in this subsystem
biotypeIx Preexisting BiotypeIndex defined in the Biotype class

Referenced by P12::P12(), ZincIon::setAmberLikeParameters(), ChlorideIon::setAmberLikeParameters(), CalciumIon::setAmberLikeParameters(), LithiumIon::setAmberLikeParameters(), CesiumIon::setAmberLikeParameters(), MagnesiumIon::setAmberLikeParameters(), PotassiumIon::setAmberLikeParameters(), SodiumIon::setAmberLikeParameters(), RubidiumIon::setAmberLikeParameters(), and Water::Water().

void setBondBendGlobalScaleFactor ( Real   )  [inherited]

scale all built-in bond bending terms

void setBondStretchGlobalScaleFactor ( Real   )  [inherited]

scale all built-in bond stretch terms

void setBondTorsionGlobalScaleFactor ( Real   )  [inherited]

scale all built-in bond torsion terms

void setChargedAtomTypeCharge ( DuMM::ChargedAtomTypeIndex  ,
Real  charge 
) [protected, inherited]
void setChargedAtomTypeCharge_KA ( DuMM::ChargedAtomTypeIndex  chargedAtomTypeIx,
Real  charge 
) [inline, protected, inherited]
void setCoulomb12ScaleFactor ( Real   )  [inherited]

default 0

void setCoulomb13ScaleFactor ( Real   )  [inherited]

default 0

void setCoulomb14ScaleFactor ( Real   )  [inherited]

default 1

void setCoulomb15ScaleFactor ( Real   )  [inherited]

default 1

void setCoulombGlobalScaleFactor ( Real   )  [inherited]

scale all Coulomb terms

void setCustomBondBendGlobalScaleFactor ( Real   )  [inherited]

scale all custom bond bending terms

void setCustomBondStretchGlobalScaleFactor ( Real   )  [inherited]

scale all custom bond stretch terms

void setCustomBondTorsionGlobalScaleFactor ( Real   )  [inherited]

scale all custom bond torsion terms

void setGbsaGlobalScaleFactor ( Real   )  [inherited]

scale all GBSA terms

void setGbsaIncludeAceApproximation ( bool   )  [inherited]
void setGbsaIncludeAceApproximationOff (  )  [inline, inherited]
void setGbsaIncludeAceApproximationOn (  )  [inline, inherited]
void setNumThreadsRequested ( int   )  [inherited]

Request how many threads we want to use if multithreading is enabled; zero (the default) means let DuMM choose, which will likely be one thread per processor.

DuMM will choose to use single threaded code if there is only one processor, but requesting 1 thread here explicitly will cause it to use multithreaded code but with a single thread (useful for debugging sometimes).

void setSoluteDielectric ( Real   )  [inherited]
void setSolventDielectric ( Real   )  [inherited]
void setTraceOpenMM ( bool  shouldTrace  )  [inline, inherited]

OBSOLETE NAME -- just use setTracing().

void setTracing ( bool   )  [inherited]

For debugging, you can ask DuMM to dump some information to std::clog about its attempts to figure out the best way to use the available hardware.

void setUseMultithreadedComputation ( bool   )  [inherited]

Enable or disable multithreaded computation (enabled by default).

void setUseOpenMMAcceleration ( bool   )  [inherited]

This determines whether we use OpenMM GPU acceleration if it is available.

By default, this is set false because OpenMM will compute only to single precision. Note that even if you set this flag, we won't use OpenMM unless (a) it is installed correctly on your machine, and (b) it can run with GPU acceleration. If you want to allow use of the non-accelerated Reference Platform provided by OpenMM, use setAllowOpenMMReference().

void setVdw12ScaleFactor ( Real   )  [inherited]

default 0

void setVdw13ScaleFactor ( Real   )  [inherited]

default 0

void setVdw14ScaleFactor ( Real   )  [inherited]

default 1

void setVdw15ScaleFactor ( Real   )  [inherited]

default 1

void setVdwGlobalScaleFactor ( Real   )  [inherited]

scale all van der Waals terms

void setVdwMixingRule ( VdwMixingRule   )  [inherited]

Set the van der Waals mixing rule -- our default is Waldman-Hagler.

SimTK::DuMM::SimTK_DEFINE_UNIQUE_INDEX_TYPE ( ClusterIndex   ) 
SimTK::DuMM::SimTK_DEFINE_UNIQUE_INDEX_TYPE ( BondIndex   ) 
SimTK::DuMM::SimTK_DEFINE_UNIQUE_INDEX_TYPE ( ChargedAtomTypeIndex   ) 
SimTK::DuMM::SimTK_DEFINE_UNIQUE_INDEX_TYPE ( AtomClassIndex   ) 
SimTK_PIMPL_DOWNCAST ( DuMMForceFieldSubsystem  ,
Subsystem   
) [protected, inherited]

Reimplemented from ForceSubsystem.

SimTK_PIMPL_DOWNCAST ( DuMMForceFieldSubsystem  ,
ForceSubsystem   
) [protected, inherited]

Friends

friend class MolecularMechanicsSystem [friend, inherited]

Generated on Thu Aug 12 16:37:47 2010 for SimTKcore by  doxygen 1.6.1