Data Structures

The primary building blocks in an mBuild hierarchy inherit from the mbuild.Compound class. Compounds maintain an ordered set of children which are other Compounds. In addition, an independent, ordered dictionary of labels is maintained through which users can reference any other Compound in the hierarchy via descriptive strings. Every Compound knows its parent Compound, one step up in the hierarchy, and knows which Compounds reference it in their labels. mbuild.Port is a special type of Compound which are used internally to connect different Compounds using the equivalence transformations described below.

Compounds at the bottom of an mBuild hierarchy, the leaves of the tree, are referred to as Particles and can be instantiated as foo = mbuild.Particle(name='bar'). Note however, that this merely serves to illustrate that this Compound is at the bottom of the hierarchy; Particle is simply an alias for Compound which can be used to clarify the intended role of an object you are creating. The method mbuild.Compound.particles() traverses the hierarchy to the bottom and yields those Compounds. mbuild.Compound.root() returns the compound at the top of the hierarchy.

Compound

class mbuild.Compound(subcompounds=None, name=None, pos=None, mass=0.0, charge=0.0, periodicity=None, box=None, element=None, port_particle=False)[source]

A building block in the mBuild hierarchy.

Compound is the superclass of all composite building blocks in the mBuild hierarchy. That is, all composite building blocks must inherit from compound, either directly or indirectly. The design of Compound follows the Composite design pattern:

@book{DesignPatterns,
    author = "Gamma, Erich and Helm, Richard and Johnson, Ralph and
    Vlissides, John M.",
    title = "Design Patterns",
    subtitle = "Elements of Reusable Object-Oriented Software",
    year = "1995",
    publisher = "Addison-Wesley",
    note = "p. 395",
    ISBN = "0-201-63361-2",
}

with Compound being the composite, and Particle playing the role of the primitive (leaf) part, where Particle is in fact simply an alias to the Compound class.

Compound maintains a list of children (other Compounds contained within), and provides a means to tag the children with labels, so that the compounds can be easily looked up later. Labels may also point to objects outside the Compound’s containment hierarchy. Compound has built-in support for copying and deepcopying Compound hierarchies, enumerating particles or bonds in the hierarchy, proximity based searches, visualization, I/O operations, and a number of other convenience methods.

Parameters
subcompoundsmb.Compound or list of mb.Compound, optional, default=None

One or more compounds to be added to self.

namestr, optional, default=self.__class__.__name__

The type of Compound.

posnp.ndarray, shape=(3,), dtype=float, optional, default=[0, 0, 0]

The position of the Compound in Cartestian space

massfloat, optional, default=0.0

The mass of the compound. If none is set, then will try to infer the mass from a compound’s element attribute. If neither mass or element are specified, then the mass will be zero.

chargefloat, optional, default=0.0

Currently not used. Likely removed in next release.

periodicitytuple of bools, length=3, optional, default=None

Whether the Compound is periodic in the x, y, and z directions. If None is provided, the periodicity is set to (False, False, False) which is non-periodic in all directions.

port_particlebool, optional, default=False

Whether or not this Compound is part of a Port

boxmb.Box, optional

The simulation box containing the compound. Also accounts for the periodicity. Defaults to None which is treated as non-periodic.

element: str, optional, default=None

The one or two character element symbol

Attributes
bond_graphmb.BondGraph

Graph-like object that stores bond information for this Compound

childrenOrderedSet

Contains all children (other Compounds).

labelsOrderedDict

Labels to Compound/Atom mappings. These do not necessarily need not be in self.children.

parentmb.Compound

The parent Compound that contains this part. Can be None if this compound is the root of the containment hierarchy.

referrersset

Other compounds that reference this part with labels.

rigid_idint, default=None

Get the rigid_id of the Compound.

boundingboxmb.Box

The bounds (xmin, xmax, ymin, ymax, zmin, zmax) of particles in Compound

center

Get the cartesian center of the Compound based on its Particles.

contains_rigid

Return True if the Compound contains rigid bodies.

mass

Return the total mass of a compound.

max_rigid_id

Return the maximum rigid body ID contained in the Compound.

n_particles

Return the number of Particles in the Compound.

n_bonds

Return the number of bonds in the Compound.

root

Get the Compound at the top of self’s hierarchy.

xyz

Return all particle coordinates in this compound.

xyz_with_ports

Return all particle coordinates in this compound including ports.

add(new_child, label=None, containment=True, replace=False, inherit_periodicity=None, inherit_box=False, reset_rigid_ids=True)[source]

Add a part to the Compound.

Note:

This does not necessarily add the part to self.children but may instead be used to add a reference to the part to self.labels. See ‘containment’ argument.

Parameters
new_childmb.Compound or list-like of mb.Compound

The object(s) to be added to this Compound.

labelstr, optional, default None

A descriptive string for the part.

containmentbool, optional, default=True

Add the part to self.children.

replacebool, optional, default=True

Replace the label if it already exists.

inherit_periodicitybool, optional, default=True

Replace the periodicity of self with the periodicity of the Compound being added

inherit_box: bool, optional, default=False

Replace the box of self with the box of the Compound being added

reset_rigid_idsbool, optional, default=True

If the Compound to be added contains rigid bodies, reset the rigid_ids such that values remain distinct from rigid_ids already present in self. Can be set to False if attempting to add Compounds to an existing rigid body.

