How to get a bonded_structure from analysis.chemenv?

[fxcoudert]

When one wants to analyse a materials' structure in terms of a molecular graph / connectivity matrix / list of bonded atoms, one approach is to use the classes in analysis.local_env, like CrystallNN or IsayevNN. You can then get the bonding information through get_bonded_structure() as a nice StructureGraph object, which internally uses a networkx graph.

On the other hand, there is the approach of analysis.chemenv, which is a different technique but has a related aim. I am puzzled however about the objects that it returns. I think the closest analogue to get the structure graph is get_structure_connectivity(): that returns a StructureConnectivity, which internally uses a networkx graph. But it appears to be different from StructureGraph, although they must store pretty much the same information.

I have two questions:

1. How can I convert between StructureConnectivity and StructureGraph? or why can't I? am I missing something?
2. How am I suppose to manipulate/query a StructureConnectivity object? Ask it to return a list of bonds, for example? (Something like the get_connected_sites() function of StructureGraph.) I don't find a lot of doc or examples, and the only function I have found is a print_links() function to display stuff. Again, I am probably missing something, but I don't know what.

kind poke to @JaGeo because I think she wrote that code

[mkhorton]

Hi @fxcoudert, your confusion is warranted. Essentially, there was parallel development of quite similar functionality in pymatgen by two different groups.

They are not directly equivalent, and you can use or try either the chemenv module or the local_env module, and both have different algorithms for determining chemical environments.

The StructureGraph object is agnostic, it is simply a data structure for recording information about edges between sites in a Structure. It is convenient to use in part because it is also supported by Crystal Toolkit for visualization purposes, i.e. the visualisations on the Materials Project website are all showing StructureGraph objects. The get_bonded_structure method inside local_env is a convenience method to create one of these StructureGraph objects based off its coordination analysis, but you could equally create one using one of the chemenv algorithms too.

In an ideal world the API would be a little cleaner, but sometimes there is a historical path dependence, and overall I think it’s healthy to explore different approaches to solve a problem so I don’t think there’s much harm in it.

[mkhorton]

To add, both approaches also had their own publications, so I would refer to those if you want to understand their merits better.

To answer your question directly, I would iterate over the StructureConnectivity object to create the StructureGraph object manually. Indeed, this would be a good feature to add via PR, and maybe there is still an opportunity to clean up this API a bit to make it easier to understand.

[fxcoudert]

For reference, after some more reading, this is the code I currently have to convert from LightStructureEnvironments to StructureGraph:

def bonded_structure_from_lse(lse: LightStructureEnvironments) -> StructureGraph:
    graph = StructureGraph.from_empty_graph(lse.structure, name="bonds")

    for nset in lse.neighbors_sets:
        # Only works if there is a unique environment
        assert len(nset) == 1
        for neighbor in nset[0].neighb_indices_and_images:
            graph.add_edge(
                from_index=nset[0].isite,
                from_jimage=(0, 0, 0),
                to_index=neighbor["index"],
                to_jimage=neighbor["image_cell"],
                weight=1,
                edge_properties=None,
                warn_duplicates=False,
            )

    graph.set_node_attributes()
    return graph

It has restrictions that make it unsuitable for a PR at this stage: it requires a single well-defined environment to be present. I'm going to test it across the MP database.

The scientific context is: I am trying to compare the bond determination from ChemEnv to CrystalNN and others. I don't think it has been benchmarked before.

[JaGeo]

@fxcoudert
Yes, @mkhorton said: These developments happened in parallel (if i remember it correctly, chemenv started earlier). I thought to implement a converter at the time when I worked more on ChemEnv but never got to it. Most important for us were the coordination environments at the time. This is also what we are showing on the MP website.

Regarding your code: it should be possible to generalize it for all strategies (simple, multiweight), if you always use the first environment set.

I am tagging @kaueltzen from my group as she nowadays works mostly on ChemEnv and might want to comment.

[mkhorton]

Yes, I believe chemenv did start quite a bit earlier! I think the local_env code started in part as a benchmarking effort (which is why it contains implementations of a lot of already existing bonding/coordination schemes), so there was perhaps less overlap in the past, but then CrystalNN was implemented there too.