Molmodel
Public Types | Public Member Functions | Protected Member Functions | Friends

SimTK::DuMMForceFieldSubsystem Class Reference

This is a concrete subsystem that provides basic molecular mechanics functionality for coarse-grained molecules built in the SimTK framework. More...

#include <DuMMForceFieldSubsystem.h>

Inheritance diagram for SimTK::DuMMForceFieldSubsystem:

List of all members.

Public Types

enum  VdwMixingRule {
  WaldmanHagler = 1, HalgrenHHG = 2, Jorgensen = 3, LorentzBerthelot = 4,
  Kong = 5
}
 

These are the van der Waals mixing rules supported by DuMM.

More...

Public Member Functions

 DuMMForceFieldSubsystem ()
 DuMMForceFieldSubsystem (MolecularMechanicsSystem &)
void setTraceOpenMM (bool shouldTrace)
 OBSOLETE NAME -- just use setTracing().
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
 Calculate 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.
Atom/bond exclusion methods (advanced users only!)

Methods in this group give you fine control over which atoms in the defined molecules are actually used in the calculation of nonbonded forces, and which bonds are included in the calculation of bonded forces.

By default, any atom or bond is included if it can contribute to forces that affect motion.

Warning:
You can cause strange effects with these methods if you're not careful -- for example excluding atoms while using GBSA will cause odd surface area calculations to be performed. Don't use these methods unless you know exactly what you're doing!

During realizeTopology() DuMM studies the defined molecules and force field terms to determine (for example) which bonds cross bodies and which atoms can generate forces. Then the list of included atoms and bonds is finalized based on the instructions you give by calling these methods.

void clearIncludedNonbondAtomList ()
 Clear the list of atoms to be included in nonbonded force calculations. This is usually called prior to including selected groups of atoms.
void clearIncludedBondList ()
 Clear the list of bonds to be included in bonded force calculations. This is usually called prior to including selected bonds.
void resetIncludedNonbondAtomListToDefault ()
 Restore the included nonbond atoms list to its DuMM-selected default. This will include all the atoms that can affect interbody nonbonded forces, including Coulomb, Van der Waals, and solvent (GBSA) forces.
void resetIncludedBondListToDefault ()
 Request that DuMM set the included bonds list to all bonds that can generate forces that affect motion, meaning bonds whose bonded force terms have a non-zero amplitude and that involve atoms from two or more bodies.
void includeNonbondAtom (DuMM::AtomIndex atom)
 Add the given atom to the included nonbond atom list if it isn't already there. This does not include bonded terms that involve this atom; you have to request that explicitly with includeAllInterbodyBondsForOneAtom().
void includeAllNonbondAtomsForOneBody (MobilizedBodyIndex mobod)
 Add to the included nonbond atom list all the atoms that are fixed to the given body.
void includeAllInterbodyBondsForOneAtom (DuMM::AtomIndex atom)
 Given an atom, include all the bonded force terms that (a) include this atom, and (b) involve a body other than the one to which this atom is fixed. This does not include the atom in nonbonded force calculations; you have to request that explicitly with includeNonbondAtom().
void includeAllInterbodyBondsWithBothAtoms (DuMM::AtomIndex atom1, DuMM::AtomIndex atom2)
 Given two atoms that may appear together in some bonded force term, include all the bonded terms that involve both of them.
void includeAllInterbodyBondsWithBothAtoms (DuMM::BondIndex bond)
 Given a bond (that is, a 1-2 connection between atoms), extract the two atoms from it and then make sure all interbody bond terms that involve both of those atoms (in any position) are included.
void includeAllInterbodyBondsForOneBody (MobilizedBodyIndex mobod)
 Include any bonds that are needed to calculate interbody bonded forces that involve any atoms fixed to the given mobilized body. This will involve a few atoms on neighboring bodies.
void includeAllInterbodyBondsBetweenTwoBodies (MobilizedBodyIndex mobod1, MobilizedBodyIndex mobod2)
 Include any bonds that are needed to calculate interbody bonded forces that act between the given pair of mobilized bodies. Any bonded term that involves both an atom from mobod1 and an atom from mobod2 will be included.
int getNumIncludedAtoms () const
 Atoms to be included in force calculations are numbered from 0 to getNumIncludedAtoms()-1 after realizeTopology() has been called.
DuMM::AtomIndex getAtomIndexOfIncludedAtom (DuMM::IncludedAtomIndex incAtomIx) const
 Given a DuMM::IncludedAtomIndex, return the corresponding DuMM::AtomIndex. You must already have called realizeTopology().
int getNumNonbondAtoms () const
 A subset of the included atoms are used in nonbond calculations. These are numbered from 0 to getNumNonbondAtoms()-1 after realizeTopology() has been called.
