QMMM Submodule¶

`janus.qmmm.QMMM`(hl_wrapper, ll_wrapper, sys_info)	QMMM class for QMMM computations
`janus.qmmm.AQMMM`(hl_wrapper, ll_wrapper, …)	AQMMM super class for adaptive QMMM computations.
`janus.qmmm.OniomXS`(hl_wrapper, ll_wrapper, …)	Class for the Oniom-XS adaptive QM/MM method.
`janus.qmmm.HotSpot`(hl_wrapper, ll_wrapper, …)	Class for the Hot Spot adaptive QM/MM method.
`janus.qmmm.PAP`(hl_wrapper, ll_wrapper, sys_info)	Class for the PAP adaptive QM/MM method.
`janus.qmmm.SAP`(hl_wrapper, ll_wrapper, sys_info)	Class for the SAP adaptive QM/MM method.
`janus.qmmm.DAS`(hl_wrapper, ll_wrapper, sys_info)	Class for the Oniom-XS adaptive QM/MM method.

QMMM¶

class janus.qmmm.QMMM(hl_wrapper, ll_wrapper, sys_info, sys_info_format='pdb', qm_atoms=[], qmmm_scheme='subtractive', embedding_method='Mechanical', boundary_treatment='link_atom', link_atom_element='H')¶

Bases: object

QMMM class for QMMM computations

Parameters:

hl_wrapper (MMWrapper subclass or QMWrapper subclass) – Wrapper for performing the high-level computation. Traditionally QM but user can define MM.
ll_wrapper (MMWrapper subclass) – Wrapper for performing the low-level computation
sys_info (str) – A string with the filename or a list with multiple filenames that contain position and topology information.
sys_info_format (str) – Describes what kind of input is contained in sys_info. Default is pdb.
qm_atoms (list) – Indicies that define the QM region. Only static in traditional QM/MM
qmmm_scheme (str) – Scheme for computing QM/MM energies and gradients, only substractive available(default)
embedding_method (str) – Embedding method to use for QM/MM. Mechanical(default) and Electrostatic available
boundary_treatment (str) – Method for treating dangling bonds in the QM region, link_atom(default), RC, and RCD available
link_atom_element – Element to use for link atom, default is H. Beware of using others (not all functionality tested)

compute_gradients(system)¶

Computes the QM/MM gradients

Note

RCD gradients currently not implemented

Parameters:	system (`System`) – The system in which to save qmmm energy and forces

convert_input(fil, form)¶

Converts a set of input files into a MD trajectory

Parameters:	fil (str) – A string with the filename or a list with multiple filenames that contain position and topology information. form (str) – Describes what kind of input is contained in sys_info. Default is pdb.
Returns:
Return type:	MDtraj object

electrostatic(system, main_info)¶

Gets energies of needed components and computes a QM/MM energy with a subtractive electrostatic embedding scheme

E(QM/MM) = E(MM no coulomb)_entire_sys - E(MM no coulomb)_primary_subsys

E(QM)_primary_subsys + E(MM just coulomb)_secondary_subsys

Parameters:	system (`System`) – The system in which to save qmmm energy and forces main_info (dict) – contains the energy and forces for the whole system

find_boundary_bonds(qm_atoms=None)¶

Identified any covalent bonds that the QM/MM boundary cuts across

Parameters:	qm_atoms (list) – Indicies corresponding to the atoms in the primary subsystem. Default is None and uses self.qm_atoms

Examples

>>> find_boundary_bonds()
find_boundary_bonds(qm_atoms=[0,1,2,3])

get_external_charges(system)¶

Gets the point charges of atoms from secondary subsystem for electrostatic embedding

Parameters:	system (`System`) – The system in which to save qmmm energy and forces
Returns:	charges and corresponding positions in angstroms as xyz coordinates
Return type:	list

get_forces(run_ID=None)¶

Function to return qmmm forces

Parameters:	run_ID (int) – identifies which step to get forces from
Returns:	qmmm forces in au/bohr
Return type:	dict

