/* * Project: MoleCuilder * Description: creates and alters molecular systems * Copyright (C) 2021 Frederik Heber. All rights reserved. * * * This file is part of MoleCuilder. * * MoleCuilder is free software: you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free Software Foundation, either version 2 of the License, or * (at your option) any later version. * * MoleCuilder is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with MoleCuilder. If not, see . */ /* * Graph6Writer.cpp * * Created on: Apr 2, 2021 * Author: heber */ // include config.h #ifdef HAVE_CONFIG_H #include #endif #include "Graph6Writer.hpp" #include "CodePatterns/Assert.hpp" #include "CodePatterns/Log.hpp" #include #include #include #include "Atom/atom.hpp" #include "Descriptors/AtomIdDescriptor.hpp" #include "Element/element.hpp" #include "Graph/BoostGraphCreator.hpp" #include "Graph/BreadthFirstSearchGatherer.hpp" #include "World.hpp" //#include "CodePatterns/MemDebug.hpp" Graph6Writer::Graph6Writer(const std::vector atoms): _atoms(atoms) {} void Graph6Writer::write_n(std::ostream& out) { const unsigned long n = _atoms.size(); if (n<62) { out << ((unsigned char)(n+63)); return; } out << ((unsigned char)126); int num_bytes = 2; if (n> 258047) { out << ((unsigned char)126); num_bytes = 3; } for(int value=num_bytes; value>=0; value--) { unsigned char c = 0; int n_pos = 6*(value+1)-1; for(int c_pos=5; c_pos>=0; n_pos--, c_pos--) { c += (n & (1<>((int)n_pos/6); } out << (c+63); } } /* Given an iterator over the adjacency matrix in the order (0,1),(0,2),(1,2),(0,3),(1,3),(2,3),...,(n-1,n) this writes a graph6 representation to out. */ void Graph6Writer::write_graph6(std::ostream& out) { write_n(out); const unsigned long n = _atoms.size(); unsigned char value = 0; int byte_pos = 5; unsigned int bytes_written = 0; for (size_t j=0; jIsBondedTo(_atoms[j]); LOG(2, "DEBUG: (" << i << "," << j << ") = " << bit << "," << value << " | " << bit << " << " << byte_pos << " = " << (unsigned int)value << " | " << (bit << byte_pos)); value = value | (bit << byte_pos--); if (byte_pos < 0) { LOG(2, "DEBUG: Writing byte " << value << " into range [" << (unsigned char)63 << "," << (unsigned char)126 << "]"); ASSERT( (value+63) <= 126, "Graph6Writer::write_graph6() - char to write is outside "+toString((unsigned char)63) +" and "+toString((unsigned char)126)); out << (unsigned char)(value+63); bytes_written++; value = 0; byte_pos = 5; } } if (byte_pos!=5) { ASSERT( (value+63) <= 126, "Graph6Writer::write_graph6() - char to write is outside "+toString((unsigned char)63) +" and "+toString((unsigned char)126)); LOG(2, "DEBUG: Writing byte " << value << " into range [" << (unsigned char)63 << "," << (unsigned char)126 << "]"); out << (unsigned char)(value+63); bytes_written++; value=0; } ASSERT( value==0, "Graph6Writer::write_graph6() - byte is not null, i.e. chars left to write?"); ASSERT( bytes_written == (unsigned int)ceil(n*(n-1)/12.0f), "Graph6Writer::write_graph6() - unexpected number of bytes written"); } /** * Picks a non-hydrogen from all hydrogen atoms in the current set * of atoms. * * Returns -1 if none could be found. */ atomId_t Graph6Writer::getBoundaryNonHydrogen() const { atomId_t start_atom_id = -1; int bond_degree = 16; for(std::vector::const_iterator iter = _atoms.begin(); iter != _atoms.end(); ++iter) { const atom *walker = *iter; 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; } } } return start_atom_id; } bool OnlyNonHydrogens(const bond &_bond) { return _bond.HydrogenBond == 0; } void Graph6Writer::write_elementlist(std::ostream& out) { /** Execute a Breadth-First Search discovery from one terminal atom (e.g., * pick random hydrogen and then it's bond-neighbor if it is non-hydrogen). * Then return the element list in that ordering. * * The graph6 string does not account for the inherent graph symmetries * (e.g., BW having 123<->321 but not 123<->132 symmetry). */ // pick bond neighbor of a hydrogen atom atomId_t start_atom_id = getBoundaryNonHydrogen(); if (start_atom_id == (unsigned int)-1) { // fall back to first atom in list start_atom_id = _atoms.front()->getId(); } // do an unlimited BFS and get set of nodes, ordered by discovery level BoostGraphCreator graphCreator; graphCreator.createFromAtoms(_atoms, OnlyNonHydrogens); BreadthFirstSearchGatherer gatherer(graphCreator); const std::vector gathered_atoms = gatherer(start_atom_id); // print all nodes const World& world = World::getConstInstance(); for (std::vector::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(); } }