add_bond(particle_pair)[source]

Add a bond between two Particles.

Parameters
particle_pairindexable object, length=2, dtype=mb.Compound

The pair of Particles to add a bond between

all_ports()[source]

Return all Ports referenced by this Compound and its successors.

Returns
list of mb.Compound

A list of all Ports referenced by this Compound and its successors

ancestors()[source]

Generate all ancestors of the Compound recursively.

Yields
mb.Compound

The next Compound above self in the hierarchy

available_ports()[source]

Return all unoccupied Ports referenced by this Compound.

Returns
list of mb.Compound

A list of all unoccupied ports referenced by the Compound

bonds()[source]

Return all bonds in the Compound and sub-Compounds.

Yields
tuple of mb.Compound

The next bond in the Compound

See also

bond_graph.edges_iter

Iterates over all edges in a BondGraph

property box

Get the box of the Compound.

Ports cannot have a box.

property center

Get the cartesian center of the Compound based on its Particles.

Returns
np.ndarray, shape=(3,), dtype=float

The cartesian center of the Compound based on its Particles

property charge

Get the charge of the Compound.

property contains_rigid

Return True if the Compound contains rigid bodies.

If the Compound contains any particle with a rigid_id != None then contains_rigid will return True. If the Compound has no children (i.e. the Compound resides at the bottom of the containment hierarchy) then contains_rigid will return False.

Returns
bool,

True if the Compound contains any particle with a rigid_id != None

Notes

The private variable ‘_check_if_contains_rigid_bodies’ is used to help cache the status of ‘contains_rigid’. If ‘_check_if_contains_rigid_bodies’ is False, then the rigid body containment of the Compound has not changed, and the particle tree is not traversed, boosting performance.

property element

Get the element of the Compound.

energy_minimize(forcefield='UFF', steps=1000, **kwargs)[source]

Perform an energy minimization on a Compound.

Default behavior utilizes Open Babel to perform an energy minimization/geometry optimization on a Compound by applying a generic force field

Can also utilize OpenMM to energy minimize after atomtyping a Compound using Foyer to apply a forcefield XML file that contains valid SMARTS strings.

This function is primarily intended to be used on smaller components, with sizes on the order of 10’s to 100’s of particles, as the energy minimization scales poorly with the number of particles.

Parameters
stepsint, optional, default=1000

The number of optimization iterations

forcefieldstr, optional, default=’UFF’

The generic force field to apply to the Compound for minimization. Valid options are ‘MMFF94’, ‘MMFF94s’, ‘’UFF’, ‘GAFF’, ‘Ghemical’. Please refer to the Open Babel documentation when considering your choice of force field. Utilizing OpenMM for energy minimization requires a forcefield XML file with valid SMARTS strings. Please refer to OpenMM docs for more information.

Other Parameters
algorithmstr, optional, default=’cg’

The energy minimization algorithm. Valid options are ‘steep’, ‘cg’, and ‘md’, corresponding to steepest descent, conjugate gradient, and equilibrium molecular dynamics respectively. For _energy_minimize_openbabel

scale_bondsfloat, optional, default=1

Scales the bond force constant (1 is completely on). For _energy_minimize_openmm

scale_anglesfloat, optional, default=1

Scales the angle force constant (1 is completely on) For _energy_minimize_openmm

scale_torsionsfloat, optional, default=1

Scales the torsional force constants (1 is completely on) For _energy_minimize_openmm Note: Only Ryckaert-Bellemans style torsions are currently supported

scale_nonbondedfloat, optional, default=1

Scales epsilon (1 is completely on) For _energy_minimize_openmm

constraintsstr, optional, default=”AllBonds”

Specify constraints on the molecule to minimize, options are: None, “HBonds”, “AllBonds”, “HAngles” For _energy_minimize_openmm

References

If using _energy_minimize_openmm(), please cite:

Eastman2013

P. Eastman, M. S. Friedrichs, J. D. Chodera, R. J. Radmer, C. M. Bruns, J. P. Ku, K. A. Beauchamp, T. J. Lane, L.-P. Wang, D. Shukla, T. Tye, M. Houston, T. Stich, C. Klein, M. R. Shirts, and V. S. Pande. “OpenMM 4: A Reusable, Extensible, Hardware Independent Library for High Performance Molecular Simulation.” J. Chem. Theor. Comput. 9(1): 461-469. (2013).

If using _energy_minimize_openbabel(), please cite:

OBoyle2011

O’Boyle, N.M.; Banck, M.; James, C.A.; Morley, C.; Vandermeersch, T.; Hutchison, G.R. “Open Babel: An open chemical toolbox.” (2011) J. Cheminf. 3, 33

OpenBabel

Open Babel, version X.X.X http://openbabel.org, (installed Month Year)

If using the ‘MMFF94’ force field please also cite the following:

Halgren1996a

T.A. Halgren, “Merck molecular force field. I. Basis, form, scope, parameterization, and performance of MMFF94.” (1996) J. Comput. Chem. 17, 490-519

Halgren1996b

T.A. Halgren, “Merck molecular force field. II. MMFF94 van der Waals and electrostatic parameters for intermolecular interactions.” (1996) J. Comput. Chem. 17, 520-552

Halgren1996c

T.A. Halgren, “Merck molecular force field. III. Molecular geometries and vibrational frequencies for MMFF94.” (1996) J. Comput. Chem. 17, 553-586