Examples

>>> forces = get_forces()

get_redistributed_positions(positions, bonds, mm)¶

Gets the positions for the redistributed point charges in the RC and RCD schemes

Parameters:	positions (list) – bonds (list) – indices of all atoms (in secondary subsystem) bonded to M1 mm (int) – the index of M1
Returns:	positions for the redistributed charges
Return type:	list

make_primary_subsys_trajectory(qm_atoms=None)¶

Creates a MDtraj trajectory object with just the primary subsystem, and adds in any link atoms

Note

Currently adds just H as link atom, need to expand. Also, H is added as a very specific H1 atom, or else the connectivity in the create_new_residue_templates function in openmm gets messed up and gives a “set of atoms match but bonds are different error”

Parameters:	qm_atoms (list) – atom indicies corresponding to the atoms in the primary subsystem. Default is None and uses self.qm_atoms
Returns:	MDtraj trajectory object list – The link atom indices in traj

Examples

>>> make_primary_subsys_trajectory([0,1,2])
make_primary_subsys_trajectory()

make_second_subsys_trajectory(qm_atoms=None)¶

Creates a MDtraj trajectory object with just the secondary subsystem

Parameters:	qm_atoms (list) – atom indicies corresponding to the atoms in the primary subsystem. Default is None and uses self.qm_atoms
Returns:
Return type:	MDtraj trajectory object

Examples

>>> make_second_subsys_trajectory([0,1,2])
make_second_subsys_trajectory()

mechanical(system, main_info)¶

Gets energies of needed components and computes a QM/MM energy with a subtractive mechanical embedding scheme using the formula

E(QM/MM) = E(MM)_entire_sys - E(MM)_primary_subsys + E(QM)_primary_subsys

Parameters:	system (`System`) – The system in which to save qmmm energy and forces main_info (dict) – contains the energy and forces for the whole system

prepare_link_atom()¶: Identifies where to put link atom. Saves the qm and mm atom associated with the bond being cut across the QM/MM boundary and computes the g scaling factor for the link atom.

run_qmmm(main_info, wrapper_type)¶

Drives QM/MM computation. Updates the positions and topology given in main_info, and determines the QM/MM energy and gradients

Parameters:	main_info (dict) – Contains the energy, forces, topology, and position information for the whole system wrapper_type (str) – Defines the program used to obtain main_info

update_traj(position, topology, wrapper_type)¶

Updates the positions and topology of self.traj, a MDtraj trajectory object

Parameters:	positions (list) – positions in nm topology (OpenMM topology object) – If the mm program is OpenMM wrapper_type (str) – Defines the program used to obtain topology and positions

AQMMM¶

class janus.qmmm.AQMMM(hl_wrapper, ll_wrapper, sys_info, sys_info_format, qmmm_param, class_type)¶

Bases: abc.ABC, janus.qmmm.qmmm.QMMM

AQMMM super class for adaptive QMMM computations. Inherits from QMMM class.

Note

Since AQMMM is a super class and has abstract methods cannot actually instantiate AQMMM object, but only its child objects

`compute_COM`(atoms)	Computes the center of mass of a specified group
`compute_gradients`(system)	Computes the QM/MM gradients
`compute_lamda_i`(r_i)	Computes the switching function and the derivative of the switching function defined as a 5th order spline:
`compute_zero_energy`()	Compute the energy of the isolated groups at their minimum geometry
`convert_input`(fil, form)	Converts a set of input files into a MD trajectory
`define_buffer_zone`(qm_center)	Determines buffer group atoms.
`edit_atoms`(atoms, res_idx[, remove, add])	Edits a given list of atoms based on give parameters.
`electrostatic`(system, main_info)	Gets energies of needed components and computes a QM/MM energy with a subtractive electrostatic embedding scheme
`find_boundary_bonds`([qm_atoms])	Identified any covalent bonds that the QM/MM boundary cuts across
`find_buffer_atoms`(qm_center)	Find the buffer groups whose COM falls in between Rmin and Rmax
`get_Rmax`()	Function to return self.Rmin
`get_Rmin`()	Function to return self.Rmin
`get_external_charges`(system)	Gets the point charges of atoms from secondary subsystem for electrostatic embedding
`get_forces`([run_ID])	Function to return qmmm forces
`get_redistributed_positions`(positions, bonds, mm)	Gets the positions for the redistributed point charges in the RC and RCD schemes
`get_residue_info`(idx[, qm_center_xyz])	Gets the COM information and distance from the qm_center for a give residue.
`get_zero_energy`()	Incorporates the zero energy of groups to the total qmmm energy
`make_primary_subsys_trajectory`([qm_atoms])	Creates a MDtraj trajectory object with just the primary subsystem, and adds in any link atoms
`make_second_subsys_trajectory`([qm_atoms])	Creates a MDtraj trajectory object with just the secondary subsystem
`mechanical`(system, main_info)	Gets energies of needed components and computes a QM/MM energy with a subtractive mechanical embedding scheme using the formula
`partition`(info)	Function implemented in individual child classes
`prepare_link_atom`()	Identifies where to put link atom.
`run_aqmmm`()	Function implemented in individual child classes
`run_qmmm`(main_info, wrapper_type)	Drives QM/MM computation.
`set_Rmax`(Rmax)	Function to set self.Rmax
`set_Rmin`(Rmin)	Function to set self.Rmin
`update_traj`(position, topology, wrapper_type)	Updates the positions and topology of self.traj, a MDtraj trajectory object

Attributes

nm_to_angstrom

compute_COM(atoms)¶

Computes the center of mass of a specified group

Parameters:	atoms (list) – indices defining the group to compute the COM for
Returns:	numpy array – COM xyz coordinates dict – the indices of each atom in the group and the atomic weight for each atom dict – the indices of each atom in the group and the weight ratio of each atom to the total weight of the group (sum of all atomic weights)

compute_lamda_i(r_i)¶

Computes the switching function and the derivative of the switching function defined as a 5th order spline:

\lambda_i = -6x^5 + 15x^4 - 10x^3 + 1 d_{\lambda_i} = -30x^4 + 60x^3 - 30x^2

where x is the reduced distance (r_i - rmin)/(rmax - rmin) of buffer group i

Parameters:	r_i (float) – the distance between the qm center and the COM (in angstroms)
Returns:	float – lamda_i, unitless float – derivative of lamda_i, unitless

compute_zero_energy()¶: Compute the energy of the isolated groups at their minimum geometry

define_buffer_zone(qm_center)¶

Determines buffer group atoms. Gets the buffer groups in the buffer zone based on a distance partitioning scheme, saves each buffer group as a Buffer object, and saves all buffer groups in the dictionary self.buffer_groups. For water as a solvent, considers the whole water molecule as a buffer group.

Note

Currently only worked with explicit solvent based systems. Cannot treat buffer atoms that are part of large molecular structures (e.g., proteins).

Parameters:	qm_center (list) – the indicies that make up the qm center

edit_atoms(atoms, res_idx, remove=False, add=False)¶

Edits a given list of atoms based on give parameters.

Parameters:	atoms (list) – List of atom indicies to performed the desired action on res_idx (int) – Index of the residue remove (bool) – Whether to remove the atoms of residue res_idx from atoms. Default is False. add (bool) – Whether to add the atoms of residue res_idx to atoms Default is False.
Returns:	List of edited atoms
Return type:	list

Examples

>>> atoms = edit_qm_atoms(atoms=[0,1,2], res_idx=0, remove=True)

find_buffer_atoms(qm_center)¶

Find the buffer groups whose COM falls in between Rmin and Rmax

Parameters:	qm_center (list) – the indicies that make up the qm center

get_Rmax()¶

Function to return self.Rmin

Returns:	the distance from qm center to outer limit of buffer zone in angstroms
Return type:	float

get_Rmin()¶

Function to return self.Rmin

Returns:	the distance from qm center to inner limit of buffer zone in angstroms
Return type:	float

get_residue_info(idx, qm_center_xyz=None)¶

Gets the COM information and distance from the qm_center for a give residue. Saves the information in a Buffer object.

Parameters:	idx (int) – index of the residue qm_center_xyz (list) – XYZ coordinates of the qm_center as a list
Returns:
Return type:	`Buffer`

get_zero_energy()¶: Incorporates the zero energy of groups to the total qmmm energy

nm_to_angstrom = 10.0¶

partition(info)¶: Function implemented in individual child classes

run_aqmmm()¶: Function implemented in individual child classes

run_qmmm(main_info, wrapper_type)¶

Drives QM/MM computation. Updates the positions and topology given in main_info, determines the partitions to be computed, for each partition, determines the QM/MM energy and gradients, then interpolates all partitions using specified adaptive QM/MM scheme

Parameters:	main_info (dict) – Contains the energy, forces, topology, and position information for the whole system wrapper_type (str) – Defines the program used to obtain main_info

set_Rmax(Rmax)¶

Function to set self.Rmax

Parameters:	Rmax (float) – the distance from qm center to outer limit of buffer zone in angstroms

set_Rmin(Rmin)¶

Function to set self.Rmin

Parameters:	Rmin (float) – the distance from qm center to inner limit of buffer zone in angstroms

ONIOM-XS¶

class janus.qmmm.OniomXS(hl_wrapper, ll_wrapper, sys_info, sys_info_format='pdb', qm_center=[0], partition_scheme='distance', Rmin=3.8, Rmax=4.5, qmmm_param={}, **kwargs)¶

Bases: janus.qmmm.aqmmm.AQMMM

Class for the Oniom-XS adaptive QM/MM method. Inherits from AQMMM class

Parameters:

hl_wrapper (MMWrapper subclass or QMWrapper subclass) – Wrapper for performing the high-level computation. Traditionally QM but user can define MM.
ll_wrapper (MMWrapper subclass) – Wrapper for performing the low-level computation
sys_info (str) – A string with the filename or a list with multiple filenames that contain position and topology information.
sys_info_format (str) – Describes what kind of input is contained in sys_info. Default is pdb.
qm_center (list) – Atoms that define the qm center, default is [0]. If more than one index is given, COM is used as qm_center
partition_scheme (str) – Scheme to use to define buffer groups, default is distance (only scheme available as of now)
Rmin (float) – Inner radius for distance partition in angstroms, default is 3.8
Rmax (float) – Outer radius for distance partition in angstroms, default is 4.5
qmmm_param – A dictionary with any parameters to pass into the QMMM class. See QMMM class for specifics

get_switching_function()¶

Averages the individual switching functions of each buffer group

Returns:	The average of the switching functions
Return type:	float

partition(qm_center=None)¶

Finds the partitions as required by the ONIOM-XS method and saves each partition as a system object. Saves all systems in the dictionary self.systems

Parameters:	qm_center (list) – atoms that define the qm center, default is None

run_aqmmm()¶: Interpolates the energy and gradients from each partition according to the ONIOM-XS method

Hot Spot¶

class janus.qmmm.HotSpot(hl_wrapper, ll_wrapper, sys_info, sys_info_format='pdb', qm_center=[0], partition_scheme='distance', Rmin=3.8, Rmax=4.5, qmmm_param={}, **kwargs)¶

Bases: janus.qmmm.aqmmm.AQMMM

Class for the Hot Spot adaptive QM/MM method. Inherits from AQMMM class

Parameters:

hl_wrapper (MMWrapper subclass or QMWrapper subclass) – Wrapper for performing the high-level computation. Traditionally QM but user can define MM.
ll_wrapper (MMWrapper subclass) – Wrapper for performing the low-level computation
sys_info (str) – A string with the filename or a list with multiple filenames that contain position and topology information.
sys_info_format (str) – Describes what kind of input is contained in sys_info. Default is pdb.
qm_center (list) – Atoms that define the qm center, default is [0]. If more than one index is given, COM is used as qm_center
partition_scheme (str) – Scheme to use to define buffer groups, default is distance (only scheme available as of now)
Rmin (float) – Inner radius for distance partition in angstroms, default is 3.8
Rmax (float) – Outer radius for distance partition in angstroms, default is 4.5
qmmm_param – A dictionary with any parameters to pass into the QMMM class. See QMMM class for specifics

compute_lamda_i(r_i)¶

Computes the switching function of the Hot-Spot method and overrides the compute_lamda_i function from the AQMMM class

Parameters:	r_i (float) – the distance between the qm center and the COM in angstroms
Returns:	float – lamda_i, unitless None – placeholder for the derivative of r_i

partition(qm_center=None)¶

Finds the partitions as required by the Hot-Spot method and saves each partition as a system object. Saves all systems in the dictionary self.systems

Parameters:	qm_center (list) – atoms that define the qm center, default is None

run_aqmmm()¶: Interpolates the energy and gradients from each partition according to the Hot-Spot method

PAP¶

class janus.qmmm.PAP(hl_wrapper, ll_wrapper, sys_info, sys_info_format='pdb', modified_variant=False, qm_center=[0], partition_scheme='distance', Rmin=3.8, Rmax=4.5, qmmm_param={}, **kwargs)¶

Bases: janus.qmmm.aqmmm.AQMMM

Class for the PAP adaptive QM/MM method. Inherits from AQMMM class

Parameters:

hl_wrapper (MMWrapper subclass or QMWrapper subclass) – Wrapper for performing the high-level computation. Traditionally QM but user can define MM.
ll_wrapper (MMWrapper subclass) – Wrapper for performing the low-level computation
sys_info (str) – A string with the filename or a list with multiple filenames that contain position and topology information.
sys_info_format (str) – Describes what kind of input is contained in sys_info. Default is pdb.
modified_variant (bool) – Whether to use the modified version mPAP, which disregards gradient terms that come from the switching function. Default is False.
qm_center (list) – Atoms that define the qm center, default is [0]. If more than one index is given, COM is used as qm_center
partition_scheme (str) – Scheme to use to define buffer groups, default is distance (only scheme available as of now)
Rmin (float) – Inner radius for distance partition in angstroms, default is 3.8
Rmax (float) – Outer radius for distance partition in angstroms, default is 4.5
qmmm_param – A dictionary with any parameters to pass into the QMMM class. See QMMM class for specifics

compute_sf_gradient()¶

Computes forces due to the gradient of the switching function

Returns:	Forces due to gradient of switching function
Return type:	dict

get_combos(items=None)¶

Gets all combinations of a given list of indices

Parameters:	items (list) – indices to get combinations for
Returns:	all possible combinations
Return type:	list

Examples

>>> combos = get_combos([1,2])

In this case, combos will return >>> [(1), (2), (1,2)]

partition(qm_center=None)¶

Finds the partitions as required by the PAP method and saves each partition as a system object. Saves all systems in the dictionary self.systems

Parameters:	qm_center (list) – atoms that define the qm center, default is None

run_aqmmm()¶: Interpolates the energy and gradients from each partition according to the PAP method

SAP¶

class janus.qmmm.SAP(hl_wrapper, ll_wrapper, sys_info, sys_info_format='pdb', modified_variant=False, qm_center=[0], partition_scheme='distance', Rmin=3.8, Rmax=4.5, qmmm_param={}, **kwargs)¶

Bases: janus.qmmm.aqmmm.AQMMM

Class for the SAP adaptive QM/MM method. Inherits from AQMMM class

Parameters:

hl_wrapper (MMWrapper subclass or QMWrapper subclass) – Wrapper for performing the high-level computation. Traditionally QM but user can define MM.
ll_wrapper (MMWrapper subclass) – Wrapper for performing the low-level computation
sys_info (str) – A string with the filename or a list with multiple filenames that contain position and topology information.
sys_info_format (str) – Describes what kind of input is contained in sys_info. Default is pdb.
modified_variant (bool) – Whether to use the modified version mPAP, which disregards gradient terms that come from the switching function. Default is False.
qm_center (list) – Atoms that define the qm center, default is [0]. If more than one index is given, COM is used as qm_center
partition_scheme (str) – Scheme to use to define buffer groups, default is distance (only scheme available as of now)
Rmin (float) – Inner radius for distance partition in angstroms, default is 3.8
Rmax (float) – Outer radius for distance partition in angstroms, default is 4.5
qmmm_param – A dictionary with any parameters to pass into the QMMM class. See QMMM class for specifics

compute_sf_gradient()¶

Computes forces due to the gradient of the switching function

Returns:	Forces due to gradient of switching function
Return type:	dict

get_combos(items=None, buffer_distance=None)¶

A given list of indices is sorted by distance and the combinations is found based on distance order. See example

Parameters:	items (list) – indices to get combinations for buffer_distance (dict) – indices from items and their distances from the qm center
Returns:	combinations
Return type:	list

Examples

>>> combos = get_combos([1,2,3])

If 1 was closer to the center than 2, and 2 closer to the center than 3, functions returns following: >>> [(1), (1,2), (1,2,3)]

get_switching_functions()¶: Computes switching function for SAP computations and saves to each buffer group object

partition(qm_center=None)¶

Finds the partitions as required by the SAP method and saves each partition as a system object. Saves all systems in the dictionary self.systems

Parameters:	qm_center (list) – atoms that define the qm center, default is None

run_aqmmm()¶: Interpolates the energy and gradients from each partition according to the SAP method

DAS¶

class janus.qmmm.DAS(hl_wrapper, ll_wrapper, sys_info, sys_info_format='pdb', qm_center=[0], partition_scheme='distance', Rmin=3.8, Rmax=4.5, qmmm_param={}, **kwargs)¶

Bases: janus.qmmm.aqmmm.AQMMM

Class for the Oniom-XS adaptive QM/MM method. Inherits from AQMMM class

Parameters:

hl_wrapper (MMWrapper subclass or QMWrapper subclass) – Wrapper for performing the high-level computation. Traditionally QM but user can define MM.
ll_wrapper (MMWrapper subclass) – Wrapper for performing the low-level computation
sys_info (str) – A string with the filename or a list with multiple filenames that contain position and topology information.
sys_info_format (str) – Describes what kind of input is contained in sys_info. Default is pdb.
qm_center (list) – Atoms that define the qm center, default is [0]. If more than one index is given, COM is used as qm_center
partition_scheme (str) – Scheme to use to define buffer groups, default is distance (only scheme available as of now)
Rmin (float) – Inner radius for distance partition in angstroms, default is 3.8
Rmax (float) – Outer radius for distance partition in angstroms, default is 4.5
qmmm_param – A dictionary with any parameters to pass into the QMMM class. See QMMM class for specifics

compute_lamda_i(r_i)¶

Computes the switching function of the DAS method and overrides the compute_lamda_i function from the AQMMM class

Parameters:	r_i (float) – the distance between the qm center and the COM in angstroms
Returns:	float – lamda_i, unitless None – placeholder for the derivative of r_i

get_combos(items=None)¶

Gets all combinations of a given list of indices according to the DAS formulation

Parameters:	items (list) – indices to get combinations for
Returns:	combinations
Return type:	list

partition(qm_center=None)¶

Finds the partitions as required by the PAP method and saves each partition as a system object. Saves all systems in the dictionary self.systems

Parameters:	qm_center (list) – atoms that define the qm center, default is None

run_aqmmm()¶: Interpolates the energy and gradients from each partition according to the PAP method