The MultiMolecule Class¶
The API of the MultiMolecule
class.
API FOX.MultiMolecule¶
-
class
FOX.classes.multi_mol.
MultiMolecule
(coords: numpy.ndarray, atoms: Optional[Dict[str, List[int]]] = None, bonds: Optional[numpy.ndarray] = None, properties: Optional[Dict[str, Any]] = None)[source]¶ A class designed for handling a and manipulating large numbers of molecules.
More specifically, different conformations of a single molecule as derived from, for example, an intrinsic reaction coordinate calculation (IRC) or a molecular dymanics trajectory (MD). The class has access to four attributes (further details are provided under parameters):
Parameters: - coords (\(m*n*3\) np.ndarray [np.float64]) – A 3D array with the cartesian coordinates of \(m\) molecules with \(n\) atoms.
- atoms (dict [str, list [int]]) – A dictionary with atomic symbols as keys and matching atomic indices as values.
Stored in the
MultiMolecule.atoms
attribute. - bonds (\(k*3\) np.ndarray [np.int64]) – A 2D array with indices of the atoms defining all \(k\) bonds
(columns 1 & 2) and their respective bond orders multiplied by 10 (column 3).
Stored in the
MultiMolecule.bonds
attribute. - properties (dict) – A Settings instance for storing miscellaneous user-defined (meta-)data.
Is devoid of keys by default.
Stored in the
MultiMolecule.properties
attribute.
-
atoms
¶ A dictionary with atomic symbols as keys and matching atomic indices as values.
Type: dict [str, list [int]]
-
bonds
¶ A 2D array with indices of the atoms defining all \(k\) bonds (columns 1 & 2) and their respective bond orders multiplied by 10 (column 3).
Type: \(k*3\) np.ndarray [np.int64]
-
properties
¶ A Settings instance for storing miscellaneous user-defined (meta-)data. Is devoid of keys by default.
Type: plams.Settings
-
round
(decimals=0, inplace=True)[source]¶ Round the Cartesian coordinates of this instance to a given number of decimals.
Parameters: Return type:
-
delete_atoms
(atom_subset)[source]¶ Create a copy of this instance with all atoms in atom_subset removed.
Parameters: atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as determined by their atomic index or atomic symbol. Include all \(n\) atoms per molecule in this instance if None
.Returns: A new MultiMolecule
instance with all atoms in atom_subset removed.Return type: FOX.MultiMolecule Raises: TypeError – Raised if atom_subset is None
.
-
add_atoms
(coords, symbols='Xx')[source]¶ Create a copy of this instance with all atoms in atom_subset appended.
Examples
>>> import numpy as np >>> from FOX import MultiMolecule, example_xyz >>> mol = MultiMolecule.from_xyz(example_xyz) >>> coords: np.ndarray = np.random.rand(73, 3) # Add 73 new atoms with random coords >>> symbols = 'Br' >>> mol_new: MultiMolecule = mol.add_atoms(coords, symbols) >>> print(repr(mol)) MultiMolecule(..., shape=(4905, 227, 3), dtype='float64') >>> print(repr(mol_new)) MultiMolecule(..., shape=(4905, 300, 3), dtype='float64')
Parameters: Returns: A new
MultiMolecule
instance with all atoms in atom_subset appended.Return type:
-
guess_bonds
(atom_subset=None)[source]¶ Guess bonds within the molecules based on atom type and inter-atomic distances.
Bonds are guessed based on the first molecule in this instance Performs an inplace modification of self.bonds
Parameters: atom_subset (Sequence) – A tuple of atomic symbols. Bonds are guessed between all atoms whose atomic symbol is in atom_subset. If None
, guess bonds for all atoms in this instance.Return type: None
-
random_slice
(start=0, stop=None, p=0.5, inplace=False)[source]¶ Construct a new
MultiMolecule
instance by randomly slicing this instance.The probability of including a particular element is equivalent to p.
Parameters: - start (int) – Start of the interval.
- stop (int) – End of the interval.
- p (float) – The probability of including each particular molecule in this instance.
Values must be between
0.0
(0%) and1.0
(100%). - inplace (bool) – Instead of returning the new coordinates, perform an inplace update of this instance.
Returns: If inplace is
True
, return a newMultiMolecule
instance.Return type: Raises: ValueError – Raised if p is smaller than
0.0
or larger than1.0
.
-
reset_origin
(mol_subset=None, atom_subset=None, inplace=True)[source]¶ Reallign all molecules in this instance.
All molecules in this instance are rotating and translating, by performing a partial partial Procrustes superimposition with respect to the first molecule in this instance.
The superimposition is carried out with respect to the first molecule in this instance.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - inplace (bool) – Instead of returning the new coordinates, perform an inplace update of this instance.
Returns: If inplace is
True
, return a newMultiMolecule
instance.Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
sort
(sort_by='symbol', reverse=False, inplace=True)[source]¶ Sort the atoms in this instance and self.atoms, performing in inplace update.
Parameters: - sort_by (str or Sequence [int]) – The property which is to be used for sorting.
Accepted values:
"symbol"
(i.e. alphabetical),"atnum"
,"mass"
,"radius"
or"connectors"
. See the plams.PeriodicTable module for more details. Alternatively, a user-specified sequence of indices can be provided for sorting. - reverse (bool) – Sort in reversed order.
- inplace (bool) – Instead of returning the new coordinates, perform an inplace update of this instance.
Returns: If inplace is
True
, return a newMultiMolecule
instance.Return type: - sort_by (str or Sequence [int]) – The property which is to be used for sorting.
Accepted values:
-
residue_argsort
(concatenate=True)[source]¶ Return the indices that would sort this instance by residue number.
Residues are defined based on moleculair fragments based on self.bonds.
Parameters: concatenate (bool) – If False
, returned a nested list with atomic indices. Each sublist contains the indices of a single residue.Returns: A 1D array of indices that would sort \(n\) atoms this instance. Return type: \(n\) np.ndarray [np.int64]
-
get_center_of_mass
(mol_subset=None, atom_subset=None)[source]¶ Get the center of mass.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A 2D array with the centres of mass of \(m\) molecules with \(n\) atoms.
Return type: \(m*3\) np.ndarray [np.float64]
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_bonds_per_atom
(atom_subset=None)[source]¶ Get the number of bonds per atom in this instance.
Parameters: atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as determined by their atomic index or atomic symbol. Include all \(n\) atoms per molecule in this instance if None
.Returns: A 1D array with the number of bonds per atom, for all \(n\) atoms in this instance. Return type: \(n\) np.ndarray [np.int64]
-
init_average_velocity
(timestep=1.0, rms=False, mol_subset=None, atom_subset=None)[source]¶ Calculate the average atomic velocty.
The average velocity (in fs/A) is calculated for all atoms in atom_subset over the course of a trajectory.
The velocity is averaged over all atoms in a particular atom subset.
Parameters: - timestep (float) – The stepsize, in femtoseconds, between subsequent frames.
- rms (bool) – Calculate the root-mean squared average velocity instead.
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A dataframe holding \(m-1\) velocities averaged over one or more atom subsets.
Return type:
-
init_time_averaged_velocity
(timestep=1.0, rms=False, mol_subset=None, atom_subset=None)[source]¶ Calculate the time-averaged velocty.
The time-averaged velocity (in fs/A) is calculated for all atoms in atom_subset over the course of a trajectory.
Parameters: - timestep (float) – The stepsize, in femtoseconds, between subsequent frames.
- rms (bool) – Calculate the root-mean squared time-averaged velocity instead.
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A dataframe holding \(m-1\) time-averaged velocities.
Return type:
-
init_rmsd
(mol_subset=None, atom_subset=None, reset_origin=True)[source]¶ Initialize the RMSD calculation, returning a dataframe.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - reset_origin (bool) – Reset the origin of each molecule in this instance by means of a partial Procrustes superimposition, translating and rotating the molecules.
Returns: A dataframe of RMSDs with one column for every string or list of ints in atom_subset. Keys consist of atomic symbols (e.g.
"Cd"
) if atom_subset contains strings, otherwise a more generic ‘series ‘ + str(int) scheme is adopted (e.g."series 2"
). Molecular indices are used as index.Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
init_rmsf
(mol_subset=None, atom_subset=None, reset_origin=True)[source]¶ Initialize the RMSF calculation, returning a dataframe.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - reset_origin (bool) – Reset the origin of each molecule in this instance by means of a partial Procrustes superimposition, translating and rotating the molecules.
Returns: A dataframe of RMSFs with one column for every string or list of ints in atom_subset. Keys consist of atomic symbols (e.g.
"Cd"
) if atom_subset contains strings, otherwise a more generic ‘series ‘ + str(int) scheme is adopted (e.g."series 2"
). Molecular indices are used as indices.Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_average_velocity
(timestep=1.0, rms=False, mol_subset=None, atom_subset=None)[source]¶ Return the mean or root-mean squared velocity.
Parameters: - timestep (float) – The stepsize, in femtoseconds, between subsequent frames.
- rms (bool) – Calculate the root-mean squared average velocity instead.
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A 1D array holding \(m-1\) velocities averaged over one or more atom subsets.
Return type: \(m-1\) np.ndarray [np.float64]
-
get_time_averaged_velocity
(timestep=1.0, rms=False, mol_subset=None, atom_subset=None)[source]¶ Return the mean or root-mean squared velocity (mean = time-averaged).
Parameters: - timestep (float) – The stepsize, in femtoseconds, between subsequent frames.
- rms (bool) – Calculate the root-mean squared average velocity instead.
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A 1D array holding \(n\) time-averaged velocities.
Return type: \(n\) np.ndarray [np.float64]
-
get_velocity
(timestep=1.0, norm=True, mol_subset=None, atom_subset=None)[source]¶ Return the atomic velocties.
The velocity (in fs/A) is calculated for all atoms in atom_subset over the course of a trajectory.
Parameters: - timestep (float) – The stepsize, in femtoseconds, between subsequent frames.
- norm (bool) – If
True
return the norm of the \(x\), \(y\) and \(z\) velocity components. - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A 2D or 3D array of atomic velocities, the number of dimensions depending on the value of norm (
True
= 2D;False
= 3D).Return type: \(m*n\) or \(m*n*3\) np.ndarray [np.float64]
-
get_rmsd
(mol_subset=None, atom_subset=None)[source]¶ Calculate the root mean square displacement (RMSD).
The RMSD is calculated with respect to the first molecule in this instance.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A dataframe with the RMSD as a function of the XYZ frame numbers.
Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_rmsf
(mol_subset=None, atom_subset=None)[source]¶ Calculate the root mean square fluctuation (RMSF).
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A dataframe with the RMSF as a function of atomic indices.
Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
init_shell_search
(mol_subset=None, atom_subset=None, rdf_cutoff=0.5)[source]¶ Calculate and return properties which can help determining shell structures.
The following two properties are calculated and returned:
- The mean distance (per atom) with respect to the center of mass (i.e. a modified RMSF).
- A series mapping abritrary atomic indices in the RMSF to the actual atomic indices.
- The radial distribution function (RDF) with respect to the center of mass.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - rdf_cutoff (float) – Remove all values in the RDF below this value (Angstrom). Usefull for dealing with divergence as the “inter-atomic” distance approaches 0.0 A.
Returns: - Returns the following items:
- A dataframe holding the mean distance of all atoms with respect the to center of mass.
- A series mapping the indices from 1. to the actual atomic indices.
- A dataframe holding the RDF with respect to the center of mass.
Return type:
-
static
get_at_idx
(rmsf, idx_series, dist_dict)[source]¶ Create subsets of atomic indices.
The subset is created (using rmsf and idx_series) based on distance criteria in dist_dict.
For example,
dist_dict = {'Cd': [3.0, 6.5]}
will create and return a dictionary with three keys: One for all atoms whose RMSF is smaller than 3.0, one where the RMSF is between 3.0 and 6.5, and finally one where the RMSF is larger than 6.5.Examples
>>> dist_dict = {'Cd': [3.0, 6.5]} >>> idx_series = pd.Series(np.arange(12)) >>> rmsf = pd.DataFrame({'Cd': np.arange(12, dtype=float)}) >>> get_at_idx(rmsf, idx_series, dist_dict) {'Cd_1': [0, 1, 2], 'Cd_2': [3, 4, 5], 'Cd_3': [7, 8, 9, 10, 11] }
Parameters: Returns: A dictionary with atomic symbols as keys, and matching atomic indices as values.
Return type: Raises: KeyError – Raised if a key in dist_dict is absent from rmsf.
-
init_rdf
(mol_subset=None, atom_subset=None, dr=0.05, r_max=12.0, mem_level=2)[source]¶ Initialize the calculation of radial distribution functions (RDFs).
RDFs are calculated for all possible atom-pairs in atom_subset and returned as a dataframe.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - dr (float) – The integration step-size in Ångström, i.e. the distance between concentric spheres.
- r_max (float) – The maximum to be evaluated interatomic distance in Ångström.
- mem_level (int) –
Set the level of to-be consumed memory and, by extension, the execution speed. Given a molecule subset of size \(m\), atom subsets of (up to) size \(n\) and the resulting RDF with \(p\) points (
p = r_max / dr
), the mem_level values can be interpreted as following:0
: Slow; memory scaling: \(n^2\)1
: Medium; memory scaling: \(n^2 + m * p\)2
: Fast; memory scaling: \(n^2 * m\)
Returns: A dataframe of radial distribution functions, averaged over all conformations in xyz_array. Keys are of the form: at_symbol1 + ‘ ‘ + at_symbol2 (e.g.
"Cd Cd"
). Radii are used as index.Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_dist_mat
(mol_subset=None, atom_subset=(None, None))[source]¶ Create and return a distance matrix for all molecules and atoms in this instance.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A 3D distance matrix of \(m\) molecules, created out of two sets of \(n\) and \(k\) atoms.
Return type: \(m*n*k\) np.ndarray [np.float64]
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_pair_dict
(atom_subset, r=2)[source]¶ Take a subset of atoms and return a dictionary.
Parameters: Return type:
-
init_power_spectrum
(mol_subset=None, atom_subset=None, freq_max=4000)[source]¶ Calculate and return the power spectrum associated with this instance.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - freq_max (int) – The maximum to be returned wavenumber (cm**-1).
Returns: A DataFrame containing the power spectrum for each set of atoms in atom_subset.
Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_vacf
(mol_subset=None, atom_subset=None)[source]¶ Calculate and return the velocity autocorrelation function (VACF).
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A DataFrame containing the power spectrum for each set of atoms in atom_subset.
Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
init_adf
(mol_subset=None, atom_subset=None, r_max=8.0, weight=<function neg_exp>)[source]¶ Initialize the calculation of distance-weighted angular distribution functions (ADFs).
ADFs are calculated for all possible atom-pairs in atom_subset and returned as a dataframe.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - r_max (float or str) – The maximum inter-atomic distance (in Angstrom) for which angles are constructed.
The distance cuttoff can be disabled by settings this value to
np.inf
,"np.inf"
or"inf"
. - weight (Callable[[np.ndarray], np.ndarray], optional) – A callable for creating a weighting factor from inter-atomic distances.
The callable should take an array as input and return an array.
Given an angle \(\phi_{ijk}\), to the distance \(r_{ijk}\) is defined
as \(max[r_{ij}, r_{jk}]\).
Set to
None
to disable distance weighting.
Returns: A dataframe of angular distribution functions, averaged over all conformations in this instance.
Return type: Note
Disabling the distance cuttoff is strongly recommended (i.e. it is faster) for large values of r_max. As a rough guideline,
r_max="inf"
is roughly as fast asr_max=15.0
(though this is, of course, system dependant).Note
The ADF construction will be conducted in parralel if the DASK package is installed. DASK can be installed, via anaconda, with the following command:
conda install -n FOX -y -c conda-forge dask
.- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
get_angle_mat
(mol_subset=0, atom_subset=(None, None, None), get_r_max=False)[source]¶ Create and return an angle matrix for all molecules and atoms in this instance.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - get_r_max (bool) – Whether or not the maximum distance should be returned or not.
Returns: A 4D angle matrix of \(m\) molecules, created out of three sets of \(n\), \(k\) and \(l\) atoms. If get_r_max =
True
, also return the maximum distance.Return type: \(m*n*k*l\) np.ndarray [np.float64] and (optionally) float
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
as_pdb
(filename, mol_subset=0)[source]¶ Convert a MultiMolecule object into one or more Protein DataBank files (.pdb).
Utilizes the plams.Molecule.write method.
Parameters: Return type:
-
as_mol2
(filename, mol_subset=0)[source]¶ Convert a MultiMolecule object into one or more .mol2 files.
Utilizes the plams.Molecule.write method.
Parameters: Return type:
-
as_mol
(filename, mol_subset=0)[source]¶ Convert a MultiMolecule object into one or more .mol files.
Utilizes the plams.Molecule.write method.
Parameters: Return type:
-
as_xyz
(filename, mol_subset=None)[source]¶ Create an .xyz file out of this instance.
Comments will be constructed by iteration through
MultiMolecule.properties["comments"]
if the following two conditions are fulfilled:- The
"comments"
key is actually present inMultiMolecule.properties
. MultiMolecule.properties["comments"]
is an iterable.
Parameters: Return type: - The
-
as_mass_weighted
(mol_subset=None, atom_subset=None, inplace=False)[source]¶ Transform the Cartesian of this instance into mass-weighted Cartesian coordinates.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
. - inplace (bool) – Instead of returning the new coordinates, perform an inplace update of this instance.
Returns: if inplace =
False
return a newMultiMolecule
instance with the mass-weighted Cartesian coordinates of \(m\) molecules with \(n\) atoms.Return type: \(m*n*3\) np.ndarray [np.float64] or None
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
from_mass_weighted
(mol_subset=None, atom_subset=None)[source]¶ Transform this instance from mass-weighted Cartesian into Cartesian coordinates.
Performs an inplace update of this instance.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Return type: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
as_Molecule
(mol_subset=None, atom_subset=None)[source]¶ Convert this instance into a list of plams.Molecule.
Parameters: - mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
None
. - atom_subset (Sequence) – Perform the calculation on a subset of atoms in this instance, as
determined by their atomic index or atomic symbol.
Include all \(n\) atoms per molecule in this instance if
None
.
Returns: A list of \(m\) PLAMS molecules constructed from this instance.
Return type: \(m\) list [plams.Molecule]
- mol_subset (slice) – Perform the calculation on a subset of molecules in this instance, as
determined by their moleculair index.
Include all \(m\) molecules in this instance if
-
classmethod
from_Molecule
(mol_list, subset='atoms')[source]¶ Construct a
MultiMolecule
instance from one or more PLAMS molecules.Parameters: - mol_list (plams.Molecule or list [plams.Molecule]) – A PLAMS molecule or list of PLAMS molecules.
- subset (Sequence [str]) – Transfer a subset of plams.Molecule attributes to this instance.
If
None
, transfer all attributes. Accepts one or more of the following values as strings:"properties"
,"atoms"
and/or"bonds"
.
Returns: A
MultiMolecule
instance constructed from mol_list.Return type:
-
classmethod
from_xyz
(filename, bonds=None, properties=None)[source]¶ Construct a
MultiMolecule
instance from a (multi) .xyz file.Comment lines extracted from the .xyz file are stored, as array, under
MultiMolecule.properties["comments"]
.Parameters: - filename (str) – The path+filename of an .xyz file.
- bonds (\(k*3\) np.ndarray [np.int64]) – An optional 2D array with indices of the atoms defining all \(k\) bonds (columns 1 & 2) and their respective bond orders multiplied by 10 (column 3). Stored in the MultieMolecule.bonds attribute.
- properties (dict) – A Settings object (subclass of dictionary) intended for storing miscellaneous user-defined (meta-)data. Is devoid of keys by default. Stored in the MultiMolecule.properties attribute.
Returns: A
MultiMolecule
instance constructed from filename.Return type:
-
classmethod
from_kf
(filename, bonds=None, properties=None)[source]¶ Construct a
MultiMolecule
instance from a KF binary file.Parameters: - filename (str) – The path+filename of an KF binary file.
- bonds (\(k*3\) np.ndarray [np.int64]) – An optional 2D array with indices of the atoms defining all \(k\) bonds (columns 1 & 2) and their respective bond orders multiplied by 10 (column 3). Stored in the MultieMolecule.bonds attribute.
- properties (dict) – A Settings object (subclass of dictionary) intended for storing miscellaneous user-defined (meta-)data. Is devoid of keys by default. Stored in the MultiMolecule.properties attribute.
Returns: A
MultiMolecule
instance constructed from filename.Return type:
API FOX._MultiMolecule¶
-
class
FOX.classes.multi_mol_magic.
_MultiMolecule
(coords: numpy.ndarray, atoms: Optional[Dict[str, List[int]]] = None, bonds: Optional[numpy.ndarray] = None, properties: Optional[Dict[str, Any]] = None)[source]¶ Private superclass of
MultiMolecule
.Handles all magic methods and @property decorated methods.
-
property
loc
¶ A getter and setter for atom-type-based slicing.
Get, set and del operations are performed using the list(s) of atomic indices associated with the provided atomic symbol(s). Accepts either one or more atomic symbols.
Examples
>>> mol = MultiMolecule(...) >>> mol.atoms['Cd'] = [0, 1, 2, 3, 4, 5] >>> mol.atoms['Se'] = [6, 7, 8, 9, 10, 11] >>> mol.atoms['O'] = [12, 13, 14] >>> (mol.loc['Cd'] == mol[mol.atoms['Cd']]).all() True >>> idx = mol.atoms['Cd'] + mol.atoms['Se'] + mol.atoms['O'] >>> (mol.loc['Cd', 'Se', 'O'] == mol[idx]).all() True >>> mol.loc['Cd'] = 1 >>> print((mol.loc['Cd'] == 1).all()) True >>> del mol.loc['Cd'] ValueError: cannot delete array elements
Parameters: mol ( MultiMolecule
) – A MultiMolecule instance; seeAtGetter.atoms
.-
mol
¶ A MultiMolecule instance.
Type: MultiMolecule
Return type: LocGetter
-
-
property
atom12
¶ Get or set the indices of the atoms for all bonds in
MultiMolecule.bonds
as 2D array.Return type: _MultiMolecule
-
property
atom1
¶ Get or set the indices of the first atoms in all bonds of
MultiMolecule.bonds
as 1D array.Return type: _MultiMolecule
-
property
atom2
¶ Get or set the indices of the second atoms in all bonds of
MultiMolecule.bonds
as 1D array.Return type: ndarray
-
property
order
¶ Get or set the bond orders for all bonds in
MultiMolecule.bonds
as 1D array.Return type: ndarray
-
property
x
¶ Get or set the x coordinates for all atoms in instance as 2D array.
Return type: _MultiMolecule
-
property
y
¶ Get or set the y coordinates for all atoms in this instance as 2D array.
Return type: _MultiMolecule
-
property
z
¶ Get or set the z coordinates for all atoms in this instance as 2D array.
Return type: _MultiMolecule
-
property
symbol
¶ Get the atomic symbols of all atoms in
MultiMolecule.atoms
as 1D array.Return type: ndarray
-
property
atnum
¶ Get the atomic numbers of all atoms in
MultiMolecule.atoms
as 1D array.Return type: ndarray
-
property
mass
¶ Get the atomic masses of all atoms in
MultiMolecule.atoms
as 1D array.Return type: ndarray
-
property
radius
¶ Get the atomic radii of all atoms in
MultiMolecule.atoms
as 1d array.Return type: ndarray
-
property
connectors
¶ Get the atomic connectors of all atoms in
MultiMolecule.atoms
as 1D array.Return type: ndarray
-
copy
(order='C', deep=True)[source]¶ Create a copy of this instance.
Parameters: - order (str) – Controls the memory layout of the copy. See np.ndarray.copy for details.
- copy_attr (bool) – Whether or not the attributes of this instance should be returned as copies or views.
Returns: A copy of this instance.
Return type:
-
property