1 /* 2 Contributors: Nicola Zonta 3 Copyright Schrodinger, LLC. All rights reserved 4 */ 5 6 #ifndef COORDGEN_MACROCYCLE_BUILDER_H 7 #define COORDGEN_MACROCYCLE_BUILDER_H 8 9 #include "CoordgenConfig.hpp" 10 #include <cstdlib> 11 #include <iostream> 12 #include <vector> 13 14 class sketcherMinimizerAtom; 15 class sketcherMinimizerRing; 16 class sketcherMinimizerBond; 17 class sketcherMinimizerPointF; 18 19 struct doubleBondConstraint { 20 bool trans; 21 int previousAtom, atom1, atom2, followingAtom; 22 }; 23 24 struct ringConstraint { ringConstraintringConstraint25 ringConstraint(int a, sketcherMinimizerRing* r, bool fo) 26 { 27 ring = r; 28 atom = a; 29 forceOutside = fo; 30 } 31 bool forceOutside; 32 int atom; 33 sketcherMinimizerRing* ring; 34 }; 35 36 struct vertexCoords { 37 bool operator!=(const vertexCoords& rhs) const 38 { 39 if (x != rhs.x) { 40 return true; 41 } 42 if (y != rhs.y) { 43 return true; 44 } 45 if (z != rhs.z) { 46 return true; 47 } 48 return false; 49 } 50 bool operator<(const vertexCoords& rhs) const 51 { 52 if (x < rhs.x) { 53 return true; 54 } 55 if (y < rhs.y) { 56 return true; 57 } 58 if (z < rhs.z) { 59 return true; 60 } 61 return false; 62 } 63 friend const vertexCoords operator+(const vertexCoords& v1, 64 const vertexCoords& v2) 65 { 66 return vertexCoords(v1.x + v2.x, v1.y + v2.y, v1.z + v2.z); 67 } 68 friend const vertexCoords operator-(const vertexCoords& v1, 69 const vertexCoords& v2) 70 { 71 return vertexCoords(v1.x - v2.x, v1.y - v2.y, v1.z - v2.z); 72 } 73 74 bool operator==(const vertexCoords& rhs) const 75 { 76 return x == rhs.x && y == rhs.y && z == rhs.z; 77 } 78 vertexCoordsvertexCoords79 vertexCoords(int ix, int iy, int iz) 80 { 81 x = ix; 82 y = iy; 83 z = iz; 84 } 85 int x, y, z; 86 87 private: 88 friend std::ostream& operator<<(std::ostream& os, const vertexCoords& v); 89 }; 90 91 struct pathRestraints { 92 std::vector<int> heteroAtoms; 93 std::vector<std::pair<int, int>> substitutedAtoms; 94 }; 95 96 struct pathConstraints { 97 std::vector<doubleBondConstraint> doubleBonds; 98 std::vector<ringConstraint> ringConstraints; 99 std::vector<int> forceOutside; 100 }; 101 102 struct hexCoords { hexCoordshexCoords103 hexCoords(int ix, int iy) 104 { 105 x = ix; 106 y = iy; 107 } 108 bool operator==(const hexCoords& rhs) const 109 { 110 return x == rhs.x && y == rhs.y; 111 } 112 int x, y; distanceFromhexCoords113 int distanceFrom(const hexCoords& origin) const 114 { 115 int dx = abs(x - origin.x); 116 int dy = abs(y - origin.y); 117 int dz = abs((-x - y) - (-origin.x - origin.y)); 118 int max = (dx > dy) ? dx : dy; 119 return (dz > max) ? dz : max; 120 } rotate30DegreeshexCoords121 hexCoords rotate30Degrees() 122 { 123 int z = -x - y; 124 return {-z, -x}; 125 } toVertexCoordshexCoords126 vertexCoords toVertexCoords() const { return {x, y, -x - y}; } 127 128 private: 129 friend std::ostream& operator<<(std::ostream& os, const hexCoords& h); 130 }; 131 132 /* 133 unit hexagon used to build polyominoes for macrocycle shapes 134 */ 135 struct Hex { HexHex136 Hex(hexCoords coords) : m_coords(coords) {} setCoordsHex137 void setCoords(hexCoords coords) { m_coords = coords; } xHex138 int x() const { return m_coords.x; } yHex139 int y() const { return m_coords.y; } zHex140 int z() const { return -x() - y(); } coordsHex141 hexCoords coords() const { return m_coords; } 142 hexCoords m_coords; 143 std::vector<hexCoords> neighbors() const; 144 static std::vector<hexCoords> neighboringPositions(hexCoords h); 145 vertexCoords followingVertex(vertexCoords v) const; 146 }; 147 148 /* 149 hex polyomino (geometrical figure built on a hexagon lattice). All functions 150 assume that the polyomino has no holes 151 */ 152 class EXPORT_COORDGEN Polyomino 153 { 154 public: 155 Polyomino(); 156 Polyomino(const Polyomino& p); 157 ~Polyomino(); 158 Polyomino& operator=(const Polyomino& rhv); 159 160 /* 161 explore the topology of the polyominoes and returns true if they have the 162 same. Takes into account translations, rotations and mirroring 163 */ 164 bool isTheSameAs(Polyomino& p) const; 165 166 /* returns the number of hexagons in the polyomino */ 167 size_t size() const; 168 169 /* empties the polyomino */ 170 void clear(); 171 172 /* marks one hexagon to be a pentagon */ 173 void 174 markOneVertexAsPentagon(); // to get a path with an odd number of vertices 175 176 /* returns all hexagons that share the given vertex */ 177 std::vector<Hex*> vertexNeighbors(vertexCoords v) const; 178 179 /* 180 returns neighboring positions that are not shared with other hexagons 181 */ 182 std::vector<hexCoords> freeVertexNeighborPositions(vertexCoords v) const; 183 184 /* 185 return the contour path of the polyomino 186 */ 187 std::vector<vertexCoords> getPath() const; 188 189 /* return the hexagon at the given position */ 190 Hex* getHex(hexCoords coords) const; 191 192 /* find an outer vertex. Used to start path around polyomino. */ 193 vertexCoords findOuterVertex() const; 194 195 /* number of Hexs present at the given vertex */ 196 size_t hexagonsAtVertex(vertexCoords v) const; 197 198 /* build a round-ish polyomino with the given number of vertices */ 199 void buildWithVerticesN(int totVertices); 200 201 /* build a box-like polyomino of the given size */ 202 void buildSkewedBoxShape(int x, int y, 203 bool pentagon = false); // squared shape 204 205 /* 206 build a zig-zag box polyomino of the given size 207 */ 208 void buildRaggedBoxShape(int x, int y, bool pentagon = false); 209 210 /* 211 build a zig-zag box polyomino of the given size 212 */ 213 void buildRaggedSmallerBoxShape( 214 int x, int y, 215 bool pentagon = false); // box alternating rows of length x and x-1 216 217 /* 218 build a zig-zag box polyomino of the given size 219 */ 220 void buildRaggedBiggerBoxShape( 221 int x, int y, 222 bool pentagon = false); // box alternating rows of length x and x+1 223 224 /* 225 return vector of coordinates of unoccupied hexagons bordering the polyomino 226 */ 227 std::vector<hexCoords> allFreeNeighbors() const; 228 229 /* 230 return number of neighbors of hexagon at given coordinates 231 */ 232 int countNeighbors(hexCoords) const; 233 234 /* add an hexagon at given coordinates */ 235 void addHex(hexCoords coords); 236 237 /* remove the hexagon at given coordinates */ 238 void removeHex(hexCoords coords); 239 240 /* does removing this hexagon yield another polyomino with the same number 241 of vertices? true if the hexagon has 3 neighbors all next to each other 242 */ 243 bool isEquivalentWithout(hexCoords c) const; 244 245 /* give the coordinates of an hypotetical substituent bound to an atom in 246 * position pos */ 247 vertexCoords coordinatesOfSubstituent(vertexCoords pos) const; 248 249 // holds pointers to all hexagons 250 std::vector<Hex*> m_list; 251 252 // holds coordinates to vertices marked as pentagonal 253 std::vector<vertexCoords> pentagonVertices; 254 255 private: 256 /* mark vertex as pentagon (i.e. fuse it with neighbor vertex) */ 257 void setPentagon(vertexCoords c); 258 259 /* resize the grid to new size i */ 260 void resizeGrid(int i) const; 261 262 /* reassign hexagons to the grid */ 263 void reassignHexs() const; 264 265 /* get index of hexagon at the given coordinates */ 266 int getIndexInList(hexCoords coords) const; 267 268 /* hold pointers to all positions in the 269 grid, NULL pointers for empty positions */ 270 std::vector<Hex*> mutable m_grid; 271 272 /* size of the grid */ 273 mutable int m_gridSize; 274 }; 275 276 /* 277 class that builds coordinates for macrocycles 278 */ 279 class EXPORT_COORDGEN CoordgenMacrocycleBuilder 280 { 281 public: 282 CoordgenMacrocycleBuilder() = default; 283 ~CoordgenMacrocycleBuilder() = default; 284 285 /* assign coordinates to macrocycle */ 286 std::vector<sketcherMinimizerPointF> 287 newMacrocycle(sketcherMinimizerRing* r, 288 std::vector<sketcherMinimizerAtom*> atoms) const; 289 290 /* assign coordinates by removing one bond, generating the coordinates of 291 the resulting molecule and re-adding the missing bond */ 292 bool openCycleAndGenerateCoords(sketcherMinimizerRing* ring) const; 293 294 /* find most suitable bond to be broken */ 295 sketcherMinimizerBond* findBondToOpen(sketcherMinimizerRing* ring) const; 296 297 /* Add constraints to stereoactive double bonds to avoid inversions 298 public to be accessed by tests */ 299 std::vector<doubleBondConstraint> 300 getDoubleBondConstraints(std::vector<sketcherMinimizerAtom*>& atoms) const; 301 302 /* Skip the polyomino approach and fall back to opening the macrocycle when 303 * generating coordinates */ 304 bool m_forceOpenMacrocycles = false; 305 306 float getPrecision() const; 307 308 void setPrecision(float f); 309 310 private: 311 /* 312 precision of calculation. Higher values result in better results but longer 313 calculation times 314 */ 315 float m_precision; 316 317 /* build polyominos with the given number of vertices */ 318 std::vector<Polyomino> buildSquaredShapes(int totVertices) const; 319 320 /* remove duplicate polyominos from the list */ 321 std::vector<Polyomino> removeDuplicates(std::vector<Polyomino>& pols) const; 322 323 std::vector<ringConstraint> 324 getRingConstraints(std::vector<sketcherMinimizerAtom*>& atoms) const; 325 326 /* get total number of atoms bound to a that are not on parent size */ 327 int getNumberOfChildren(sketcherMinimizerAtom* a, 328 sketcherMinimizerAtom* parent) const; 329 330 /* get restraints to try to avoid things like having heteroatoms or 331 * substituents pointing inside of the macrocycle */ 332 pathRestraints 333 getPathRestraints(std::vector<sketcherMinimizerAtom*>& atoms) const; 334 335 /* a shape is not allowed to break a constraint (i.e. invert 336 stereochemistry on a double bond or put a fused ring inside 337 the macrocycle) */ 338 pathConstraints 339 getPathConstraints(std::vector<sketcherMinimizerAtom*>& atoms) const; 340 341 /* build a list of polyominoes with the same number of vertices by removing 342 * hexagons with 3 neighbors */ 343 std::vector<Polyomino> listOfEquivalent(const Polyomino& p) const; 344 345 std::vector<Polyomino> 346 listOfEquivalents(const std::vector<Polyomino>& l) const; 347 348 /* check that no ring constraints are violated */ 349 bool checkRingConstraints(std::vector<ringConstraint>& ringConstraints, 350 Polyomino& p, std::vector<vertexCoords>& path, 351 std::vector<int>& neighborNs, int& startI) const; 352 353 /* check that no double bond constraints are violated */ 354 bool checkDoubleBoundConstraints( 355 std::vector<doubleBondConstraint>& dbConstraints, 356 std::vector<vertexCoords>& vertices, int& startI) const; 357 358 /* score the path restraints */ 359 int scorePathRestraints(pathRestraints& pr, Polyomino& p, 360 std::vector<vertexCoords>& path, 361 std::vector<int>& neighborNs, int& startI) const; 362 363 /* check that no path constraints are violated */ 364 bool scorePathConstraints(pathConstraints& pc, Polyomino& p, 365 std::vector<vertexCoords>& path, 366 std::vector<int>& neighborNs, int& startI) const; 367 368 /* score restraints and constraints of the given path */ 369 int scorePath(Polyomino& p, std::vector<vertexCoords>& path, 370 std::vector<int>& neighborNs, int& startI, 371 pathConstraints& pc, pathRestraints& pr) const; 372 373 /* acceptable score to consider a shape appropriate */ 374 int acceptableShapeScore(int numberOfAtoms) const; 375 std::vector<int> getVertexNeighborNs(Polyomino& p, 376 std::vector<vertexCoords>& path) const; 377 378 /* return lowest period after which there is a rotation symmetry */ 379 int getLowestPeriod(std::vector<int>& neighbors) const; 380 381 /* explore the given polyominoes and score them */ 382 bool matchPolyominoes(std::vector<Polyomino>& pols, pathConstraints& pc, 383 pathRestraints& pr, int& bestP, int& bestScore, 384 int& bestStart, int& checkedMacrocycles) const; 385 386 /* explore the best starting point to place atoms around the given polyomino 387 */ 388 bool matchPolyomino(Polyomino& p, pathConstraints& pc, pathRestraints& pr, 389 int& bestStart, int& bestScore) const; 390 391 /* write the coordinates of the given path */ 392 void 393 writePolyominoCoordinates(std::vector<vertexCoords>& path, 394 const std::vector<sketcherMinimizerAtom*>& atoms, 395 int startI) const; 396 397 /* return the coordinates of the given vertex */ 398 sketcherMinimizerPointF coordsOfVertex(vertexCoords& v) const; 399 }; 400 401 #endif /* defined(COORDGEN_MACROCYCLE_BUILDER_H) */ 402