DuMM::IncludedAtomIndex getIncludedAtomIndexOfNonbondAtom (DuMM::NonbondAtomIndex nonbondAtomIx) const
 Given a DuMM::NonbondAtomIndex, return the corresponding DuMM::IncludedAtomIndex. You must already have called realizeTopology().
DuMM::AtomIndex getAtomIndexOfNonbondAtom (DuMM::NonbondAtomIndex nonbondAtomIx) const
 Given a DuMM::NonbondAtomIndex, return the corresponding DuMM::AtomIndex. You must already have called realizeTopology().
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.

void defineAtomClass (DuMM::AtomClassIndex atomClassIx, const char *atomClassName, int atomicNumber, int expectedValence, Real vdwRadiusInNm, Real vdwWellDepthInKJ)
 Define a new atom class for this force field, for identifying atoms of a particular element, number of bonds, and van der Waals parameters.
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, that is, 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)
 Obsolete method -- use the other signature.
bool hasAtomClass (DuMM::AtomClassIndex) const
 Check whether an atom class has been defined using this index.
bool hasAtomClass (const String &atomClassName) const
 Check whether an atom class has been defined using this name.
DuMM::AtomClassIndex getAtomClassIndex (const String &atomClassName) const
 Obtain the atom class index corresponding to this atom class name.
DuMM::AtomClassIndex getNextUnusedAtomClassIndex () const
 Obtain an atom class index that is numerically larger than the largest currently-defined atom class index.
DuMM::AtomClassIndex getAtomClassIndex (DuMM::AtomIndex atomIx) const
 Get the index number of the atom class associated with this atom.
Real getVdwRadius (DuMM::AtomClassIndex atomClassIx) const
 Get the van der Waals radius shared by all atoms that belong to the indicated atom class. See comments for this group for a precise definition of what this means; there are ambiguities so don't assume you already know.
Real getVdwWellDepth (DuMM::AtomClassIndex atomClassIx) const
 Get the van der Waals energy well depth shared by all atoms that belong to the indicated atom class. See comments for this group for a precise definition of what this means; there are ambiguities so don't assume you already know.
void defineChargedAtomType (DuMM::ChargedAtomTypeIndex atomTypeIx, const char *atomTypeName, DuMM::AtomClassIndex atomClassIx, Real partialChargeInE)
 Define a new ChargedAtomType for this force field, for identifying atoms of a particular AtomClass that have a particular partial charge.
void defineChargedAtomType_KA (DuMM::ChargedAtomTypeIndex atomTypeIx, const char *atomTypeName, DuMM::AtomClassIndex atomClassIx, Real partialChargeInE)
 This is an alternate name for defineChargedAtomType() but since partial charge uses the same unit (e) in both the MD and KA unit systems, this is identical to defineChargedAtomType().
void defineChargedAtomType_KA (int atomTypeIx, const char *atomTypeName, int atomClassIx, Real partialChargeInE)
 Obsolete method -- use the other signature.
bool hasChargedAtomType (DuMM::ChargedAtomTypeIndex) const
 Check whether a charged atom type has been defined using this index.
bool hasChargedAtomType (const String &chargedTypeName) const
 Check whether a charged atom type has been defined using this name.
DuMM::ChargedAtomTypeIndex getChargedAtomTypeIndex (const String &chargedTypeName) const
 Obtain the charged atom type index corresponding to this charged atom type name.
DuMM::ChargedAtomTypeIndex getNextUnusedChargedAtomTypeIndex () const
 Obtain a charged atom type index that is numerically larger than the largest currently-defined charged atom type index.
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. (Only the Amber99 force field has been tested.)
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.

Protected Member Functions

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

Detailed Description

This is a concrete subsystem that provides basic molecular mechanics functionality for coarse-grained molecules built in the SimTK framework.

UNITS: This subsystem requires that the system be modeled in "MD units" of nanometers, daltons (g/mol), and picoseconds, yielding consistent energy units of kJ/mol==(Da-nm^2/ps^2), and force in kJ/mol-nm. Charge is in proton charge units e, and angles are in radians. For convenience, we allow the force field to be defined in "KA" units, that is, angstroms instead of nanometers, and energy in kcal rather than kJ, and we also allow angles to be supplied in degrees. However, these are immediately converted to the MD units described above.


Member Enumeration Documentation

These are the van der Waals mixing rules supported by DuMM.

Enumerator:
WaldmanHagler 

Our default, Waldman & Hagler, J.Comp.Chem. 14(9) 1993.

HalgrenHHG 

MMFF, AMOEBA.

Jorgensen 

OPLS.

LorentzBerthelot 

AMBER, CHARMM.

Kong 

Kong, J.Chem.Phys. 59(5) 1973.


The documentation for this class was generated from the following file:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines