OPLSPotentialBuilder¶
- class OPLSPotentialBuilder(use_bonds=None, include_bond_stretching=None, include_bond_angles=None, include_torsions=None, include_inversions=None, include_lennard_jones=None, include_electrostatic=None, lennard_jones_cutoff=None, lennard_jones_smoothing_length=None, electrostatic_pairwise_cutoff=None, electrostatic_smoothing_length=None, electrostatic_spme_cutoff=None, electrostatic_accuracy=None, electrostatic_scale=None, use_electrostatic_direct_sum=None, ignore_missing_types=None)¶
This object builds a OPLS potential using TremoloX.
- Parameters:
use_bonds (bool) – Controls if any bonding (stretching, angle, torsion and inversion) terms are to be included. Default: True.
include_bond_stretching (bool) – Controls if bonding stretching terms are included. Default: use_bonds.
include_bond_angles (bool) – Controls if angle bending terms are included. Default: use_bonds.
include_torsions (bool) – Controls if bond torsion terms are included. Default: use_bonds.
include_inversions (bool) – Controls if bond inversion torsion terms are included. Default: use_bonds.
include_lennard_jones (bool) – Controls if Lennard-Jones terms are included. Default: use_bonds.
include_electrostatic (bool) – Controls if electrostatic terms are included. Default: use_bonds.
lennard_jones_cutoff (PhysicalQuantity of type distance) – Distance beyond which the Lennard-Jones terms are truncated. Default: 10.0 * Angstrom.
lennard_jones_smoothing_length (PhysicalQuantity of type distance) – Distance over which the Lennard-Jones terms are brought to zero at the cutoff. Default: 2.0 * Angstrom.
electrostatic_pairwise_cutoff (PhysicalQuantity of type distance) – Distance at which the pairwise electrostatic terms are brought to zero in pairwise summation. Default: 12 * Angstrom
electrostatic_smoothing_length (PhysicalQuantity of type distance) – Distance over which the electrostatic terms are brought to zero at the cutoff when using the direct electrostatic sum. Default: 2.0 * Angstrom
electrostatic_spme_cutoff (PhysicalQuantity of type distance) – Distance at which the pairwise electrostatic terms are brought to zero in SPME summation. Default: 7.5 * Angstrom
electrostatic_accuracy (float) – Sets the relative accuracy of the SPME summation. Default: 0.0001
electrostatic_scale (float) – Set the relative strength of the electrostatic terms. Default: 1.0
use_electrostatic_direct_sum (bool) – Use the direct Coulomb sum to calculate the electrostatic terms in MoleculeConfigurations. Default: False
ignore_missing_types (bool) – Whether or not molecules missing OPLS types are ignored. Default: False
- assignAtomTypes(configuration, search_types=None, overwrite=False)¶
Automatically assign OPLS atom tags to a configuration.
- Parameters:
configuration (AtomicConfiguration) – The configuration on which to assign tags.
search_types (list of str | None) – List of types to assign onto the configuration. If specified only these types will be assigned, with OPLS_UNKNOWN assigned to atoms that don’t match any of the given types. The default argument None means that all available types are assigned. Default: None.
overwrite (bool) – Whether or not to assign new types if the configuration already has a valid set of types. Default: False
- checkAtomTypes(configuration, add_missing_terms=False)¶
Function that tests that each atom has only one valid OPLS type.
- Parameters:
configuration (AtomicConfiguration) – The configuration that is being tested.
add_missing_terms (bool) – Whether or not to add terms for missing OPLS types. Default: False
- createCalculator(configuration, check_types=True, verbose=False, set_improper_torsion_indices=True, atomic_charges=None, charge_averaging=None)¶
Create and return a calculator with the OPLS potential.
- Parameters:
configuration (AtomicConfiguration) – The configuration the potential is to be used with.
check_types (bool) – Whether or not the types are checked for basic consistency.
verbose (bool) – Whether or not to print the terms assigned in the potential set. Default: False
set_improper_torsion_indices (bool) – Whether or not to add improper torsions onto the input configuration. Default: True
atomic_charges (PhysicalQuantity of type charge) – Custom atomic charges to be used in the potential.
charge_averaging (NLFlag | None) – The type of averaging of charges to be used, if any. Accepted values are None, AtomConnectivity to average over atoms bonded to the same elements and AtomTypes to average over atoms of the same forcefield type. Default: None
- Returns:
The calculator containing the potential.
- Return type:
- electrostaticAccuracy()¶
- Returns:
The relative accuracy of the electrostatic summation for SPME.
- Return type:
float
- electrostaticPairwiseCutoff()¶
- Returns:
The cutoff used in pairwise electrostatic summation.
- Return type:
PhysicalQuantity of type length
- electrostaticScale()¶
- Returns:
Scaling factor for the electrostatic terms.
- Return type:
float
- electrostaticSmoothingLength()¶
- Returns:
Distance over which electrostatic terms are splined when using direct summation.
- Return type:
PhysicalQuantity of type length
- electrostaticSpmeCutoff()¶
- Returns:
The cutoff used in SPME electrostatic summation.
- Return type:
PhysicalQuantity of type length
- ignoreMissingTypes()¶
- Returns:
Whether or not missing types are ignored on molecules.
- Return type:
bool
- includeBondAngles()¶
- Returns:
Whether or not bond angle terms are included in the potential set.
- Return type:
bool
- includeBondStretching()¶
- Returns:
Whether or not bond stretching terms are included in the potential set.
- Return type:
bool
- includeElectrostatic()¶
- Returns:
Whether or not electrostatic terms are included in the potential set.
- Return type:
bool
- includeInversions()¶
- Returns:
Whether or not bond inversion terms are included in the potential set.
- Return type:
bool
- includeLennardJones()¶
- Returns:
Whether or not Lennard Jones terms are included in the potential set.
- Return type:
bool
- includeTorsions()¶
- Returns:
Whether or not bond torsion terms are included in the potential set.
- Return type:
bool
- isOplsType(type_label)¶
Function that tests if a type is in the potential.
- Parameters:
type_label (str) – The type being tested.
- Returns:
Whether or not the tag is defined in the potential.
- Return type:
bool
- lennardJonesCutoff()¶
- Returns:
Distance beyond which Lennard Jones terms are truncated.
- Return type:
PhysicalQuantity of type length
- lennardJonesSmoothingLength()¶
- Returns:
Distance over which Lennard Jones terms are splined.
- Return type:
PhysicalQuantity of type length
- oplsType(index, configuration, add_missing_terms=False)¶
Gives the OPLS type on an atom in a configuration, ignoring other tags.
- Parameters:
index (int) – Index of the atom for which the type is sought.
configuration (AtomicConfiguration) – The configuration containing the atom.
add_missing_terms (bool) – Whether or not to add terms for missing OPLS types. Default: False
- Returns:
The type of the atom.
- Return type:
str
- printTotalTypes(type_class=None)¶
Print a table of available atom types. This provides a list of types and descriptions which can be useful when doing manual type assignment.
- Parameters:
type_class (str | None) – Display only atom types with this class label. None displays all types. Default: None
- printUsedTypes(configuration)¶
Print a table of the types used on the configuration.
- Parameters:
configuration (AtomicConfiguration) – The configuration whose types are being displayed.
- totalCharge(configuration, verbose=False)¶
Return the total charge of the configuration according to the atom types on the molecule.
- Parameters:
configuration (AtomicConfiguration) – The configuration whose charges are being reported.
verbose (bool) – Flag as to whether or not additional information about charges is printed. Default: False
- Returns:
The total charge of the configuration according to the atom types on the molecule.
- Return type:
PhysicalQuantity with units of charge.
- typePartialCharge(type_label)¶
Return the atomic partial charge associated with the given OPLS type label.
- Parameters:
type_label (str) – The type label for the OPLS atom type.
- Returns:
The atom charge associated with the type label.
- Return type:
PhysicalQuantity of type charge
- uniqueString()¶
Return a unique string representing the state of the object.
- useBonds()¶
- Returns:
Whether or not bond terms are included in the potential set.
- Return type:
bool
- useElectrostaticDirectSum()¶
- Returns:
Whether or not direct Coulomb sums are used for molecules.
- Return type:
bool
Usage Examples¶
Build a bulk configuration of poly(vinyl chloride) plasticized with dimethyl sebacate using the OPLS potential.
# ---------------------------------------------
# Monomer And Included Molecules
# ---------------------------------------------
PVC_monomer = nlread('PVC_Monomer.hdf5')[-1]
dimthyl_sebacate = nlread('Dimethyl_Sebacate.hdf5')[-1]
# -------------------------------------------------------------
# OPLS Potential Definition
# -------------------------------------------------------------
opls_potential_builder = OPLSPotentialBuilder()
# Display the tags and charges placed on the dimethyl sebacate molecule
opls_potential_builder.printUsedTypes(dimthyl_sebacate)
opls_potential_builder.totalCharge(dimthyl_sebacate, verbose=True)
# ---------------------------------------------
# Polymer Sequence
# ---------------------------------------------
polymer_sequence = PolymerSequence(
monomer_configurations=[PVC_monomer],
number_of_monomers=20,
number_of_chains=20,
tactic_ratio=0.5,
)
# -------------------------------------------------------------
# Polymer Monte Carlo Builder
# -------------------------------------------------------------
polymer_builder = PolymerMonteCarloBuilder(
polymer_sequence=[polymer_sequence],
density=1.4000*gram/cm**3,
included_molecules=[(dimthyl_sebacate, 4)],
monte_carlo_temperature=300.0*Kelvin,
angle_sampling_points=20,
repack_molecules=False,
)
bulk_configuration = polymer_builder.buildPolymerConfiguration()
# -------------------------------------------------------------
# Add OPLS Potential To Polymer Configuration
# -------------------------------------------------------------
calculator = opls_potential_builder.createCalculator(bulk_configuration)
bulk_configuration.setCalculator(calculator)
nlsave('PVC_DMS.hdf5', bulk_configuration)
# -------------------------------------------------------------
# Force Capped Equilibration
# -------------------------------------------------------------
equilibration_method = ForceCappedEquilibration(
temperature=300.00*Kelvin,
md_steps_per_force_capped_simulation=40000,
md_time_step=0.50*fs,
force_capped_simulations=4,
starting_factor=1.040,
ending_factor=0.800,
fixed_indices=None,
)
bulk_configuration = equilibration_method.runEquilibration(bulk_configuration)
nlsave('PVC_DMS.hdf5', bulk_configuration)
OPLS_PVC_DMS_Example.py
PVC_Monomer.hdf5
Dimethyl_Sebacate.hdf5
The given example above demonstrates both how to create a calculator that uses the OPLS-AA potential and also how to assign tags for the dimethyl sebacate molecule that is included in the polymer melt. Both the atom types and charges are printed. This shows both the appropriate types that have been selected and also that the overall charge of the molecule is zero.
Notes¶
The OPLSPotentialBuilder class enables the creation of a TremoloXCalculator
that
implements the OPLS-AA potential of Mayo, Olafson and Goddard[1]. This is a bonded
valence forcefield that represents the potential energy of a molecule as a collection of simple
functions based on bond lengths, bond angles, torsion angles, inversion angles and inter-atomic
distances. The overall mathematical form of the OPLS-AA potential can be given as:
Here \(r\), \(\theta\), \(\phi\) and \(\chi\) represent inter-atomic and bond distances, angles, torsions and inversions respectively in the configuration.
To create a OPLS-AA calculator, an instance of the OPLSPotentialBuilder
must first
be created. The constructor for this class has a number of options that controls the way in which
the calculator built. Individual types of interactions can be excluded by setting the include_
options to False
By default all possible terms are included. All valence terms can also be
excluded by setting the argument use_bonds
to False. This is appropriate for configurations
consisting of single atoms that are not bonded to each other.
As the OPLSPotentialBuilder
class uses the bonding in the configuration to determine
appropriate potentials, bonds should be set on the configuration. This can be done with the
findBonds()
method in either the MoleculeConfiguration
or
BulkConfiguration
class. This automatically assigns bonds between atoms that are
within the covalent radii of each other. Specific bonds can also be assigned with the
setBonds()
, which takes a list of the indices of pairs of atoms that are
to be considered bonded.
There are then a number of options that control the non-bonding interactions. The parameter
lennard_jones_cutoff
sets the distance beyond which Lennard-Jones interactions are truncated. To
avoid integration problems in molecular dynamics, Lennard-Jones interactions are brought to zero at
the cutoff using a spline function. The width of this spline can be set using the
lennard_jones_smoothing_length
parameter. There are no specific hydrogen hydrogen bonding terms
defined within the OPLS-AA potential.
Electrostatic interactions are a little more complicated. By default the
OPLSPotentialBuilder
adds a CoulombSPME
electrostatic solver for
calculators to be used with a bulk configuration, and a CoulombDSF
electrostatic solver
to calculators to be used with a molecule configuration. When adding an CoulombSPME
solver the cutoff used for calculating the real-space interactions can be set using the
electrostatic_spme_cutoff
argument. Likewise the accuracy of the smoothed particle mesh Ewald
solver can be set using the electrostatic_accuracy
argument. In the case where a
CoulombDSF
solver is added to the potential, the cutoff used can be set using the
electrostatic_pairwise_cutoff
argument. In cases where molecular electrostatics are to be
compared with values calculated in a bulk configuration, an explicit splined Coulomb potential
can be used for molecules instead of the damped shifted force potential. In this case a
CoulombN2Spline
electrostatic solver can be added by using
use_electrostatic_direct_sum
argument. As with the damped shifted force potential the cutoff is
controlled using the electrostatic_pairwise_cutoff
argument, and like the Lennard-Jones
potential is brought to zero with a spline whose length is defined by the
electrostatic_smoothing_length
argument. For all Coulomb solvers the total electrostatic
energies and forces can be scaled up or down using the electrostatic_scale
argument.
Once the OPLSPotentialBuilder
object is instantiated, a calculator is created using
the createCalculator
method. This takes a configuration and returns the appropriate OPLS-AA
calculator for that configuration. To use electrostatic potentials in the calculator, atomic
partial charges have to be assigned to each atom. The OPLS-AA defines appropriate charges for each
atom type, which by default are added to the created calculator. It is also possible to add
pre-calculated charges from different sources using the atomic_charges
argument. This argument
takes a list of charges that have the same ordering as the atoms in the configuration. The partial
atomic charges can also be averaged according to a few different schemes. This is controlled by the
charge_averaging
parameter, which takes an NLFlag which indicates the type of averaging. The
flag AtomConnectivity
averages charges on each atom based on the element and the other elements
it is bonded to. The flag AtomTypes
averages according to the assigned OPLS-AA type. For the
calculator to recognize appropriate inversion terms, improper torsion indices also need to be set on
the configuration. By default this is done automatically, however existing improper torsions can be
preserved by setting the set_improper_torsion_indices
to False
. As the OPLS-AA potential
contains a number of different potentials, details of each potential used to create the calculator
can be printed out using the verbose
argument, which is False
by default.
When creating a OPLS-AA calculator for a configuration, the configuration must have the appropriate
tags for each atom type. Atom types in the OPLS-AA forcefield are grouped into classes. Atom types
in a class all have the same valence potential terms and individual Lennard-Jones and electrostatic
potentials. This reduces the number of parameters required to construct a transferable forcefield
while still allowing specific tailoring of the potentials to different molecules. Tags used by the
OPLS-AA forcefield all begin with the prefix OPLS_
. The next section of tag gives the label of
the class of the atomic type, and the final section gives the number of the type in its class.
Appropriate OPLS tags can be assigned automatically using the method assignAtomTypes
. This
method uses the bonding environment of each atom to estimate appropriate types. These tags are
then added to the input configuration. Type definitions have been added for atom types common to
polymers, although not all atom types have a definition. In cases where the
OPLSPotentialBuilder
is unable to recognize the type of an atom, the tag OPLS_UNKNOWN
is added to that atom.
The OPLSPotentialBuilder contains a few methods to assist in manually assigning the correct types
to configurations in cases where the types cannot be automatically assigned. Atom types in the
OPLS-AA potential are generally defined for specific functional groups. As a result, atomic partial
charges for a given functional group usually add to zero or an integer charge for charged groups.
To discover the appropriate types for a configuration, the method printTotalTypes
prints a
table of the available atomic types, with a short description of the kinds of atoms this type
represents. As they are defined together, types for different functional groups are often grouped
together. Types for just one class can also be printed using the type_class
argument. The types
assigned on a specific configuration can also be displayed using the printUsedTypes
method. As
functional groups generally have zero overall charge, correctly typed configurations should also
have zero overall charge. Calculating the total assigned charge to the configuration is therefore a
useful check to make sure the correct types have been added to the configuration. The total assigned
charge to a configuration using the OPLS forcefield can be calculated using the totalCharge
method.