Halgren1996d

T.A. Halgren and R.B. Nachbar, “Merck molecular force field. IV. Conformational energies and geometries for MMFF94.” (1996) J. Comput. Chem. 17, 587-615

Halgren1996e

T.A. Halgren, “Merck molecular force field. V. Extension of MMFF94 using experimental data, additional computational data, and empirical rules.” (1996) J. Comput. Chem. 17, 616-641

If using the ‘MMFF94s’ force field please cite the above along with:

Halgren1999

T.A. Halgren, “MMFF VI. MMFF94s option for energy minimization studies.” (1999) J. Comput. Chem. 20, 720-729

If using the ‘UFF’ force field please cite the following:

Rappe1992

Rappe, A.K., Casewit, C.J., Colwell, K.S., Goddard, W.A. III, Skiff, W.M. “UFF, a full periodic table force field for molecular mechanics and molecular dynamics simulations.” (1992) J. Am. Chem. Soc. 114, 10024-10039

If using the ‘GAFF’ force field please cite the following:

Wang2004

Wang, J., Wolf, R.M., Caldwell, J.W., Kollman, P.A., Case, D.A. “Development and testing of a general AMBER force field” (2004) J. Comput. Chem. 25, 1157-1174

If using the ‘Ghemical’ force field please cite the following:

Hassinen2001

T. Hassinen and M. Perakyla, “New energy terms for reduced protein models implemented in an off-lattice force field” (2001) J. Comput. Chem. 22, 1229-1242

freud_generate_bonds(name_a, name_b, dmin, dmax, exclude_ii=True)[source]

Add Bonds between all pairs of types a/b within [dmin, dmax].

Parameters
name_astr

The name of one of the Particles to be in each bond

name_bstr

The name of the other Particle to be in each bond

dminfloat

The minimum distance (in nm) between Particles for considering a bond

dmaxfloat

The maximum distance (in nm) between Particles for considering a bond

exclude_iibool, optional, default=True

Whether or not to include neighbors with the same index.

Notes

This is an experimental feature and some behavior might change out of step of a standard development release.

from_gmso(topology, coords_only=False, infer_hierarchy=True)[source]

Convert a GMSO Topology to mBuild Compound.

Returns
compoundmb.Compound
from_parmed(structure, coords_only=False, infer_hierarchy=True)[source]

Extract atoms and bonds from a pmd.Structure.

Will create sub-compounds for every chain if there is more than one and sub-sub-compounds for every residue.

Parameters
structurepmd.Structure

The structure to load.

coords_onlybool

Set preexisting atoms in compound to coordinates given by structure.

infer_hierarchybool, optional, default=True

If true, infer compound hierarchy from chains and residues

from_pybel(pybel_mol, use_element=True, coords_only=False, infer_hierarchy=True, ignore_box_warn=False)[source]

Create a Compound from a Pybel.Molecule.

Parameters
pybel_mol: pybel.Molecule
use_elementbool, default True

If True, construct mb Particles based on the pybel Atom’s element. If False, construcs mb Particles based on the pybel Atom’s type

coords_onlybool, default False

Set preexisting atoms in compound to coordinates given by structure. Note: Not yet implemented, included only for parity with other conversion functions

infer_hierarchybool, optional, default=True

If True, infer hierarchy from residues

ignore_box_warnbool, optional, default=False

If True, ignore warning if no box is present.

See also

mbuild.conversion.from_pybel
from_trajectory(traj, frame=- 1, coords_only=False, infer_hierarchy=True)[source]

Extract atoms and bonds from a md.Trajectory.

Will create sub-compounds for every chain if there is more than one and sub-sub-compounds for every residue.

Parameters
trajmdtraj.Trajectory

The trajectory to load.

frameint, optional, default=-1 (last)

The frame to take coordinates from.

coords_onlybool, optional, default=False

Only read coordinate information

infer_hierarchybool, optional, default=True

If True, infer compound hierarchy from chains and residues

See also

mbuild.conversion.from_trajectory
generate_bonds(name_a, name_b, dmin, dmax)[source]

Add Bonds between all pairs of types a/b within [dmin, dmax].

Parameters
name_astr

The name of one of the Particles to be in each bond

name_bstr

The name of the other Particle to be in each bond

dminfloat

The minimum distance (in nm) between Particles for considering a bond

dmaxfloat

The maximum distance (in nm) between Particles for considering a bond

get_boundingbox(pad_box=None)[source]

Compute the bounding box of the compound.

Compute and store the rectangular bounding box of the Compound.

Parameters
pad_box: Sequence, optional, default=None

Pad all lengths or a list of lengths by a specified amount in nm. Acceptable values are:

  • A single float: apply this pad value to all 3 box lengths.

  • A sequence of length 1: apply this pad value to all 3 box lengths.

  • A sequence of length 3: apply these pad values to the a, b, c box lengths.

Returns
mb.Box

The bounding box for this Compound.

Notes

Triclinic bounding boxes are supported, but only for Compounds that are generated from mb.Lattice’s and the resulting mb.Lattice.populate method

get_smiles()[source]

Get SMILES string for compound.

Bond order is guessed with pybel and may lead to incorrect SMILES strings.

Returns
smiles_string: str
label_rigid_bodies(discrete_bodies=None, rigid_particles=None)[source]

Designate which Compounds should be treated as rigid bodies.

If no arguments are provided, this function will treat the compound as a single rigid body by providing all particles in self with the same rigid_id. If discrete_bodies is not None, each instance of a Compound with a name found in discrete_bodies will be treated as a unique rigid body. If rigid_particles is not None, only Particles (Compounds at the bottom of the containment hierarchy) matching this name will be considered part of the rigid body.

Parameters
discrete_bodiesstr or list of str, optional, default=None

Name(s) of Compound instances to be treated as unique rigid bodies. Compound instances matching this (these) name(s) will be provided with unique rigid_ids

rigid_particlesstr or list of str, optional, default=None

Name(s) of Compound instances at the bottom of the containment hierarchy (Particles) to be included in rigid bodies. Only Particles matching this (these) name(s) will have their rigid_ids altered to match the rigid body number.

Examples

Creating a rigid benzene

>>> import mbuild as mb
>>> from mbuild.utils.io import get_fn
>>> benzene = mb.load(get_fn('benzene.mol2'))
>>> benzene.label_rigid_bodies()

Creating a semi-rigid benzene, where only the carbons are treated as a rigid body

>>> import mbuild as mb
>>> from mbuild.utils.io import get_fn
>>> benzene = mb.load(get_fn('benzene.mol2'))
>>> benzene.label_rigid_bodies(rigid_particles='C')

Create a box of rigid benzenes, where each benzene has a unique rigid body ID.

>>> import mbuild as mb
>>> from mbuild.utils.io import get_fn
>>> benzene = mb.load(get_fn('benzene.mol2'))
>>> benzene.name = 'Benzene'
>>> filled = mb.fill_box(benzene,
...                      n_compounds=10,
...                      box=[0, 0, 0, 4, 4, 4])
>>> filled.label_rigid_bodies(distinct_bodies='Benzene')

Create a box of semi-rigid benzenes, where each benzene has a unique rigid body ID and only the carbon portion is treated as rigid.

>>> import mbuild as mb
>>> from mbuild.utils.io import get_fn
>>> benzene = mb.load(get_fn('benzene.mol2'))
>>> benzene.name = 'Benzene'
>>> filled = mb.fill_box(benzene,
...                      n_compounds=10,
...                      box=[0, 0, 0, 4, 4, 4])
>>> filled.label_rigid_bodies(distinct_bodies='Benzene',
...                           rigid_particles='C')
property mass

Return the total mass of a compound.

If the compound contains children compouds, the total mass of all children compounds is returned. If the compound contains element information (Compound.element) then the mass is inferred from the elemental mass. If Compound.mass has been set explicitly, then it will override the mass inferred from Compound.element. If neither of a Compound’s element or mass attributes have been set, then a mass of zero is returned.

property max_rigid_id

Return the maximum rigid body ID contained in the Compound.

This is usually used by compound.root to determine the maximum rigid_id in the containment hierarchy.

Returns
int or None

The maximum rigid body ID contained in the Compound. If no rigid body IDs are found, None is returned

property maxs

Return the maximum x, y, z coordinate of any particle in this compound.

min_periodic_distance(xyz0, xyz1)[source]

Vectorized distance calculation considering minimum image.

Only implemented for orthorhombic simulation boxes.

Parameters
xyz0np.ndarray, shape=(3,), dtype=float

Coordinates of first point

xyz1np.ndarray, shape=(3,), dtype=float

Coordinates of second point

Returns
float

Vectorized distance between the two points following minimum image convention

property mins

Return the mimimum x, y, z coordinate of any particle in this compound.

property n_bonds

Return the number of bonds in the Compound.

Returns
int

The number of bonds in the Compound

property n_particles

Return the number of Particles in the Compound.

Returns
int,

The number of Particles in the Compound

particles(include_ports=False)[source]

Return all Particles of the Compound.

Parameters
include_portsbool, optional, default=False

Include port particles

Yields
mb.Compound

The next Particle in the Compound

particles_by_element(element)[source]

Return all Particles of the Compound with a specific element.

Parameters
namestr or ele.Element

element abbreviation or element

Yields
mb.Compound

The next Particle in the Compound with the user-specified element

particles_by_name(name)[source]

Return all Particles of the Compound with a specific name.

Parameters
namestr

Only particles with this name are returned

Yields
mb.Compound

The next Particle in the Compound with the user-specified name

particles_in_range(compound, dmax, max_particles=20, particle_kdtree=None, particle_array=None)[source]

Find particles within a specified range of another particle.

Parameters
compoundmb.Compound

Reference particle to find other particles in range of

dmaxfloat

Maximum distance from ‘compound’ to look for Particles

max_particlesint, optional, default=20

Maximum number of Particles to return

particle_kdtreemb.PeriodicKDTree, optional

KD-tree for looking up nearest neighbors. If not provided, a KD- tree will be generated from all Particles in self

particle_arraynp.ndarray, shape=(n,), dtype=mb.Compound, optional

Array of possible particles to consider for return. If not provided, this defaults to all Particles in self

Returns
np.ndarray, shape=(n,), dtype=mb.Compound

Particles in range of compound according to user-defined limits

See also

periodic_kdtree.PerioidicKDTree

mBuild implementation of kd-trees

scipy.spatial.kdtree

Further details on kd-trees

