3

As per Should I ask my question about which data structure to use here? I'm asking this here. Hopefully this isn't too implementation specific.

I'm currently developing a program that will represent chemical structures (and eventually allow manipulations/reactions of them). I've hashed out two different systems of storing this data but I'm not sure how to pick which one.

Option 1: Nested matrices

Chemicals have (more or less) definite structures and relative positions of items, and it seems as though it would be easier to represent something like that with a matrix - then I can assign definite 'coordinates' to each atom or bond. Unfortunately, you end up with some overlap if you have a structure that branches out a fair bit. To solve this I ended up working out a solution that used nested, 5x5 matrices.

[[None,  None, [CH3], None, None ],
 [None,  None, Bond,  None, None ],
 [[CH3], Bond, C,     Bond, [CH3]],
 [None,  None, Bond,  None, None ],
 [None,  None, [CH3], None, None ]]

Here C represents the element carbon (this could be represented a couple different ways, and isn't particularly relevant here). Bond represents a bond, and [CH3] represents another 5x5 matrix. This is the uppermost ([0][2]) CH3

[[None, None, H,    None, None],
 [None, None, Bond, None, None],
 [H,    Bond, C,    Bond, H   ],
 [None, None, Bond, None, None],
 [None, None, Root, None, None]]

Here the lowermost Bond is a reference to the same object as the uppermost Bond in the first structure, while Root is (some undetermined) way of relating that this is where the matrix connects to the next level up. If this were a more complex structure, one of the Hs would be replaced by another matrix, such as an additional [CH3].

Advantages of this method:

  • Everything is explicitly laid out, including structure, which is good according to The Zen of Python.
  • When these chemicals are (eventually) represented graphically it will (theoretically) be much easier to make sure everything is in its place

Disadvantages:

  • This isn't readable at all if you have a more complex structure
  • This is nested, not flat
  • This is going to require a lot of recursion to walk through

Option 2: Dictionaries

Another way I could represent these with significantly less nesting and better readability would be using dictionaries. Each molecule would be a key in a dictionary {Molecule1 : Data, Molecule2 : Data} and the data would be dictionaries as well Atoms: Data, Bonds: Data. This requires quite a bit less nesting, but would require a little more work when I eventually want to display them.

Closing thoughts

Before asking this question, and indeed while writing it, I've been pretty set on dictionaries. It seems like a cleaner and easier method of storing this data, and won't require figuring out a really nasty recursive technique to walk through a molecule. However I wouldn't be surprised if I missed something which is why I'm asking here.

If this seems too implementation specific please let me know what you think I should change and I'll do my best to make it more general

  • 1
    Think about how you would represent Ferricyanide or Sulfur hexafloride. – user40980 Aug 18 '14 at 1:03
  • 2
    Disregarding the question. I must say that that is a very inefficient way to store chemical representations. I personally would recommend you store them similar to how SMILES does. Use a new object for each atom (ion, etc.) and have each object contain an array for its bonds where each item in the array is a reference to the object it bonds with (if that makes sense). Something like a linked list that works both ways and each "item" is either an array of separate lists or an atom (ion, etc.) – Cole Johnson Aug 18 '14 at 1:37
  • @ColeJohnson That is what I had in mind for the individual atoms, but I wanted something to represent the structure as well. The answers to this question (and that SMILES link) have started to persuade me that it would be entirely unnecessary – Dannnno Aug 18 '14 at 2:30
  • @MichaelT One of the many reasons I didn't like the matrix method – Dannnno Aug 18 '14 at 2:30
3

The "matrix method" has several problems:

  • Why only two dimensions? Very few molecules are flat, two dimensional structures (unless you work exclusively with graphene). Pretty much any program deals with chemical structures will want to work in 3D (excluding ones whose sole purpose is to draw structural formulas). By collapsing everything into 2D, you are losing valuable information.

  • What benefit would this method give? For one thing, finding adjacent atoms is rather complex: you have to look at adjacent neighbors in the matrix, of which you only get a maximum of 4 before you have to branch off using "Root". The artificial restriction of 4 neighbors in a matrix will make your code very ugly.

My only guess as to why you want to use this is because you want to eventually display the molecules. The display is the complicated part (flattening a 3D molecule into a 2D formula is a bit of an art even for a human), but you should seek to separate the display logic from the internal data that stores the molecular structure.


I'm not entirely clear what your "dictionary method" means: for example, what does Atoms denote? Are they just arbitrary labels, IDs, or the element's symbol? (Keep in mind that dictionary keys must be unique so element symbol won't work.)

You can take advantage of the fact that Python allows you to have cyclic references and construct graphs of the atoms. This allows you easily figure out the adjacent atoms.

Here's an example:

class Atom(object):
    def __init__(self, element, position):
        self.element = element
        self.position = position
        # each one for a different bond order, starting at 1
        self.bonds = [set(), set(), set(), set()]
    def bond_with(self, other, order=1):
        self.bonds[order].add(other)
        other.bonds[order].add(self)

# construct hydrogen fluoride
H = Atom("H", (0.0, 0.0, 0.0))
F = Atom("F", (0.0, 0.0, 1.0))
H.bond_with(F)
HF = set((H, F))
  • This answer. The best way to represent chemicals is as linked list of atoms or ions like how SMILES does. – Cole Johnson Aug 18 '14 at 1:39
  • The reason for only two dimensions was that it seemed harder to represent a structure in precise three dimensions - the Bond objects would contain information as to R/S of the bond and I'd have to write some functions that determined the approximate shapes of different structures (trigonal planar and all that). I apologize for being unclear with the dictionary method. I meant that each atom would have a unique id such that {id:atom_object}. All in all, this is along the lines of how I was considering making this work. Thanks! – Dannnno Aug 18 '14 at 2:33
2

Dictionaries are what Python's primary graph library NetworkX uses. Each graph structure is a dictionary of dictionaries. So it's nested in a sense but not by much. Molecules are very much representable by graphs, being connected groups of nodes.

If you want to skip all that and go straight to the modelling work, there's a full toolkit for that. The Molecular Modelling toolkit is a set of libraries for dealing with molecular modelling in Python. It is mostly aimed at Bioscience but I'm sure you could use it for other kinds of chemical work say Geochemistry or Materials research. It can even output to various visualizers.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.