Changeset 78e5cf for src/Graph/Graph6Writer.cpp

Timestamp:

Oct 28, 2022, 7:23:01 AM (3 years ago)

Author:

Frederik Heber <frederik.heber@…>

Branches:

Candidate_v1.7.0, stable

Children:

80bf6e

Parents:

1d9586

git-author:

Frederik Heber <frederik.heber@…> (10/25/22 16:00:09)

git-committer:

Frederik Heber <frederik.heber@…> (10/28/22 07:23:01)

Message:

Graph6Writer::write_elementlist uses least connected non-hydrogen.

• this should enforce outermost atom in linear chain.

File:

• src/Graph/Graph6Writer.cpp (modified) (5 diffs)

src/Graph/Graph6Writer.cpp

@@ -125 +125 @@
 
 /**
- * Picks a non-hydrogen from all hydrogen atoms in the current set
- * of atoms.
+ * Picks a non-hydrogen from all atoms in the current set of atoms
+ * with lowest non-hydrogen bonds.
  *
  * Returns -1 if none could be found.
@@ -132 +132 @@
 atomId_t Graph6Writer::getBoundaryNonHydrogen() const {
   atomId_t start_atom_id = -1;
-  int bond_degree = 16;
+  int lowest_non_hydrogen_count = 16;
   for(std::vector<const atom *>::const_iterator iter = _atoms.begin();
       iter != _atoms.end(); ++iter) {
     const atom *walker = *iter;
-    if (walker->getElement().getSymbol() == "H") {
+    if (walker->getElement().getSymbol() != "H") {
       const BondList& bond_list = walker->getListOfBonds();
-      assert(bond_list.size() == 1);
-      const atom *other = bond_list.front().get()->GetOtherAtom(walker);
-      const int number_of_bonds = other->getListOfBonds().size();
-      if ((other->getElement().getSymbol() != "H") && (bond_degree > number_of_bonds)) {
-        start_atom_id = other->getId();
-        bond_degree = number_of_bonds;
+      int number_of_non_hydrogen_bonds = 0;
+      for (BondList::const_iterator iter = bond_list.begin();
+        iter != bond_list.end(); ++iter) {
+        number_of_non_hydrogen_bonds += (*iter)->GetOtherAtom(walker)->getElement().getSymbol() != "H";
       }
-    }
+      if (lowest_non_hydrogen_count > number_of_non_hydrogen_bonds) {
+        start_atom_id = walker->getId();
+        lowest_non_hydrogen_count = number_of_non_hydrogen_bonds;
+      }
+    }
+  }
+  if ((start_atom_id == -1) && (!_atoms.empty())) {
+    // we only have hydrogens, just pick the first
+    start_atom_id = (*_atoms.begin())->getId();
   }
   return start_atom_id;
@@ -162 +168 @@
    * (e.g., BW having 123<->321 but not 123<->132 symmetry).
    */
+  const World& world = World::getConstInstance();
   // pick bond neighbor of a hydrogen atom
   atomId_t start_atom_id = getBoundaryNonHydrogen();
@@ -168 +175 @@
     start_atom_id = _atoms.front()->getId();
   }
+  const atom* start_atom = world.getAtom(AtomById(start_atom_id));
+  LOG(1, "INFO: Start atom is " << *start_atom << ".");
 
   // do an unlimited BFS and get set of nodes, ordered by discovery level
@@ -173 +182 @@
   graphCreator.createFromAtoms(_atoms, OnlyNonHydrogens);
   BreadthFirstSearchGatherer gatherer(graphCreator);
-  const std::vector<atomId_t> gathered_atoms = gatherer(start_atom_id);
-
-  // print all nodes
-  const World& world = World::getConstInstance();
-  for (std::vector<atomId_t>::const_iterator iter = gathered_atoms.begin();
-      iter != gathered_atoms.end(); ++iter) {
-    if (iter != gathered_atoms.begin())
-      out << ' ';
-    const atom* walker = world.getAtom(AtomById(*iter));
-    assert(walker != NULL);
-    out << walker->getElement().getSymbol();
-  }
-}
-
+  gatherer(start_atom_id);
+
+  // go through distance map and print sorted by discovery level
+  const BreadthFirstSearchGatherer::distance_map_t &distances = gatherer.getDistances();
+  using pairtype = std::pair<atomId_t, size_t>;
+  const size_t max_distance = std::max_element(distances.begin(), distances.end(), [] (const pairtype & p1, const pairtype & p2) {
+    return p1.second < p2.second;
+  })->second;
+  bool isFirst = true;
+  /**
+   * This is O(N^2) and a stupid implementation. However, we only intend to
+   * use this for small molecules, so I don't care at the moment. The better
+   * approach is to revert the map into a multimap and then traverse that.
+   */
+  for (size_t i=0; i<= max_distance; ++i) {
+    for (BreadthFirstSearchGatherer::distance_map_t::const_iterator iter = distances.begin();
+        iter != distances.end(); ++iter) {
+      if (iter->second != i)
+        continue;
+      const atom* walker = world.getAtom(AtomById(iter->first));
+      assert(walker != NULL);
+      LOG(1, "INFO: Gathered atom " << *walker);
+      if (!isFirst)
+        out << ' ';
+      isFirst = false;
+      out << walker->getElement().getSymbol();
+    }
+  }
+}
+