property periodicity

Get the periodicity of the Compound.

property pos

Get the position of the Compound.

If the Compound contains children, returns the center.

The position of a Compound containing children can’t be set.

referenced_ports()[source]

Return all Ports referenced by this Compound.

Returns
list of mb.Compound

A list of all ports referenced by the Compound

remove(objs_to_remove)[source]

Remove children from the Compound cleanly.

Parameters
objs_to_removemb.Compound or list of mb.Compound

The Compound(s) to be removed from self

remove_bond(particle_pair)[source]

Delete a bond between a pair of Particles.

Parameters
particle_pairindexable object, length=2, dtype=mb.Compound

The pair of Particles to remove the bond between

property rigid_id

Get the rigid_id of the Compound.

rigid_particles(rigid_id=None)[source]

Generate all particles in rigid bodies.

If a rigid_id is specified, then this function will only yield particles with a matching rigid_id.

Parameters
rigid_idint, optional

Include only particles with this rigid body ID

Yields
mb.Compound

The next particle with a rigid_id that is not None, or the next particle with a matching rigid_id if specified

property root

Get the Compound at the top of self’s hierarchy.

Returns
mb.Compound

The Compound at the top of self’s hierarchy

rotate(theta, around)[source]

Rotate Compound around an arbitrary vector.

Parameters
thetafloat

The angle by which to rotate the Compound, in radians.

aroundnp.ndarray, shape=(3,), dtype=float

The vector about which to rotate the Compound.

save(filename, show_ports=False, forcefield_name=None, forcefield_files=None, forcefield_debug=False, box=None, overwrite=False, residues=None, combining_rule='lorentz', foyer_kwargs=None, **kwargs)[source]

Save the Compound to a file.

Parameters
filenamestr

Filesystem path in which to save the trajectory. The extension or prefix will be parsed and control the format. Supported extensions: ‘hoomdxml’, ‘gsd’, ‘gro’, ‘top’, ‘lammps’, ‘lmp’, ‘mcf’

show_portsbool, optional, default=False

Save ports contained within the compound.

forcefield_filesstr, optional, default=None

Apply a forcefield to the output file using a forcefield provided by the foyer package.

forcefield_namestr, optional, default=None

Apply a named forcefield to the output file using the foyer package, e.g. ‘oplsaa’. Foyer forcefields

forcefield_debugbool, optional, default=False

Choose verbosity level when applying a forcefield through foyer. Specifically, when missing atom types in the forcefield xml file, determine if the warning is condensed or verbose.

boxmb.Box, optional, default=self.boundingbox (with buffer)

Box information to be written to the output file. If ‘None’, a bounding box is used with 0.25nm buffers at each face to avoid overlapping atoms.

overwritebool, optional, default=False

Overwrite if the filename already exists

residuesstr of list of str

Labels of residues in the Compound. Residues are assigned by checking against Compound.name.

combining_rulestr, optional, default=’lorentz’

Specify the combining rule for nonbonded interactions. Only relevant when the foyer package is used to apply a forcefield. Valid options are ‘lorentz’ and ‘geometric’, specifying Lorentz-Berthelot and geometric combining rules respectively.

foyer_kwargsdict, optional, default=None

Keyword arguments to provide to foyer.Forcefield.apply. Depending on the file extension these will be passed to either write_gsd, write_hoomdxml, write_lammpsdata, write_mcf, or parmed.Structure.save. See parmed structure documentation

Other Parameters
ref_distancefloat, optional, default=1.0

Normalization factor used when saving to .gsd and .hoomdxml formats for converting distance values to reduced units.

ref_energyfloat, optional, default=1.0

Normalization factor used when saving to .gsd and .hoomdxml formats for converting energy values to reduced units.

ref_massfloat, optional, default=1.0

Normalization factor used when saving to .gsd and .hoomdxml formats for converting mass values to reduced units.

atom_style: str, default=’full’

Defines the style of atoms to be saved in a LAMMPS data file. The following atom styles are currently supported: ‘full’, ‘atomic’, ‘charge’, ‘molecular’ See LAMMPS atom style documentation for more information.

unit_style: str, default=’real’

Defines to unit style to be save in a LAMMPS data file. Defaults to ‘real’ units. Current styles are supported: ‘real’, ‘lj’. See LAMMPS unit style documentation_ for more information.

See also

conversion.save

Main saver logic

formats.gsdwrite.write_gsd

Write to GSD format

formats.hoomdxml.write_hoomdxml

Write to Hoomd XML format

formats.xyzwriter.write_xyz

Write to XYZ format

formats.lammpsdata.write_lammpsdata

Write to LAMMPS data format

formats.cassandramcf.write_mcf

Write to Cassandra MCF format

formats.json_formats.compound_to_json

Write to a json file

Notes

When saving the compound as a json, only the following arguments are used: * filename * show_ports

spin(theta, around)[source]

Rotate Compound in place around an arbitrary vector.

Parameters
thetafloat

The angle by which to rotate the Compound, in radians.

aroundnp.ndarray, shape=(3,), dtype=float

The axis about which to spin the Compound.

successors()[source]

Yield Compounds below self in the hierarchy.

Yields
mb.Compound

The next Particle below self in the hierarchy

to_gmso()[source]

Create a GMSO Topology from a mBuild Compound.

Parameters
compoundmb.Compound

The mb.Compound to be converted.

Returns
topologygmso.Topology

The converted gmso Topology

to_intermol(molecule_types=None)[source]

Create an InterMol system from a Compound.

Parameters
molecule_typeslist or tuple of subclasses of Compound
Returns
intermol_systemintermol.system.System

See also

mbuild.conversion.to_intermol
to_networkx(names_only=False)[source]

Create a NetworkX graph representing the hierarchy of a Compound.

Parameters
names_onlybool, optional, default=False

Store only the names of the compounds in the graph, appended with their IDs, for distinction even if they have the same name. When set to False, the default behavior, the nodes are the compounds themselves.

Returns
Gnetworkx.DiGraph

See also

mbuild.conversion.to_networkx
mbuild.bond_graph

Notes

This digraph is not the bondgraph of the compound.

to_parmed(box=None, title='', residues=None, show_ports=False, infer_residues=False)[source]

Create a ParmEd Structure from a Compound.

Parameters
boxmb.Box, optional, default=self.boundingbox (with buffer)

Box information to be used when converting to a Structure. If ‘None’, self.box is used. If self.box is None, a bounding box is used with 0.5 nm buffer in each dimension to avoid overlapping atoms.

titlestr, optional, default=self.name

Title/name of the ParmEd Structure

residuesstr of list of str

Labels of residues in the Compound. Residues are assigned by checking against Compound.name.

show_portsboolean, optional, default=False

Include all port atoms when converting to a Structure.

infer_residuesbool, optional, default=False

Attempt to assign residues based on names of children.

Returns
parmed.structure.Structure

ParmEd Structure object converted from self

See also

mbuild.conversion.to_parmed
parmed.structure.Structure

Details on the ParmEd Structure object

to_pybel(box=None, title='', residues=None, show_ports=False, infer_residues=False)[source]

Create a pybel.Molecule from a Compound.

Parameters
boxmb.Box, def None
titlestr, optional, default=self.name

Title/name of the ParmEd Structure

residuesstr of list of str

Labels of residues in the Compound. Residues are assigned by checking against Compound.name.

show_portsboolean, optional, default=False

Include all port atoms when converting to a Structure.

infer_residuesbool, optional, default=False

Attempt to assign residues based on names of children

Returns
pybel.Molecule

See also

mbuild.conversion.to_pybel

Notes

Most of the mb.Compound is first converted to openbabel.OBMol And then pybel creates a pybel.Molecule from the OBMol Bond orders are assumed to be 1 OBMol atom indexing starts at 1, with spatial dimension Angstrom

to_trajectory(show_ports=False, chains=None, residues=None, box=None)[source]

Convert to an md.Trajectory and flatten the compound.

Parameters
show_portsbool, optional, default=False

Include all port atoms when converting to trajectory.

chainsmb.Compound or list of mb.Compound

Chain types to add to the topology

residuesstr of list of str

Labels of residues in the Compound. Residues are assigned by checking against Compound.name.

boxmb.Box, optional, default=self.boundingbox (with buffer)

Box information to be used when converting to a Trajectory. If ‘None’, self.box is used. If self.box is None, a bounding box is used with a 0.5 nm buffer in each dimension to avoid overlapping atoms.

Returns
trajectorymd.Trajectory

See also

_to_topology
translate(by)[source]

Translate the Compound by a vector.

Parameters
bynp.ndarray, shape=(3,), dtype=float
translate_to(pos)[source]

Translate the Compound to a specific position.

Parameters
posnp.ndarray, shape=3(,), dtype=float
unlabel_rigid_bodies()[source]

Remove all rigid body labels from the Compound.

update_coordinates(filename, update_port_locations=True)[source]

Update the coordinates of this Compound from a file.

Parameters
filenamestr

Name of file from which to load coordinates. Supported file types are the same as those supported by load()

update_port_locationsbool, optional, default=True

Update the locations of Ports so that they are shifted along with their anchor particles. Note: This conserves the location of Ports with respect to the anchor Particle, but does not conserve the orientation of Ports with respect to the molecule as a whole.

See also

load

Load coordinates from a file

visualize(show_ports=False, backend='py3dmol', color_scheme={})[source]

Visualize the Compound using py3dmol (default) or nglview.

Allows for visualization of a Compound within a Jupyter Notebook.

Parameters
show_portsbool, optional, default=False

Visualize Ports in addition to Particles

backendstr, optional, default=’py3dmol’

Specify the backend package to visualize compounds Currently supported: py3dmol, nglview

color_schemedict, optional

Specify coloring for non-elemental particles keys are strings of the particle names values are strings of the colors i.e. {‘_CGBEAD’: ‘blue’}

property xyz

Return all particle coordinates in this compound.

Returns
posnp.ndarray, shape=(n, 3), dtype=float

Array with the positions of all particles.

property xyz_with_ports

Return all particle coordinates in this compound including ports.

Returns
posnp.ndarray, shape=(n, 3), dtype=float

Array with the positions of all particles and ports.

Box

class mbuild.Box(lengths, angles=None, precision=None)[source]

A box representing the bounds of the system.

Parameters
lengthslist-like, shape=(3,), dtype=float

Lengths of the edges of the box.

angleslist-like, shape=(3,), dtype=float, default=None

