Forcefield Class

The primary data structure in foyer is the Forcefield class, which inherits from the OpenMM class of the same name. The primary operation on this class is the .apply() function, which takes a chemical topology and returns a parametrized ParmEd Structure. The user may pass some otions to this function based on a particular use case.

class foyer.forcefield.Forcefield(forcefield_files=None, name=None, validation=True, debug=False)[source]

Specialization of OpenMM’s Forcefield allowing SMARTS based atomtyping.

Parameters
  • forcefield_files (list of str, optional, default=None) – List of forcefield files to load.

  • name (str, optional, default=None) – Name of a forcefield to load that is packaged within foyer.

apply(structure, references_file=None, use_residue_map=True, assert_bond_params=True, assert_angle_params=True, assert_dihedral_params=True, assert_improper_params=False, verbose=False, *args, **kwargs)[source]

Apply the force field to a molecular structure.

Parameters
  • structure (parmed.Structure or mbuild.Compound) – Molecular structure to apply the force field to

  • references_file (str, optional, default=None) – Name of file where force field references will be written (in Bibtex format)

  • use_residue_map (boolean, optional, default=True) – Store atomtyped topologies of residues to a dictionary that maps them to residue names. Each topology, including atomtypes, will be copied to other residues with the same name. This avoids repeatedly calling the subgraph isomorphism on idential residues and should result in better performance for systems with many identical residues, i.e. a box of water. Note that for this to be applied to independent molecules, they must each be saved as different residues in the topology.

  • assert_bond_params (bool, optional, default=True) – If True, Foyer will exit if parameters are not found for all system bonds.

  • assert_angle_params (bool, optional, default=True) – If True, Foyer will exit if parameters are not found for all system angles.

  • assert_dihedral_params (bool, optional, default=True) – If True, Foyer will exit if parameters are not found for all system proper dihedrals.

  • assert_improper_params (bool, optional, default=False) – If True, Foyer will exit if parameters are not found for all system improper dihedrals.

  • verbose (bool, optional, default=False) – If True, Foyer will print debug-level information about notable or potentially problematic details it encounters.

property combining_rule

Return the combining rule of this force field.

property coulomb14scale

Get Coulombic 1-4 scale for this forcefield.

createSystem(topology, nonbondedMethod=NoCutoff, nonbondedCutoff=Quantity(value=1.0, unit=nanometer), constraints=None, rigidWater=True, removeCMMotion=True, hydrogenMass=None, switchDistance=None, **args)[source]

Construct an OpenMM System representing a Topology with this force field.

Parameters
  • topology (Topology) – The Topology for which to create a System

  • nonbondedMethod (object=NoCutoff) – The method to use for nonbonded interactions. Allowed values are NoCutoff, CutoffNonPeriodic, CutoffPeriodic, Ewald, or PME.

  • nonbondedCutoff (distance=1*nanometer) – The cutoff distance to use for nonbonded interactions

  • constraints (object=None) – Specifies which bonds and angles should be implemented with constraints. Allowed values are None, HBonds, AllBonds, or HAngles.

  • rigidWater (boolean=True) – If true, water molecules will be fully rigid regardless of the value passed for the constraints argument

  • removeCMMotion (boolean=True) – If true, a CMMotionRemover will be added to the System

  • hydrogenMass (mass=None) – The mass to use for hydrogen atoms bound to heavy atoms. Any mass added to a hydrogen is subtracted from the heavy atom to keep their total mass the same.

  • switchDistance (float=None) – The distance at which the potential energy switching function is turned on for

  • args – Arbitrary additional keyword arguments may also be specified. This allows extra parameters to be specified that are specific to particular force fields.

Returns

the newly created System

Return type

system

static get_generator(ff, gen_type)[source]

Return a specific force generator for this Forcefield.

Parameters
  • ff (foyer.Forcefield) – The forcefield to return generator for

  • gen_type (Type) – The generator type to return

Returns

The instance of gen_type for this Forcefield

Return type

instance of gen_type

Raises

MissingForceError – If the Forcefield doesn’t have a generator of type gen_type

get_parameters(group, key, keys_are_atom_classes=False)[source]

Get parameters for a specific group of Forces in this Forcefield.

Parameters
  • group (str) – One of {“atoms”, “harmonic_bonds”, “harmonic_angles”, “periodic_propers”, “periodic_impropers”, “rb_propers”, “rb_impropers”}. Note that these entries are case insensitive

  • key (str, list of str) – The atom type(s)/class(es) to extract parameters for

  • keys_are_atom_classes (bool, default=False) – If True, the entries in key are considered to be atom classes rather than atom types

Examples

Following shows some example usage of this function:
>>> from foyer import Forcefield
>>> ff = Forcefield(name='oplsaa')
>>> ff.get_parameters('atoms', key='opls_235')  # NonBonded Parameters
{'charge': -0.18, 'sigma': 0.35, 'epsilon': 0.276144}
>>> ff.get_parameters('rb_propers', ['opls_137', 'opls_209', 'opls_193', 'opls_154'])  # RBPropers
{'c0': 2.87441, 'c1': 0.58158, 'c2': 2.092, 'c3': -5.54799, 'c4': 0.0, 'c5': 0.0}
Returns

A dictionary (or an item in the list of dictionaries) with parameter names as key and the parameter value as value as it pertains to specific force. Note that the keys of the dictionary can be different for different forces.

Return type

dict or list of dict

Raises
  • MissingParametersError – Raised when parameters are missing from the forcefield or no matching parameters found

  • MissingForceError – Raised when a particular force generator is missing from the Forcefield

property included_forcefields

Return a dictionary mappying for all included forcefields.

property lj14scale

Get LJ 1-4 scale for this forcefield.

map_atom_classes_to_types(atom_classes_keys, strict=False)[source]

Replace the atom_classes in the provided list by atom_types.

property name

Return the name of the force field XML.

parametrize_system(structure=None, references_file=None, assert_bond_params=True, assert_angle_params=True, assert_dihedral_params=True, assert_improper_params=False, verbose=False, *args, **kwargs)[source]

Create system based on resulting typemapping.

registerAtomType(parameters)[source]

Register a new atom type.

run_atomtyping(structure, use_residue_map=True, **kwargs)[source]

Atomtype the topology.

Parameters
  • structure (parmed.structure.Structure) – Molecular structure to find atom types of

  • use_residue_map (boolean, optional, default=True) – Store atomtyped topologies of residues to a dictionary that maps them to residue names. Each topology, including atomtypes, will be copied to other residues with the same name. This avoids repeatedly calling the subgraph isomorphism on idential residues and should result in better performance for systems with many identical residues, i.e. a box of water. Note that for this to be applied to independent molecules, they must each be saved as different residues in the topology.

static substitute_wildcards(atom_types, wildcard)[source]

Return possible wildcard options.

property version

Return version number of the force field XML file.