Modifying a SMIRNOFF force field

In this example, we’ll parameterize a ligand automatically, and then play with its parameters to demonstrate the Toolkit’s ability to facilitate force field optimization. For each modification, we’ll calculate the energy for the original conformation, and then minimize the energy and visualize the result.

from copy import deepcopy

from openff.units import unit
from openff.units.openmm import from_openmm, to_openmm
from openmm import LangevinIntegrator, app
from openmm import unit as openmm_unit

import openff.toolkit.typing.engines.smirnoff.parameters as offtk_parameters
from openff.toolkit import ForceField, Molecule

We’re going to do a lot of changing a parameter and then visualising what happened, so let’s define a convenience function to do just that

def minimize_and_visualize(molecule, forcefield):
    # Sort out our input data
    mol_topology = molecule.to_topology()
    mol_system = forcefield.create_openmm_system(
        mol_topology,
        charge_from_molecules=[molecule],
    )

    # Set up the minimization and point calculation
    integrator = LangevinIntegrator(
        300 * openmm_unit.kelvin,
        1 / openmm_unit.picosecond,
        0.002 * openmm_unit.picoseconds,
    )
    simulation = app.Simulation(mol_topology.to_openmm(), mol_system, integrator)
    simulation.context.setPositions(to_openmm(molecule.conformers[0]))

    # Get the initial energy
    initial_potential = simulation.context.getState(getEnergy=True).getPotentialEnergy()

    # Energy minimize
    simulation.minimizeEnergy()
    minimized_state = simulation.context.getState(getPositions=True, getEnergy=True)
    minimized_potential = minimized_state.getPotentialEnergy()
    minimized_coords = from_openmm(minimized_state.getPositions(asNumpy=True))

    # Visualize
    vis_mol = deepcopy(molecule)
    vis_mol.conformers[0] = minimized_coords
    view = vis_mol.visualize(backend="nglview")
    print(
        f"Initial energy is {initial_potential.format('%0.1F')};",
        f"Minimized energy is {minimized_potential.format('%0.1F')}",
    )
    return view

Getting to know you — the molecule

This “ligand” is a modified version of the molecule we introduced in the Toolkit Showcase. It’s just been altered for a slightly more exciting example here. This also lets us demonstrate constructing a Molecule from a SMILES string!

ligand_smiles = "CC(C)(C)c1c(O)c(O)c2c(c1O)[C@H]1OCCC[C@H]1[C@H](c1cc(O)c(O)c(F)c1)N2"
ligand = Molecule.from_smiles(ligand_smiles)
ligand.generate_conformers(n_conformers=1)
force_field = ForceField("openff-2.2.0.offxml")

Computing charges is expensive, and we’re going to be changing the force field a lot, so we can save time by computing them just once and caching them.

ℹ️ Note that when we call create_openmm_system above we pass in the charges with the charge_from_molecules argument!
ligand.partial_charges = force_field.get_partial_charges(ligand)

The ligand, visualised

Let’s take a close look at the ligand and decide what we want to modify. We’ll label the atoms with their indexes so we can identify them later

view = minimize_and_visualize(ligand, force_field)
view.add_label(label_type="atomindex", color="black", attachment="middle-center")
view
Initial energy is 39.6 kJ/mol; Minimized energy is -313.3 kJ/mol

Getting to know all about you — Investigating assigned parameters

Let’s start with something simple — lengthening the bond to the fluorine. We can use ForceField.label_molecules to identify the appropriate parameters:

ff_applied_parameters = force_field.label_molecules(ligand.to_topology())
ff_applied_parameters
[{'Constraints': ValenceDict(
    {(0, 30): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (0, 31): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (0, 32): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (2, 33): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (2, 34): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (2, 35): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (3, 36): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (3, 37): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (3, 38): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (6, 39): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (8, 40): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (12, 41): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (13, 42): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (15, 43): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (15, 44): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (16, 45): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (16, 46): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (17, 47): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (17, 48): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (18, 49): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (19, 50): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (21, 51): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (23, 52): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (25, 53): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (28, 54): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >,
     (29, 55): <ConstraintType with smirks: [#1:1]-[*:2]  id: c1  >}),
  'Bonds': ValenceDict(
    {(0,
      1): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (0,
      30): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (0,
      31): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (0,
      32): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (1,
      2): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (1,
      3): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (1,
      4): <BondType with smirks: [#6X4:1]-[#6X3:2]  id: b2  length: 1.50959128588 angstrom  k: 478.7391355181 kilocalorie_per_mole / angstrom ** 2  >,
     (2,
      33): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (2,
      34): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (2,
      35): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (3,
      36): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (3,
      37): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (3,
      38): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (4,
      5): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (4,
      11): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (5,
      6): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
     (5,
      7): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (6,
      39): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
     (7,
      8): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
     (7,
      9): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (8,
      40): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
     (9,
      10): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (9,
      29): <BondType with smirks: [#6X3:1]-[#7X3:2]  id: b8  length: 1.389502949007 angstrom  k: 658.5261427735 kilocalorie_per_mole / angstrom ** 2  >,
     (10,
      11): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (10,
      13): <BondType with smirks: [#6X4:1]-[#6X3:2]  id: b2  length: 1.50959128588 angstrom  k: 478.7391355181 kilocalorie_per_mole / angstrom ** 2  >,
     (11,
      12): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
     (12,
      41): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
     (13,
      14): <BondType with smirks: [#6X4:1]-[#8X2H0:2]  id: b16  length: 1.436323905911 angstrom  k: 454.0742933374 kilocalorie_per_mole / angstrom ** 2  >,
     (13,
      18): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (13,
      42): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (14,
      15): <BondType with smirks: [#6X4:1]-[#8X2H0:2]  id: b16  length: 1.436323905911 angstrom  k: 454.0742933374 kilocalorie_per_mole / angstrom ** 2  >,
     (15,
      16): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (15,
      43): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (15,
      44): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (16,
      17): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (16,
      45): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (16,
      46): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (17,
      18): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (17,
      47): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (17,
      48): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (18,
      19): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
     (18,
      49): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (19,
      20): <BondType with smirks: [#6X4:1]-[#6X3:2]  id: b2  length: 1.50959128588 angstrom  k: 478.7391355181 kilocalorie_per_mole / angstrom ** 2  >,
     (19,
      29): <BondType with smirks: [#6:1]-[#7:2]  id: b7  length: 1.477375101044 angstrom  k: 451.9422814679 kilocalorie_per_mole / angstrom ** 2  >,
     (19,
      50): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
     (20,
      21): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (20,
      28): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (21,
      22): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (21,
      51): <BondType with smirks: [#6X3:1]-[#1:2]  id: b85  length: 1.08603610657 angstrom  k: 772.8798736582 kilocalorie_per_mole / angstrom ** 2  >,
     (22,
      23): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
     (22,
      24): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (23,
      52): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
     (24,
      25): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
     (24,
      26): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (25,
      53): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
     (26,
      27): <BondType with smirks: [#6:1]-[#9:2]  id: b68  length: 1.353861210984 angstrom  k: 710.5975770066 kilocalorie_per_mole / angstrom ** 2  >,
     (26,
      28): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
     (28,
      54): <BondType with smirks: [#6X3:1]-[#1:2]  id: b85  length: 1.08603610657 angstrom  k: 772.8798736582 kilocalorie_per_mole / angstrom ** 2  >,
     (29,
      55): <BondType with smirks: [#7:1]-[#1:2]  id: b87  length: 1.018604730378 angstrom  k: 961.9809758788 kilocalorie_per_mole / angstrom ** 2  >}),
  'Angles': ValenceDict(
    {(0,
      1,
      2): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (0,
      1,
      3): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (0,
      1,
      4): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      0,
      30): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      0,
      31): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      0,
      32): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      2,
      33): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      2,
      34): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      2,
      35): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      3,
      36): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      3,
      37): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      3,
      38): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (1,
      4,
      5): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (1,
      4,
      11): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (2,
      1,
      3): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (2,
      1,
      4): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (3,
      1,
      4): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (4,
      5,
      6): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (4,
      5,
      7): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (4,
      11,
      10): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (4,
      11,
      12): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (5,
      4,
      11): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (5,
      6,
      39): <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >,
     (5,
      7,
      8): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (5,
      7,
      9): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (6,
      5,
      7): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (7,
      8,
      40): <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >,
     (7,
      9,
      10): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (7,
      9,
      29): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (8,
      7,
      9): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (9,
      10,
      11): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (9,
      10,
      13): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (9,
      29,
      19): <AngleType with smirks: [*:1]~[#7X3$(*~[#6X3,#6X2,#7X2+0]):2]~[*:3]  angle: 121.7705343867 degree  k: 158.9566366884 kilocalorie_per_mole / radian ** 2  id: a20  >,
     (9,
      29,
      55): <AngleType with smirks: [#1:1]-[#7X3$(*~[#6X3,#6X2,#7X2+0]):2]-[*:3]  angle: 118.0636611293 degree  k: 74.20968211314 kilocalorie_per_mole / radian ** 2  id: a21  >,
     (10,
      9,
      29): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (10,
      11,
      12): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (10,
      13,
      14): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (10,
      13,
      18): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (10,
      13,
      42): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (11,
      10,
      13): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (11,
      12,
      41): <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >,
     (13,
      14,
      15): <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >,
     (13,
      18,
      17): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (13,
      18,
      19): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (13,
      18,
      49): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (14,
      13,
      18): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (14,
      13,
      42): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (14,
      15,
      16): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (14,
      15,
      43): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (14,
      15,
      44): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (15,
      16,
      17): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (15,
      16,
      45): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (15,
      16,
      46): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (16,
      15,
      43): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (16,
      15,
      44): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (16,
      17,
      18): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (16,
      17,
      47): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (16,
      17,
      48): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (17,
      16,
      45): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (17,
      16,
      46): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (17,
      18,
      19): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (17,
      18,
      49): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (18,
      13,
      42): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (18,
      17,
      47): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (18,
      17,
      48): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (18,
      19,
      20): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (18,
      19,
      29): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (18,
      19,
      50): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (19,
      18,
      49): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (19,
      20,
      21): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (19,
      20,
      28): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (19,
      29,
      55): <AngleType with smirks: [#1:1]-[#7X3$(*~[#6X3,#6X2,#7X2+0]):2]-[*:3]  angle: 118.0636611293 degree  k: 74.20968211314 kilocalorie_per_mole / radian ** 2  id: a21  >,
     (20,
      19,
      29): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (20,
      19,
      50): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (20,
      21,
      22): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (20,
      21,
      51): <AngleType with smirks: [#1:1]-[#6X3:2]~[*:3]  angle: 119.8544859917 degree  k: 65.84475468128 kilocalorie_per_mole / radian ** 2  id: a11  >,
     (20,
      28,
      26): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (20,
      28,
      54): <AngleType with smirks: [#1:1]-[#6X3:2]~[*:3]  angle: 119.8544859917 degree  k: 65.84475468128 kilocalorie_per_mole / radian ** 2  id: a11  >,
     (21,
      20,
      28): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (21,
      22,
      23): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (21,
      22,
      24): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (22,
      21,
      51): <AngleType with smirks: [#1:1]-[#6X3:2]~[*:3]  angle: 119.8544859917 degree  k: 65.84475468128 kilocalorie_per_mole / radian ** 2  id: a11  >,
     (22,
      23,
      52): <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >,
     (22,
      24,
      25): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (22,
      24,
      26): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (23,
      22,
      24): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (24,
      25,
      53): <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >,
     (24,
      26,
      27): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (24,
      26,
      28): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (25,
      24,
      26): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (26,
      28,
      54): <AngleType with smirks: [#1:1]-[#6X3:2]~[*:3]  angle: 119.8544859917 degree  k: 65.84475468128 kilocalorie_per_mole / radian ** 2  id: a11  >,
     (27,
      26,
      28): <AngleType with smirks: [*:1]~[#6X3:2]~[*:3]  angle: 120.1049996386 degree  k: 169.0953155969 kilocalorie_per_mole / radian ** 2  id: a10  >,
     (29,
      19,
      50): <AngleType with smirks: [*:1]~[#6X4:2]-[*:3]  angle: 109.9380527584 degree  k: 133.5126366319 kilocalorie_per_mole / radian ** 2  id: a1  >,
     (30,
      0,
      31): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (30,
      0,
      32): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (31,
      0,
      32): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (33,
      2,
      34): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (33,
      2,
      35): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (34,
      2,
      35): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (36,
      3,
      37): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (36,
      3,
      38): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (37,
      3,
      38): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (43,
      15,
      44): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (45,
      16,
      46): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >,
     (47,
      17,
      48): <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >}),
  'ProperTorsions': ValenceDict(
    {(0,
      1,
      2,
      33): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      2,
      34): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      2,
      35): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      3,
      36): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      3,
      37): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      3,
      38): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      4,
      5): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (0,
      1,
      4,
      11): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (1,
      4,
      5,
      6): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (1,
      4,
      5,
      7): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (1,
      4,
      11,
      10): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (1,
      4,
      11,
      12): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      0,
      30): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      0,
      31): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      0,
      32): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      3,
      36): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      3,
      37): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      3,
      38): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      4,
      5): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (2,
      1,
      4,
      11): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      0,
      30): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      0,
      31): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      0,
      32): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      2,
      33): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      2,
      34): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      2,
      35): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      4,
      5): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (3,
      1,
      4,
      11): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      0,
      30): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      0,
      31): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      0,
      32): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      2,
      33): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      2,
      34): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      2,
      35): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      3,
      36): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      3,
      37): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      1,
      3,
      38): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      5,
      6,
      39): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      5,
      7,
      8): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      5,
      7,
      9): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      11,
      10,
      9): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      11,
      10,
      13): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (4,
      11,
      12,
      41): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (5,
      4,
      11,
      10): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (5,
      4,
      11,
      12): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (5,
      7,
      8,
      40): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (5,
      7,
      9,
      10): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (5,
      7,
      9,
      29): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (6,
      5,
      4,
      11): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (6,
      5,
      7,
      8): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (6,
      5,
      7,
      9): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (7,
      5,
      4,
      11): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (7,
      5,
      6,
      39): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (7,
      9,
      10,
      11): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (7,
      9,
      10,
      13): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (7,
      9,
      29,
      19): <ProperTorsionType with smirks: [*:1]~[#7X3,#7X2-1:2]-[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t73  k1: 0.7020804385738 kilocalorie / mole  idivf1: 1.0  >,
     (7,
      9,
      29,
      55): <ProperTorsionType with smirks: [*:1]~[#7X3,#7X2-1:2]-[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t73  k1: 0.7020804385738 kilocalorie / mole  idivf1: 1.0  >,
     (8,
      7,
      9,
      10): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (8,
      7,
      9,
      29): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (9,
      7,
      8,
      40): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (9,
      10,
      11,
      12): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (9,
      10,
      13,
      14): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (9,
      10,
      13,
      18): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (9,
      10,
      13,
      42): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (9,
      29,
      19,
      18): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#7X3$(*~[#6X3,#6X2]):3]~[*:4]  periodicity1: 2  periodicity2: 3  phase1: 0.0 degree  phase2: 0.0 degree  id: t64  k1: 0.3215160932459 kilocalorie / mole  k2: 0.2510071310501 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (9,
      29,
      19,
      20): <ProperTorsionType with smirks: [#6X3:1]-[#7X3:2]-[#6X4:3]-[#6X3:4]  periodicity1: 2  periodicity2: 1  phase1: 180.0 degree  phase2: 0.0 degree  id: t66  k1: -0.5801351028897 kilocalorie / mole  k2: -0.5728076214143 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (9,
      29,
      19,
      50): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#7X3$(*~[#6X3,#6X2]):3]~[*:4]  periodicity1: 2  periodicity2: 3  phase1: 0.0 degree  phase2: 0.0 degree  id: t64  k1: 0.3215160932459 kilocalorie / mole  k2: 0.2510071310501 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (10,
      9,
      29,
      19): <ProperTorsionType with smirks: [*:1]~[#7X3,#7X2-1:2]-[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t73  k1: 0.7020804385738 kilocalorie / mole  idivf1: 1.0  >,
     (10,
      9,
      29,
      55): <ProperTorsionType with smirks: [*:1]~[#7X3,#7X2-1:2]-[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t73  k1: 0.7020804385738 kilocalorie / mole  idivf1: 1.0  >,
     (10,
      11,
      12,
      41): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (10,
      13,
      14,
      15): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#8X2H0:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t95  k1: 0.7893340970024 kilocalorie / mole  idivf1: 3.0  >,
     (10,
      13,
      18,
      17): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (10,
      13,
      18,
      19): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (10,
      13,
      18,
      49): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (11,
      10,
      9,
      29): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (11,
      10,
      13,
      14): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (11,
      10,
      13,
      18): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (11,
      10,
      13,
      42): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (12,
      11,
      10,
      13): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (13,
      10,
      9,
      29): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (13,
      14,
      15,
      16): <ProperTorsionType with smirks: [#6X4:1]-[#6X4:2]-[#8X2H0:3]-[#6X4:4]  periodicity1: 3  periodicity2: 2  phase1: 0.0 degree  phase2: 180.0 degree  id: t96  k1: 0.3102287948961 kilocalorie / mole  k2: -0.04857068407726 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (13,
      14,
      15,
      43): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#8X2H0:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t95  k1: 0.7893340970024 kilocalorie / mole  idivf1: 3.0  >,
     (13,
      14,
      15,
      44): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#8X2H0:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t95  k1: 0.7893340970024 kilocalorie / mole  idivf1: 3.0  >,
     (13,
      18,
      17,
      16): <ProperTorsionType with smirks: [#6X4:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  periodicity2: 2  periodicity3: 1  phase1: 0.0 degree  phase2: 180.0 degree  phase3: 180.0 degree  id: t2  k1: 0.3549667638131 kilocalorie / mole  k2: 0.2399307824079 kilocalorie / mole  k3: 0.9156101580834 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  idivf3: 1.0  >,
     (13,
      18,
      17,
      47): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (13,
      18,
      17,
      48): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (13,
      18,
      19,
      20): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (13,
      18,
      19,
      29): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (13,
      18,
      19,
      50): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (14,
      13,
      18,
      17): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (14,
      13,
      18,
      19): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (14,
      13,
      18,
      49): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#8X2:4]  periodicity1: 3  periodicity2: 1  phase1: 0.0 degree  phase2: 0.0 degree  id: t9  k1: 0.1311347012225 kilocalorie / mole  k2: 0.436546671918 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (14,
      15,
      16,
      17): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (14,
      15,
      16,
      45): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#8X2:4]  periodicity1: 3  periodicity2: 1  phase1: 0.0 degree  phase2: 0.0 degree  id: t9  k1: 0.1311347012225 kilocalorie / mole  k2: 0.436546671918 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (14,
      15,
      16,
      46): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#8X2:4]  periodicity1: 3  periodicity2: 1  phase1: 0.0 degree  phase2: 0.0 degree  id: t9  k1: 0.1311347012225 kilocalorie / mole  k2: 0.436546671918 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (15,
      14,
      13,
      18): <ProperTorsionType with smirks: [#6X4:1]-[#6X4:2]-[#8X2H0:3]-[#6X4:4]  periodicity1: 3  periodicity2: 2  phase1: 0.0 degree  phase2: 180.0 degree  id: t96  k1: 0.3102287948961 kilocalorie / mole  k2: -0.04857068407726 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (15,
      14,
      13,
      42): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#8X2H0:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t95  k1: 0.7893340970024 kilocalorie / mole  idivf1: 3.0  >,
     (15,
      16,
      17,
      18): <ProperTorsionType with smirks: [#6X4:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  periodicity2: 2  periodicity3: 1  phase1: 0.0 degree  phase2: 180.0 degree  phase3: 180.0 degree  id: t2  k1: 0.3549667638131 kilocalorie / mole  k2: 0.2399307824079 kilocalorie / mole  k3: 0.9156101580834 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  idivf3: 1.0  >,
     (15,
      16,
      17,
      47): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (15,
      16,
      17,
      48): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (16,
      17,
      18,
      19): <ProperTorsionType with smirks: [#6X4:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  periodicity2: 2  periodicity3: 1  phase1: 0.0 degree  phase2: 180.0 degree  phase3: 180.0 degree  id: t2  k1: 0.3549667638131 kilocalorie / mole  k2: 0.2399307824079 kilocalorie / mole  k3: 0.9156101580834 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  idivf3: 1.0  >,
     (16,
      17,
      18,
      49): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (17,
      16,
      15,
      43): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (17,
      16,
      15,
      44): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (17,
      18,
      13,
      42): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (17,
      18,
      19,
      20): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (17,
      18,
      19,
      29): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (17,
      18,
      19,
      50): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (18,
      17,
      16,
      45): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (18,
      17,
      16,
      46): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (18,
      19,
      20,
      21): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (18,
      19,
      20,
      28): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (18,
      19,
      29,
      55): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#7X3$(*~[#6X3,#6X2]):3]~[*:4]  periodicity1: 2  periodicity2: 3  phase1: 0.0 degree  phase2: 0.0 degree  id: t64  k1: 0.3215160932459 kilocalorie / mole  k2: 0.2510071310501 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (19,
      18,
      13,
      42): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (19,
      18,
      17,
      47): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (19,
      18,
      17,
      48): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#6X4:4]  periodicity1: 3  phase1: 0.0 degree  id: t4  k1: 0.0985377828486 kilocalorie / mole  idivf1: 1.0  >,
     (19,
      20,
      21,
      22): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (19,
      20,
      21,
      51): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (19,
      20,
      28,
      26): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (19,
      20,
      28,
      54): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (20,
      19,
      18,
      49): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (20,
      19,
      29,
      55): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#7X3$(*~[#6X3,#6X2]):3]~[*:4]  periodicity1: 2  periodicity2: 3  phase1: 0.0 degree  phase2: 0.0 degree  id: t64  k1: 0.3215160932459 kilocalorie / mole  k2: 0.2510071310501 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >,
     (20,
      21,
      22,
      23): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (20,
      21,
      22,
      24): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (20,
      28,
      26,
      24): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (20,
      28,
      26,
      27): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      20,
      19,
      29): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      20,
      19,
      50): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      20,
      28,
      26): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      20,
      28,
      54): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      22,
      23,
      52): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      22,
      24,
      25): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (21,
      22,
      24,
      26): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (22,
      21,
      20,
      28): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (22,
      24,
      25,
      53): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (22,
      24,
      26,
      27): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (22,
      24,
      26,
      28): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (23,
      22,
      21,
      51): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (23,
      22,
      24,
      25): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (23,
      22,
      24,
      26): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (24,
      22,
      21,
      51): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (24,
      22,
      23,
      52): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (24,
      26,
      28,
      54): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (25,
      24,
      26,
      27): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (25,
      24,
      26,
      28): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (26,
      24,
      25,
      53): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#8X2:3]-[#1:4]  periodicity1: 2  phase1: 180.0 degree  id: t106  k1: 0.9816061886151 kilocalorie / mole  idivf1: 1.0  >,
     (27,
      26,
      28,
      54): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (28,
      20,
      19,
      29): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (28,
      20,
      19,
      50): <ProperTorsionType with smirks: [*:1]~[#6X3:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t17  k1: 0.1898990064757 kilocalorie / mole  idivf1: 1.0  >,
     (28,
      20,
      21,
      51): <ProperTorsionType with smirks: [*:1]~[#6X3:2]:[#6X3:3]~[*:4]  periodicity1: 2  phase1: 180.0 degree  id: t44  k1: 3.268474846449 kilocalorie / mole  idivf1: 1.0  >,
     (29,
      19,
      18,
      49): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#6X4:3]-[*:4]  periodicity1: 3  phase1: 0.0 degree  id: t1  k1: 0.143904881748 kilocalorie / mole  idivf1: 1.0  >,
     (42,
      13,
      18,
      49): <ProperTorsionType with smirks: [*:1]~[#6x3:2](~[#7,#8,#16])-[#6x3:3]~[*:4]  periodicity1: 1  phase1: 0.0 degree  id: t141c  k1: -3.525605054758 kilocalorie / mole  idivf1: 1.0  >,
     (43,
      15,
      16,
      45): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (43,
      15,
      16,
      46): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (44,
      15,
      16,
      45): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (44,
      15,
      16,
      46): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (45,
      16,
      17,
      47): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (45,
      16,
      17,
      48): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (46,
      16,
      17,
      47): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (46,
      16,
      17,
      48): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (47,
      17,
      18,
      49): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (48,
      17,
      18,
      49): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (49,
      18,
      19,
      50): <ProperTorsionType with smirks: [#1:1]-[#6X4:2]-[#6X4:3]-[#1:4]  periodicity1: 3  phase1: 0.0 degree  id: t3  k1: 0.2279946533057 kilocalorie / mole  idivf1: 1.0  >,
     (50,
      19,
      29,
      55): <ProperTorsionType with smirks: [*:1]-[#6X4:2]-[#7X3$(*~[#6X3,#6X2]):3]~[*:4]  periodicity1: 2  periodicity2: 3  phase1: 0.0 degree  phase2: 0.0 degree  id: t64  k1: 0.3215160932459 kilocalorie / mole  k2: 0.2510071310501 kilocalorie / mole  idivf1: 1.0  idivf2: 1.0  >}),
  'ImproperTorsions': <openff.toolkit.topology.topology.ImproperDict at 0x7f20b02abee0>,
  'vdW': ValenceDict(
    {(0,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (1,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (2,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (3,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (4,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (5,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (6,): <vdWType with smirks: [#8X2H1+0:1]  epsilon: 0.2094735324129 kilocalorie / mole  id: n19  rmin_half: 1.682099169199 angstrom  >,
     (7,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (8,): <vdWType with smirks: [#8X2H1+0:1]  epsilon: 0.2094735324129 kilocalorie / mole  id: n19  rmin_half: 1.682099169199 angstrom  >,
     (9,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (10,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (11,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (12,): <vdWType with smirks: [#8X2H1+0:1]  epsilon: 0.2094735324129 kilocalorie / mole  id: n19  rmin_half: 1.682099169199 angstrom  >,
     (13,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (14,): <vdWType with smirks: [#8X2H0+0:1]  epsilon: 0.1684651402602 kilocalorie / mole  id: n18  rmin_half: 1.697783613804 angstrom  >,
     (15,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (16,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (17,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (18,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (19,): <vdWType with smirks: [#6X4:1]  epsilon: 0.1088406109251 kilocalorie / mole  id: n16  rmin_half: 1.896698071741 angstrom  >,
     (20,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (21,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (22,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (23,): <vdWType with smirks: [#8X2H1+0:1]  epsilon: 0.2094735324129 kilocalorie / mole  id: n19  rmin_half: 1.682099169199 angstrom  >,
     (24,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (25,): <vdWType with smirks: [#8X2H1+0:1]  epsilon: 0.2094735324129 kilocalorie / mole  id: n19  rmin_half: 1.682099169199 angstrom  >,
     (26,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (27,): <vdWType with smirks: [#9:1]  epsilon: 0.061 kilocalorie / mole  id: n23  rmin_half: 1.75 angstrom  >,
     (28,): <vdWType with smirks: [#6:1]  epsilon: 0.0868793154488 kilocalorie / mole  id: n14  rmin_half: 1.953447017081 angstrom  >,
     (29,): <vdWType with smirks: [#7:1]  epsilon: 0.1676915150424 kilocalorie / mole  id: n20  rmin_half: 1.799798315098 angstrom  >,
     (30,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (31,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (32,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (33,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (34,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (35,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (36,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (37,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (38,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (39,): <vdWType with smirks: [#1:1]-[#8]  epsilon: 1.232599966667e-05 kilocalorie / mole  id: n12  rmin_half: 0.2999999999997 angstrom  >,
     (40,): <vdWType with smirks: [#1:1]-[#8]  epsilon: 1.232599966667e-05 kilocalorie / mole  id: n12  rmin_half: 0.2999999999997 angstrom  >,
     (41,): <vdWType with smirks: [#1:1]-[#8]  epsilon: 1.232599966667e-05 kilocalorie / mole  id: n12  rmin_half: 0.2999999999997 angstrom  >,
     (42,): <vdWType with smirks: [#1:1]-[#6X4]-[#7,#8,#9,#16,#17,#35]  epsilon: 0.01640924602775 kilocalorie / mole  id: n3  rmin_half: 1.449786411317 angstrom  >,
     (43,): <vdWType with smirks: [#1:1]-[#6X4]-[#7,#8,#9,#16,#17,#35]  epsilon: 0.01640924602775 kilocalorie / mole  id: n3  rmin_half: 1.449786411317 angstrom  >,
     (44,): <vdWType with smirks: [#1:1]-[#6X4]-[#7,#8,#9,#16,#17,#35]  epsilon: 0.01640924602775 kilocalorie / mole  id: n3  rmin_half: 1.449786411317 angstrom  >,
     (45,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (46,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (47,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (48,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (49,): <vdWType with smirks: [#1:1]-[#6X4]  epsilon: 0.01577948280971 kilocalorie / mole  id: n2  rmin_half: 1.48419980825 angstrom  >,
     (50,): <vdWType with smirks: [#1:1]-[#6X4]-[#7,#8,#9,#16,#17,#35]  epsilon: 0.01640924602775 kilocalorie / mole  id: n3  rmin_half: 1.449786411317 angstrom  >,
     (51,): <vdWType with smirks: [#1:1]-[#6X3]  epsilon: 0.01561134320353 kilocalorie / mole  id: n7  rmin_half: 1.443812569645 angstrom  >,
     (52,): <vdWType with smirks: [#1:1]-[#8]  epsilon: 1.232599966667e-05 kilocalorie / mole  id: n12  rmin_half: 0.2999999999997 angstrom  >,
     (53,): <vdWType with smirks: [#1:1]-[#8]  epsilon: 1.232599966667e-05 kilocalorie / mole  id: n12  rmin_half: 0.2999999999997 angstrom  >,
     (54,): <vdWType with smirks: [#1:1]-[#6X3]  epsilon: 0.01561134320353 kilocalorie / mole  id: n7  rmin_half: 1.443812569645 angstrom  >,
     (55,): <vdWType with smirks: [#1:1]-[#7]  epsilon: 0.01409081474669 kilocalorie / mole  id: n11  rmin_half: 0.6192778454102 angstrom  >}),
  'Electrostatics': ValenceDict({}),
  'LibraryCharges': {},
  'ToolkitAM1BCC': ValenceDict({})}]

label_molecules() returns a list of the molecules that can be parameterized from the arguments. It returns them as dictionaries that give us access to the ForceField parameters that are used for the molecule. Since we’re only passing in a single molecule, ff_applied_parameters is a list of one element. We can see the bonds that are used for the ligand by converting the appropriate ValenceDict to a regular dictionary:

dict(ff_applied_parameters[0]["Bonds"])
{(0,
  1): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (0,
  30): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (0,
  31): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (0,
  32): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (1,
  2): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (1,
  3): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (1,
  4): <BondType with smirks: [#6X4:1]-[#6X3:2]  id: b2  length: 1.50959128588 angstrom  k: 478.7391355181 kilocalorie_per_mole / angstrom ** 2  >,
 (2,
  33): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (2,
  34): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (2,
  35): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (3,
  36): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (3,
  37): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (3,
  38): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (4,
  5): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (4,
  11): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (5,
  6): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
 (5,
  7): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (6,
  39): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
 (7,
  8): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
 (7,
  9): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (8,
  40): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
 (9,
  10): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (9,
  29): <BondType with smirks: [#6X3:1]-[#7X3:2]  id: b8  length: 1.389502949007 angstrom  k: 658.5261427735 kilocalorie_per_mole / angstrom ** 2  >,
 (10,
  11): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (10,
  13): <BondType with smirks: [#6X4:1]-[#6X3:2]  id: b2  length: 1.50959128588 angstrom  k: 478.7391355181 kilocalorie_per_mole / angstrom ** 2  >,
 (11,
  12): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
 (12,
  41): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
 (13,
  14): <BondType with smirks: [#6X4:1]-[#8X2H0:2]  id: b16  length: 1.436323905911 angstrom  k: 454.0742933374 kilocalorie_per_mole / angstrom ** 2  >,
 (13,
  18): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (13,
  42): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (14,
  15): <BondType with smirks: [#6X4:1]-[#8X2H0:2]  id: b16  length: 1.436323905911 angstrom  k: 454.0742933374 kilocalorie_per_mole / angstrom ** 2  >,
 (15,
  16): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (15,
  43): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (15,
  44): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (16,
  17): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (16,
  45): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (16,
  46): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (17,
  18): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (17,
  47): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (17,
  48): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (18,
  19): <BondType with smirks: [#6X4:1]-[#6X4:2]  id: b1  length: 1.533627603844 angstrom  k: 430.6027811279 kilocalorie_per_mole / angstrom ** 2  >,
 (18,
  49): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (19,
  20): <BondType with smirks: [#6X4:1]-[#6X3:2]  id: b2  length: 1.50959128588 angstrom  k: 478.7391355181 kilocalorie_per_mole / angstrom ** 2  >,
 (19,
  29): <BondType with smirks: [#6:1]-[#7:2]  id: b7  length: 1.477375101044 angstrom  k: 451.9422814679 kilocalorie_per_mole / angstrom ** 2  >,
 (19,
  50): <BondType with smirks: [#6X4:1]-[#1:2]  id: b84  length: 1.094206460752 angstrom  k: 715.6382948318 kilocalorie_per_mole / angstrom ** 2  >,
 (20,
  21): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (20,
  28): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (21,
  22): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (21,
  51): <BondType with smirks: [#6X3:1]-[#1:2]  id: b85  length: 1.08603610657 angstrom  k: 772.8798736582 kilocalorie_per_mole / angstrom ** 2  >,
 (22,
  23): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
 (22,
  24): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (23,
  52): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
 (24,
  25): <BondType with smirks: [#6X3:1]-[#8X2H1:2]  id: b18  length: 1.360361459598 angstrom  k: 692.682789202 kilocalorie_per_mole / angstrom ** 2  >,
 (24,
  26): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (25,
  53): <BondType with smirks: [#8:1]-[#1:2]  id: b88  length: 0.9752301008088 angstrom  k: 1076.70158803 kilocalorie_per_mole / angstrom ** 2  >,
 (26,
  27): <BondType with smirks: [#6:1]-[#9:2]  id: b68  length: 1.353861210984 angstrom  k: 710.5975770066 kilocalorie_per_mole / angstrom ** 2  >,
 (26,
  28): <BondType with smirks: [#6X3:1]:[#6X3:2]  id: b5  length: 1.400706064009 angstrom  k: 753.7514414653 kilocalorie_per_mole / angstrom ** 2  >,
 (28,
  54): <BondType with smirks: [#6X3:1]-[#1:2]  id: b85  length: 1.08603610657 angstrom  k: 772.8798736582 kilocalorie_per_mole / angstrom ** 2  >,
 (29,
  55): <BondType with smirks: [#7:1]-[#1:2]  id: b87  length: 1.018604730378 angstrom  k: 961.9809758788 kilocalorie_per_mole / angstrom ** 2  >}

Take a look at the first entry as an example. The dictionary is keyed by the atomic indices of the particles in the molecule that the parameter applies to, and the values are special types that specify the parameters from the force field. We’ll look in more detail at one of these parameters in a second.

A little goes a long way — Changing parameters

We already decided we wanted to adjust the torsion parameter to the hydroxyl group, so let’s do that. We know parameters are indexed by the atoms they’ve been applied to, so we can look at our labelled widget above and pull out exactly the parameter we care about:

fluorine = next(atom for atom in ligand.atoms if atom.symbol == "F")
carbon = next(iter(fluorine.bonded_atoms))
fluorine_carbon_bond_indices = (
    fluorine.molecule_atom_index,
    carbon.molecule_atom_index,
)
print("The fluorine-carbon bond is between atoms", fluorine_carbon_bond_indices)

ff_applied_parameters[0]["Bonds"][fluorine_carbon_bond_indices]
The fluorine-carbon bond is between atoms (27, 26)
<BondType with smirks: [#6:1]-[#9:2]  id: b68  length: 1.353861210984 angstrom  k: 710.5975770066 kilocalorie_per_mole / angstrom ** 2  >
📗 We have to get the indices of these atoms out programmatically because our SMILES code doesn't specify atom indices, and so you might have different indices to us. Take a look at the widget above and try manually specifying atom indices to get the parameters of a bond!

Let’s dig into this type a bit more. It has a few attributes in its textual representation. The first of these is maybe the most important: the smirks attribute, which tells the Toolkit which atoms this parameter applies to. SMIRKS is a chemical pattern matching format; think of it as the result of a regular expression having a baby with a SMILES string. This one is very simple: it is just a Carbon atom ([#6:1], atomic number 6) singly bonded (-) to a Fluorine atom ([#9:2], atomic number 9). The numbers after the colons just label the atoms — this is a bond, so it has two atoms labelled. This is helpful for when we want to match against atoms that aren’t a part of the parameter; we just don’t label the additional atoms. We’ll use this trick later.

The other important attributes provide the actual parameterization, and so are different for different kinds of parameters. For proper torsions, this is the periodicity of the sinewave that describes the torsion, as well as its phase and amplitude (or force constant, \(k\)). These parameter values are similar to the equivalent description in most other force field formats.

ℹ️ The SMIRKS specification is available online:

https://www.daylight.com/dayhtml/doc/theory/theory.smirks.html

It is closely related to the SMARTS molecular pattern matching language, whose specification is probably more useful for working with the Toolkit and is also available online:

https://www.daylight.com/dayhtml/doc/theory/theory.smarts.html

Modifying a parameter

Unfortunately we can’t just modify this parameter and see the results reflected in the parameterization. We need to get the appropriate parameter from the force field and modify it there.

fluorine_bond = force_field["Bonds"].parameters["[#6:1]-[#9:2]"]
fluorine_bond.length = 10 * unit.angstrom

Here, we’ve selected the proper torsion with the SMIRKS code we found earlier, and changed its force constant by an order of magnitude in the opposite direction! Let’s see what we’ve wrought:

minimize_and_visualize(ligand, force_field)
Initial energy is 111123.7 kJ/mol; Minimized energy is -298.8 kJ/mol

Turns out Pinocchio was a real molecule all along!

Parameters affecting multiple atom groups

Ok, that was fun, but it’s only one parameter; we could easily have made this change in the OpenMM System or a GROMACS ITP file or whatever. What’s the toolkit really giving us here?

Let’s mess with all the H-X-H angles. And let’s not get into SMIRKS this time, lets let the toolkit do the thinking:

ff_applied_parameters = force_field.label_molecules(ligand.to_topology())
for atoms, parameter in ff_applied_parameters[0]["Angles"].items():
    ele_1 = ligand.atoms[atoms[0]].symbol
    ele_2 = ligand.atoms[atoms[1]].symbol
    ele_3 = ligand.atoms[atoms[2]].symbol
    if ele_1 == "H" and ele_3 == "H":
        print(atoms, parameter)
(30, 0, 31) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(30, 0, 32) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(31, 0, 32) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(33, 2, 34) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(33, 2, 35) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(34, 2, 35) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(36, 3, 37) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(36, 3, 38) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(37, 3, 38) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(43, 15, 44) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(45, 16, 46) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >
(47, 17, 48) <AngleType with smirks: [#1:1]-[#6X4:2]-[#1:3]  angle: 108.3543135269 degree  k: 73.55327570151 kilocalorie_per_mole / radian ** 2  id: a2  >

Wow, there’s just one smirks parameter describing all of the H-X-H angles!

ℹ️ It's worth mentioning here that a force field might include multiple parameters with SMIRKS that match a particular group of atoms. When this happens, only the last parameter specified is applied. This allows force field authors to define general parameters first, and then override them with more specific parameters.
📗 Convince yourself that this SMIRKS code matches all of the H-X-H angles in our molecule. What H-X-H angles would it not match?
hxh_angle = force_field["Angles"].parameters["[#1:1]-[#6X4:2]-[#1:3]"]
hxh_angle.angle = 50 * unit.degree

minimize_and_visualize(ligand, force_field)
Initial energy is 113040.5 kJ/mol; Minimized energy is -121.8 kJ/mol

That is not how a methyl group is supposed to look. Guess we did it right!

Introducing a new parameter

Now let’s get really crazy. I’ve always thought heteroatoms in rings looked a bit too comfortable, haven’t you?

Since we have two cyclic heteroatoms in our molecule with very different chemistries, there’ll probably be two seperate parameters we need to modify. We can just pass the SMIRKS that come up for C-X-C angles to code that modifies the parameters directly, no need to manually copy it over.

ff_applied_parameters = force_field.label_molecules(ligand.to_topology())
for atoms, parameter in ff_applied_parameters[0]["Angles"].items():
    ele_1 = ligand.atoms[atoms[0]].symbol
    ele_2 = ligand.atoms[atoms[1]].symbol
    ele_3 = ligand.atoms[atoms[2]].symbol
    if (ele_2 not in ["C", "H"]) and ele_1 == "C" and ele_3 == "C":
        print(atoms, parameter)

        smirks = parameter.smirks
        heteroatom_angle = force_field["Angles"].parameters[smirks]
        # dihedral angles are undefined when linear, so let's not set to 180
        heteroatom_angle.angle = 179 * unit.degree

        # We'll need this later
        if ele_2 == "N":
            cyclic_nitrogen_smirks = smirks

minimize_and_visualize(ligand, force_field)
(9, 29, 19) <AngleType with smirks: [*:1]~[#7X3$(*~[#6X3,#6X2,#7X2+0]):2]~[*:3]  angle: 121.7705343867 degree  k: 158.9566366884 kilocalorie_per_mole / radian ** 2  id: a20  >
(13, 14, 15) <AngleType with smirks: [*:1]-[#8:2]-[*:3]  angle: 112.8736161145 degree  k: 239.791779107 kilocalorie_per_mole / radian ** 2  id: a28  >
Initial energy is 116332.9 kJ/mol; Minimized energy is 178.0 kJ/mol

Oops! We didn’t mean to make all the hydroxyl groups linear! One of our parameters must be applied to Oxygen both when its in a ring and in a hydroxyl group. Sure enough, the second SMIRKS code printed above ([*:1]-[#8:2]-[*:3]) matches any X-O-X angle! We’ll have to define a new, more specific SMIRKS entry so we only capture the appropriate oxygen.

The SMIRKS for our new parameter

Our modified force field is getting a bit messy, so let’s start again. We still want the nitrogen angle from before, so lets save that too

force_field = ForceField("openff-2.1.0.offxml")
cyclic_nitrogen_angle = force_field["Angles"].parameters[cyclic_nitrogen_smirks]
cyclic_nitrogen_angle.angle = 179 * unit.degree

Now we need to define a new SMIRKS entry, specific to oxygens in a ring. We can be as specific as we like, as our new parameter will be at the end of the force field and so will override any conflicting parameters. Let’s define our smirks to match any singly-bonded C-O-C angle in a six-membered pyran ring.

First, we’ll write down the SMILES for a pyran:

C1CCOCC1

The first carbon has a digit following it, C1, which allows us to close the ring later. Then we have two more carbons, CC, an oxygen, O, and the final two carbons, CC. Finally, we close the ring by repeating the same digit from before, 1. Lets check that this produces the chemical we expect before we move on:

Molecule.from_smiles("C1CCOCC1").visualize()
../../../../_images/5a3591a33c17cbec32b1889a969299944f79bd68aaa95e63fa7770ed08815950.svg

Perfect. We can ignore the hydrogens, as we haven’t specified any bond orders - our smirks will match the saturated structure or any of its derivatives that maintain the ring. Now all we have to do is label the atoms we want to extract the bond from! This way, all the atoms must be present for the SMIRKS to match, but only these three atoms define the angle for the parameter. This lets us be more specific than we could be if we could only describe the atoms in the angle itself. We do that the same way we did before, put them in square brackets and append a colon and a number:

cyclic_oxygen_smirks = "C1C[C:1]-[O:2]-[C:3]C1"
📗 Try replacing this SMIRKS code with one of your own. Hint: For this force field and ligand, there's a much simpler code that will uniquely specify the C-O-C bond!

Defining and registering the new parameter

Now, all that’s left is to define the angle parameter, add it to the force field, and see what we’ve wrought!

# Define the new angle parameter
angle_parameter = offtk_parameters.AngleHandler.AngleType(
    smirks=cyclic_oxygen_smirks,
    angle=179 * unit.degree,
    k=134.5019777341 * unit.kilocalorie / (unit.mole * unit.radians**2),
)

# Add the parameter to the force field
angles_handler = force_field["Angles"]
angles_handler.add_parameter(parameter=angle_parameter)

# Visualize the newly parameterized molecule
minimize_and_visualize(ligand, force_field)
Initial energy is 700.4 kJ/mol; Minimized energy is -89.5 kJ/mol

Perfect! The two cyclic heteroatoms are nearly linear while the hydroxyl group retains its characteristic angle.

Bonus: Adding a cosmetic attribute to the new parameter

In some cases, you might want to tag a particular attribute with extra data. This might be useful, for example, in communicating to a fitting tool which parameters should be fitted, or simply to add extra metadata… The toolkit allows for cosmetic attributes, which can be stored to individual parameters, will be included out when the ForceField object is written to disk, but will not otherwise affect the parameterization machinery. The API point for this is ParameterType.add_cosmetic_attribute.

Let’s use this to tag our modified and new parameters with a note about why we introduced it.

angles_handler.parameters[cyclic_nitrogen_smirks].add_cosmetic_attribute(
    attr_name="note",
    attr_value="funky pseudo-linear cyclic nitrogen",
)

angles_handler.parameters[cyclic_oxygen_smirks].add_cosmetic_attribute(
    attr_name="note",
    attr_value="funky pseudo-linear cyclic oxygen",
)

We we look at these parameter now, they include this attribute tacked on at the end:

angles_handler.parameters[cyclic_nitrogen_smirks]
<AngleType with smirks: [*:1]~[#7X3$(*~[#6X3,#6X2,#7X2+0]):2]~[*:3]  angle: 179 degree  k: 151.142556131 kilocalorie / mole / radian ** 2  id: a20  note: funky pseudo-linear cyclic nitrogen  >
angles_handler.parameters[cyclic_oxygen_smirks]
<AngleType with smirks: C1C[C:1]-[O:2]-[C:3]C1  angle: 179 degree  k: 134.5019777341 kilocalorie / mole / radian ** 2  note: funky pseudo-linear cyclic oxygen  >

Cometic attributes are included by default when writing to disk. (This can be turned off by flipping the argument discard_cosmetic_attributes to True.) Let’s write out our modified force field and search through the file for the note we added.

force_field.to_file("modified.offxml")
!grep funky modified.offxml
        <Angle smirks="[*:1]~[#7X3$(*~[#6X3,#6X2,#7X2+0]):2]~[*:3]" angle="179 * degree ** 1" k="151.142556131 * mole ** -1 * radian ** -2 * kilocalorie ** 1" id="a20" note="funky pseudo-linear cyclic nitrogen"></Angle>
        <Angle smirks="C1C[C:1]-[O:2]-[C:3]C1" angle="179 * degree ** 1" k="134.5019777341 * kilocalorie ** 1 * mole ** -1 * radian ** -2" note="funky pseudo-linear cyclic oxygen"></Angle>