Angles (in degrees) that define the tilt of the edges of the box. If None is given, angles are assumed to be [90.0, 90.0, 90.0]. These are also known as alpha, beta, gamma in the crystallography community.

precisionint, optional, default=None

Control the precision of the floating point representation of box attributes. If none provided, the default is 6 decimals.

Notes

Box vectors are expected to be provided in row-major format.

Attributes
vectorsnp.ndarray, shape=(3,3), dtype=float

Box representation as a 3x3 matrix.

lengthstuple, shape=(3,), dtype=float

Lengths of the box.

anglestuple, shape=(3,), dtype=float

Angles defining the tilt of the box (alpha, beta, gamma).

Lxfloat

Length in the x direction.

Lyfloat

Length in the y direction.

Lzfloat

Length in the z direction.

xyfloat

Tilt factor xy of the box.

xzfloat

Tilt factor xz of the box.

yzfloat

Tilt factor yz of the box.

precisionint

Amount of decimals to represent floating point values.

property Lx

Length in the x direction.

property Ly

Length in the y direction.

property Lz

Length in the z direction.

property angles

Angles defining the tilt of the box (alpha, beta, gamma).

property box_parameters

Lengths and tilt factors of the box.

property bravais_parameters

Return the Box representation as Bravais lattice parameters.

Based on the box vectors, return the parameters to describe the box in terms of the Bravais lattice parameters:

a,b,c = the edges of the Box alpha, beta, gamma = angles(tilt) of the parallelepiped, in degrees

Returns
parameterstuple of floats,

(a, b, c, alpha, beta, gamma)

classmethod from_lengths_angles(lengths, angles, precision=None)[source]

Generate a box from lengths and angles.

classmethod from_lengths_tilt_factors(lengths, tilt_factors=None, precision=None)[source]

Generate a box from box lengths and tilt factors.

classmethod from_lo_hi_tilt_factors(lo, hi, tilt_factors, precision=None)[source]

Generate a box from a lo, hi convention and tilt factors.

classmethod from_mins_maxs_angles(mins, maxs, angles, precision=None)[source]

Generate a box from min/max distance calculations and angles.

classmethod from_uvec_lengths(uvec, lengths, precision=None)[source]

Generate a box from unit vectors and lengths.

classmethod from_vectors(vectors, precision=None)[source]

Generate a box from box vectors.

property lengths

Lengths of the box.

property precision

Amount of decimals to represent floating point values.

property tilt_factors

Return the 3 tilt_factors (xy, xz, yz) of the box.

property vectors

Box representation as a 3x3 matrix.

property xy

Tilt factor xy of the box.

property xz

Tilt factor xz of the box.

property yz

Tilt factor yz of the box.

Lattice

class mbuild.Lattice(lattice_spacing=None, lattice_vectors=None, lattice_points=None, angles=None)[source]

Develop crystal structure from user defined inputs.

Lattice, the abstract building block of a crystal cell. Once defined by the user, the lattice can then be populated with Compounds and replicated as many cell lengths desired in 3D space.

A Lattice is defined through the Bravais lattice definitions. With edge vectors a1, a2, a3; lattice spacing a,b,c; and lattice points at unique fractional positions between 0-1 in 3 dimensions. This encapsulates distance, area, volume, depending on the parameters defined.

Parameters
lattice_spacingarray-like, shape=(3,), required, dtype=float

Array of lattice spacings a,b,c for the cell.

lattice_vectorsarray-like, shape=(3, 3), optional, default=None,

Vectors that encase the unit cell corresponding to dimension. Will only default to these values if no angles were defined as well. If None is given, assumes an identity matrix [[1,0,0], [0,1,0], [0,0,1]]

