%load_ext autoreload
%autoreload 2

import molsysmt as msm

Get bondgraph

Having the topology of the molecular system given by the covalent bonds can be very useful to, for instance, getting the number of disconnected components we have in the system. The method molsysmt.bondgraph returns the undirected network or graph where the nodes are the atoms and the undirected edges are given by the covalent bonds. The resulting graph can be encoded as a native python object from libraries such as NetworkX.

Lets load a molecular system to show how molsysmt.bondgraph works:

molecular_system = msm.systems['TcTIM']['1tcd.h5msm']
molecular_system = msm.convert(molecular_system, 'molsysmt.MolSys')

Lets get the bond graph defined by the atoms of the first molecule:

graph = msm.topology.get_bondgraph(molecular_system, selection='molecule_index==0', to_form="networkx.Graph")

n_nodes = graph.number_of_nodes()
n_edges = graph.number_of_edges()
print('The graph has {} nodes and {} edges'.format(n_nodes, n_edges))

The graph has 1906 nodes and 1942 edges

n_atoms, n_bonds = msm.get(molecular_system, element='atom', selection='molecule_index==0', n_atoms=True,
                           n_inner_bonds=True)

print('The molecule has {} atoms and {} bonds'.format(n_atoms, n_bonds))

The molecule has 1906 atoms and 1942 bonds

The output object was chosen to be a NetworkX graph. Lets see how many disconnected components are in the molecule:

from networkx import connected_components

components = connected_components(graph)

print('Number of components in the graph: {}'.format(len(list(components))))

Number of components in the graph: 1

Actually, this is the definition of the element “component” in MolSysMT:

msm.get(molecular_system, element='atom', selection='molecule_index==0', n_components=True)

1