lattice_pointsdictionary, shape={‘id’: [[nested list of positions]]

optional, default={‘default’: [[0.,0.,0.]]} Locations of all lattice points in cell using fractional coordinates.

anglesarray-like, shape=(3,), optional, dtype=float

Array of inter-planar Bravais angles in degrees.

Examples

Generating a triclinic lattice for cholesterol.

>>> import mbuild as mb
>>> from mbuild.utils.io import get_fn
>>> # reading in the lattice parameters for crystalline cholesterol
>>> angle_values = [94.64, 90.67, 96.32]
>>> spacing = [1.4172, 3.4209, 1.0481]
>>> basis = {'cholesterol':[[0., 0., 0.]]}
>>> cholesterol_lattice = mb.Lattice(spacing,
...                                  angles=angle_values,
...                                  lattice_points=basis)
>>>
>>> # The lattice based on the bravais lattice parameters of crystalline
>>> # cholesterol was generated.
>>>
>>> # Replicating the triclinic unit cell out 3 replications
>>> # in x,y,z directions.
>>>
>>> cholesterol_unit = mb.Compound()
>>> cholesterol_unit = mb.load(get_fn('cholesterol.pdb'))
>>> # associate basis vector with id 'cholesterol' to cholesterol Compound
>>> basis_dictionary = {'cholesterol' : cholesterol_unit}
>>> expanded_cell = cholesterol_lattice.populate(x=3, y=3, z=3,
...                              compound_dict=basis_dictionary)

The unit cell of cholesterol was associated with a Compound that contains the connectivity data and spatial arrangements of a cholesterol molecule. The unit cell was then expanded out in x,y,z directions and cholesterol Compounds were populated.

Generating BCC CsCl crystal structure

>>> import mbuild as mb
>>> chlorine = mb.Compound(name='Cl')
>>> # angles not needed, when not provided, defaults to 90,90,90
>>> cesium = mb.Compound(name='Cs')
>>> spacing = [.4123, .4123, .4123]
>>> basis = {'Cl' : [[0., 0., 0.]], 'Cs' : [[.5, .5, .5]]}
>>> cscl_lattice = mb.Lattice(spacing, lattice_points=basis)
>>>
>>> # Now associate id with Compounds for lattice points and replicate 3x
>>>
>>> cscl_dict = {'Cl' : chlorine, 'Cs' : cesium}
>>> cscl_compound = cscl_lattice.populate(x=3, y=3, z=3,
...                                       compound_dict=cscl_dict)

A multi-Compound basis was created and replicated. For each unique basis atom position, a separate entry must be completed for the basis_atom input.

Generating FCC Copper cell with lattice_vectors instead of angles

>>> import mbuild as mb
>>> copper = mb.Compound(name='Cu')
>>> lattice_vector = [[1, 0, 0], [0, 1, 0], [0, 0, 1]]
>>> spacing = [.36149, .36149, .36149]
>>> copper_locations = [[0., 0., 0.], [.5, .5, 0.],
...                     [.5, 0., .5], [0., .5, .5]]
>>> basis = {'Cu' : copper_locations}
>>> copper_lattice = mb.Lattice(lattice_spacing = spacing,
...                             lattice_vectors=lattice_vector,
...                             lattice_points=basis)
>>> copper_dict = {'Cu' : copper}
>>> copper_pillar = copper_lattice.populate(x=3, y=3, z=20,
...                                       compound_dict=copper_dict)

Generating the 2d Structure Graphene carbon backbone

>>> import mbuild as mb
>>> carbon = mb.Compound(name='C')
>>> angles = [90, 90, 120]
>>> carbon_locations = [[0, 0, 0], [2/3, 1/3, 0]]
>>> basis = {'C' : carbon_locations}
>>> graphene = mb.Lattice(lattice_spacing=[.2456, .2456, 0],
...                        angles=angles, lattice_points=basis)
>>> carbon_dict = {'C' : carbon}
>>> graphene_cell = graphene.populate(compound_dict=carbon_dict,
...                                   x=3, y=3, z=1)
Attributes
dimensionint, 3

Default dimensionality within mBuild. If choosing a lower dimension, pad the relevant arrays with zeroes.

lattice_spacingnumpy array, shape=(3,), required, dtype=float

Array of lattice spacings a,b,c for the cell.

lattice_vectorsnumpy array, shape=(3, 3), optional

default=[[1,0,0], [0,1,0], [0,0,1]] Vectors that encase the unit cell corresponding to dimension. Will only default to these values if no angles were defined as well.

lattice_pointsdictionary, shape={‘id’: [[nested list of positions]]

optional, default={‘default’: [[0.,0.,0.]]} Locations of all lattice points in cell using fractional coordinates.

anglesnumpy array, shape=(3,), optional, dtype=float

Array of inter-planar Bravais angles

populate(compound_dict=None, x=1, y=1, z=1)[source]

Expand lattice and create compound from lattice.

Expands lattice based on user input. The user must also pass in a dictionary that contains the keys that exist in the basis_dict. The corresponding Compound will be the full lattice returned to the user.

If no dictionary is passed to the user, Dummy Compounds will be used.

Parameters
xint, optional, default=1

How many iterations in the x direction.

yint, optional, default=1

How many iterations in the y direction.

zint, optional, default=1

How many iterations in the z direction.

compound_dictdictionary, optional, default=None

Link between basis_dict and Compounds.

Raises
ValueError

incorrect x,y, or z values.

TypeError

incorrect type for basis vector

Notes

Called after constructor by user.

Port

class mbuild.Port(anchor=None, orientation=None, separation=0)[source]

A set of four ghost Particles used to connect parts.

Parameters
anchormb.Particle, optional, default=None

A Particle associated with the port. Used to form bonds.

orientationarray-like, shape=(3,), optional, default=[0, 1, 0]

Vector along which to orient the port

separationfloat, optional, default=0

Distance to shift port along the orientation vector from the anchor particle position. If no anchor is provided, the port will be shifted from the origin.

Attributes
anchormb.Particle, optional, default=None

A Particle associated with the port. Used to form bonds.

upmb.Compound

Collection of 4 ghost particles used to perform equivalence transforms. Faces the opposite direction as self[‘down’].

downmb.Compound

Collection of 4 ghost particles used to perform equivalence transforms. Faces the opposite direction as self[‘up’].

usedbool

Status of whether a port has been occupied following an equivalence transform.

property access_labels

List labels used to access the Port.

Returns
list of str

Strings that can be used to access this Port relative to self.root

property center

Get the cartesian center of the Port.

property direction

Get the unit vector pointing in the ‘direction’ of the Port.

property separation

Get the distance between a port and its anchor particle.

If the port has no anchor particle, returns None.

update_orientation(orientation)[source]

Change the direction between a port and its anchor particle.

orientationarray-like, shape=(3,), required

Vector along which to orient the port

update_separation(separation)[source]

Change the distance between a port and its anchor particle.

separationfloat, required

Distance to shift port along the orientation vector from the anchor particle position. If no anchor is provided, the port will be shifted from the origin.