1 package Graph::CyclesDetection; 2 # 3 # $RCSfile: CyclesDetection.pm,v $ 4 # $Date: 2015/02/28 20:49:06 $ 5 # $Revision: 1.27 $ 6 # 7 # Author: Manish Sud <msud@san.rr.com> 8 # 9 # Copyright (C) 2015 Manish Sud. All rights reserved. 10 # 11 # This file is part of MayaChemTools. 12 # 13 # MayaChemTools is free software; you can redistribute it and/or modify it under 14 # the terms of the GNU Lesser General Public License as published by the Free 15 # Software Foundation; either version 3 of the License, or (at your option) any 16 # later version. 17 # 18 # MayaChemTools is distributed in the hope that it will be useful, but without 19 # any warranty; without even the implied warranty of merchantability of fitness 20 # for a particular purpose. See the GNU Lesser General Public License for more 21 # details. 22 # 23 # You should have received a copy of the GNU Lesser General Public License 24 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or 25 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330, 26 # Boston, MA, 02111-1307, USA. 27 # 28 29 use strict; 30 use Carp; 31 use Exporter; 32 use Graph; 33 use Graph::Path; 34 use Graph::PathGraph; 35 use BitVector; 36 37 use vars qw(@ISA @EXPORT @EXPORT_OK %EXPORT_TAGS); 38 39 @ISA = qw(Exporter); 40 @EXPORT = qw(); 41 @EXPORT_OK = qw(); 42 43 %EXPORT_TAGS = (all => [@EXPORT, @EXPORT_OK]); 44 45 # Setup class variables... 46 my($ClassName); 47 _InitializeClass(); 48 49 # Overload Perl functions... 50 use overload '""' => 'StringifyCyclesDetection'; 51 52 # Class constructor... 53 sub new { 54 my($Class, $Graph) = @_; 55 56 # Initialize object... 57 my $This = {}; 58 bless $This, ref($Class) || $Class; 59 $This->_InitializeCyclesDetection($Graph); 60 61 return $This; 62 } 63 64 # Initialize object data... 65 sub _InitializeCyclesDetection { 66 my($This, $Graph) = @_; 67 68 $This->{Graph} = $Graph; 69 70 # Cycles list... 71 @{$This->{AllCyclicPaths}} = (); 72 73 # Cyclic paths which are not part of any other cycle... 74 @{$This->{IndependentCyclicPaths}} = (); 75 76 return $This; 77 } 78 79 # Initialize class ... 80 sub _InitializeClass { 81 #Class name... 82 $ClassName = __PACKAGE__; 83 } 84 85 # Detect all cycles in graph... 86 # 87 sub DetectCycles { 88 my($This) = @_; 89 90 return $This->DetectCyclesUsingCollapsingPathGraphMethodology(); 91 } 92 93 # Detect all cycles in the graph using collapsing path graph [Ref 31] methodology... 94 # 95 # Note: 96 # . For topologically complex graphs containing large number of cycles, 97 # CollapseVertexAndCollectCyclicPathsDetectCycles method implemented in 98 # PathGraph can time out in which case no cycles are detected or 99 # assigned. 100 # 101 sub DetectCyclesUsingCollapsingPathGraphMethodology { 102 my($This) = @_; 103 my($PathGraph); 104 105 # Create a path graph... 106 $PathGraph = new Graph::PathGraph($This->{Graph}); 107 108 # For a vertex to be in a cycle, its degree must be >=2. So delete vertices recursively 109 # till all vertices with degree less than 2 have been deleted... 110 $PathGraph->DeleteVerticesWithDegreeLessThan(2); 111 112 # Setup a VertexID and EdgeID to index map usage during retrieval of independent cycles to 113 # avoid going over all vertices in all cylic paths later... 114 # 115 my($VertexIDsToIndicesRef, $LargestVertexIndex, $EdgeIDsToIndicesRef, $LargestEdgeIDIndex); 116 ($VertexIDsToIndicesRef, $LargestVertexIndex) = $This->_SetupVertexIDsToIndicesMap($PathGraph); 117 ($EdgeIDsToIndicesRef, $LargestEdgeIDIndex) = $This->_SetupEdgeIDsToIndicesMap($PathGraph); 118 119 # Recursively collapse vertices with lowest degree... 120 my($VertexID, $CycleVertexID); 121 while ($VertexID = $PathGraph->GetVertexWithSmallestDegree()) { 122 if (!$PathGraph->CollapseVertexAndCollectCyclicPaths($VertexID)) { 123 # Cycles detection didn't finish... 124 return undef; 125 } 126 } 127 128 # Get detected cycles and save 'em sorted by size... 129 push @{$This->{AllCyclicPaths}}, sort { $a->GetLength() <=> $b->GetLength() } $PathGraph->GetCyclicPaths(); 130 131 # Retrieve independent cyclic paths... 132 return $This->_RetrieveIndependentCycles($VertexIDsToIndicesRef, $LargestVertexIndex, $EdgeIDsToIndicesRef, $LargestEdgeIDIndex); 133 } 134 135 # Retrieve and save independent cyclic paths... 136 # 137 # Set of independent cycles identified using this method doesn't correspond to basis set of 138 # rings or smallest set of smallest rings (SSSR) [ Refs 29-30 ]; instead, set of cycles identified 139 # as independent cycles simply correspond to cycles which contain no other cycle as their 140 # subcycles and can't be described as linear combination of smaller cycles. And it also happen 141 # to contain all the rings in basis set of rings and SSSR. In other words, it's a superset of basis set 142 # of cycles and SSSR. For example, six four membered cycles are identified for cubane which is one 143 # more than the basis set of cycles. 144 # 145 sub _RetrieveIndependentCycles { 146 my($This, $VertexIDsToIndicesRef, $LargestVertexIndex, $EdgeIDsToIndicesRef, $LargestEdgeIDIndex) = @_; 147 148 # Only 1 or 0 cyclic paths... 149 if (@{$This->{AllCyclicPaths}} <= 1) { 150 push @{$This->{IndependentCyclicPaths}}, @{$This->{AllCyclicPaths}}; 151 return $This; 152 } 153 154 # Setup bit vectors for each cyclic path with size of each bit vector corresponding to 155 # maximum vertex index plus one... 156 my($CyclicPath, $BitVector, @BitNums, @CyclicPathBitVectors, @CyclicPathEdgeIDsBitVectors); 157 158 @CyclicPathBitVectors = (); @CyclicPathEdgeIDsBitVectors = (); 159 160 # Set bits corresponding to VertexIDs EdgeIDs using their indices... 161 for $CyclicPath (@{$This->{AllCyclicPaths}}) { 162 $BitVector = new BitVector($LargestVertexIndex); 163 @BitNums = map { $VertexIDsToIndicesRef->{$_} } $CyclicPath->GetVertices(); 164 $BitVector->SetBits(@BitNums); 165 push @CyclicPathBitVectors, $BitVector; 166 167 $BitVector = new BitVector($LargestEdgeIDIndex); 168 @BitNums = map { $EdgeIDsToIndicesRef->{$_} } $This->_GetPathEdgeIDs($CyclicPath); 169 $BitVector->SetBits(@BitNums); 170 push @CyclicPathEdgeIDsBitVectors, $BitVector; 171 } 172 173 # First smallest cyclic path always ends up as an independent cyclic path... 174 push @{$This->{IndependentCyclicPaths}}, $This->{AllCyclicPaths}[0]; 175 176 # Retrieve other independent cyclic paths... 177 my($CurrentIndex, $PreviousIndex, $CurrentCyclicPath, $PreviousCyclicPath, $CurrentPathLength, $PreviousPathLength, $CurrentBitVector, $PreviousBitVector, $CurrentAndPreviousBitVectpor, $AllPreviousSmallerPathsBitVector, $AllPreviousSmallerPathsEdgeIDsBitVector, $CurrentEdgeIDsBitVector, $AndBitVector, %SmallerPathAlreadyAdded, %SkipPath); 178 179 %SkipPath = (); 180 %SmallerPathAlreadyAdded = (); 181 $AllPreviousSmallerPathsBitVector = new BitVector($LargestVertexIndex); 182 $AllPreviousSmallerPathsEdgeIDsBitVector = new BitVector($LargestEdgeIDIndex); 183 184 CURRENT: for $CurrentIndex (1 .. $#{$This->{AllCyclicPaths}}) { 185 if (exists $SkipPath{$CurrentIndex}) { 186 next CURRENT; 187 } 188 $CurrentCyclicPath = $This->{AllCyclicPaths}[$CurrentIndex]; 189 $CurrentBitVector = $CyclicPathBitVectors[$CurrentIndex]; 190 $CurrentPathLength = $CurrentCyclicPath->GetLength(); 191 192 PREVIOUS: for $PreviousIndex (0 .. ($CurrentIndex - 1)) { 193 if (exists $SkipPath{$PreviousIndex}) { 194 next PREVIOUS; 195 } 196 $PreviousCyclicPath = $This->{AllCyclicPaths}[$PreviousIndex]; 197 $PreviousBitVector = $CyclicPathBitVectors[$PreviousIndex]; 198 199 # Is previous path a subset of current path? 200 $CurrentAndPreviousBitVectpor = $PreviousBitVector & $CurrentBitVector; 201 if ($PreviousBitVector->GetNumOfSetBits() == $CurrentAndPreviousBitVectpor->GetNumOfSetBits()) { 202 # Current path doesn't qualify an independent path... 203 $SkipPath{$CurrentIndex} = 1; 204 next CURRENT; 205 } 206 207 $PreviousPathLength = $PreviousCyclicPath->GetLength(); 208 if ($PreviousPathLength < $CurrentPathLength) { 209 # Mark cycle vertices seen in cyclic paths with length smaller than current path... 210 if (! exists $SmallerPathAlreadyAdded{$PreviousIndex}) { 211 $SmallerPathAlreadyAdded{$PreviousIndex} = 1; 212 $AllPreviousSmallerPathsBitVector = $AllPreviousSmallerPathsBitVector | $PreviousBitVector; 213 $AllPreviousSmallerPathsEdgeIDsBitVector = $AllPreviousSmallerPathsEdgeIDsBitVector | $CyclicPathEdgeIDsBitVectors[$PreviousIndex]; 214 } 215 } 216 } 217 if ($AllPreviousSmallerPathsBitVector->GetNumOfSetBits()) { 218 # Is current path a linear combination of smaller paths? 219 $AndBitVector = $AllPreviousSmallerPathsBitVector & $CurrentBitVector; 220 if ($CurrentBitVector->GetNumOfSetBits() == $AndBitVector->GetNumOfSetBits()) { 221 # Are all edges in the current path already present in smaller paths? 222 $CurrentEdgeIDsBitVector = $CyclicPathEdgeIDsBitVectors[$CurrentIndex]; 223 $AndBitVector = $AllPreviousSmallerPathsEdgeIDsBitVector & $CurrentEdgeIDsBitVector; 224 225 if ($CurrentEdgeIDsBitVector->GetNumOfSetBits() == $AndBitVector->GetNumOfSetBits()) { 226 # Current path doesn't qualify an independent path... 227 $SkipPath{$CurrentIndex} = 1; 228 next CURRENT; 229 } 230 } 231 } 232 # It's an independent cyclic path... 233 push @{$This->{IndependentCyclicPaths}}, $CurrentCyclicPath; 234 } 235 return $This; 236 } 237 238 # Setup a VertexID to index map... 239 # 240 sub _SetupVertexIDsToIndicesMap { 241 my($This, $PathGraph) = @_; 242 my($LargestVertexIndex, @VertexIDs, %VertexIDsMap); 243 244 %VertexIDsMap = (); @VertexIDs = (); $LargestVertexIndex = 0; 245 246 @VertexIDs = $PathGraph->GetVertices(); 247 if (! @VertexIDs) { 248 return (\%VertexIDsMap, $LargestVertexIndex); 249 } 250 @VertexIDsMap{ @VertexIDs } = (0 .. $#VertexIDs); 251 $LargestVertexIndex = scalar @VertexIDs; 252 253 return (\%VertexIDsMap, $LargestVertexIndex); 254 } 255 256 # Setup a Edge to index map using paths associated to egdes in an intial 257 # path graph... 258 # 259 sub _SetupEdgeIDsToIndicesMap { 260 my($This, $PathGraph) = @_; 261 my($Path, $LargestEdgeIndex, @EdgeIDs, %EdgeIDsMap); 262 263 %EdgeIDsMap = (); @EdgeIDs = (); $LargestEdgeIndex = 0; 264 265 @EdgeIDs = (); 266 for $Path ($PathGraph->GetPaths()) { 267 push @EdgeIDs, $This->_GetPathEdgeIDs($Path); 268 } 269 270 if (! @EdgeIDs) { 271 return (\%EdgeIDsMap, $LargestEdgeIndex); 272 } 273 274 @EdgeIDsMap{ @EdgeIDs } = (0 .. $#EdgeIDs); 275 $LargestEdgeIndex = scalar @EdgeIDs; 276 277 return (\%EdgeIDsMap, $LargestEdgeIndex); 278 } 279 280 # Get path edge IDs or number of edges. The edge IDs are generated from 281 # edge vertices and correpond to VertexID1-VertexID2 where VertexID1 <= 282 # VertexID2. 283 # 284 sub _GetPathEdgeIDs { 285 my($This, $Path) = @_; 286 my(@EdgesVertexIDs, @EdgeIDs); 287 288 @EdgesVertexIDs = (); @EdgeIDs = (); 289 @EdgesVertexIDs = $Path->GetEdges(); 290 if (!@EdgesVertexIDs) { 291 return wantarray ? @EdgeIDs : (scalar @EdgeIDs); 292 } 293 294 # Set up edge IDs... 295 my($Index, $VertexID1, $VertexID2, $EdgeID); 296 297 for ($Index = 0; $Index < $#EdgesVertexIDs; $Index += 2) { 298 $VertexID1 = $EdgesVertexIDs[$Index]; $VertexID2 = $EdgesVertexIDs[$Index + 1]; 299 $EdgeID = ($VertexID1 <= $VertexID2) ? "$VertexID1-$VertexID2" : "$VertexID2-$VertexID1"; 300 push @EdgeIDs, $EdgeID; 301 } 302 303 return wantarray ? @EdgeIDs : (scalar @EdgeIDs); 304 } 305 306 # Return an array containing references to cyclic paths or number of cylic paths... 307 # 308 sub GetAllCyclicPaths { 309 my($This) = @_; 310 311 return wantarray ? @{$This->{AllCyclicPaths}} : scalar @{$This->{AllCyclicPaths}}; 312 } 313 314 # Get cyclic paths which are independent. In otherwords, cycles which don't contain any other 315 # cycle as their subset... 316 # 317 sub GetIndependentCyclicPaths { 318 my($This) = @_; 319 320 return wantarray ? @{$This->{IndependentCyclicPaths}} : scalar @{$This->{IndependentCyclicPaths}}; 321 } 322 323 # Return a string containg data for CyclesDetection object... 324 sub StringifyCyclesDetection { 325 my($This) = @_; 326 my($CyclesDetectionString, $CyclesCount, $CyclicPath); 327 328 $CyclesCount = @{$This->{AllCyclicPaths}}; 329 $CyclesDetectionString = "AllCycles: Count - $CyclesCount"; 330 for $CyclicPath (@{$This->{AllCyclicPaths}}) { 331 $CyclesDetectionString .= "; Cycle: " . join('-', $CyclicPath->GetVertices()); 332 } 333 334 $CyclesCount = @{$This->{IndependentCyclicPaths}}; 335 $CyclesDetectionString .= "\nIndependentCycles: Count - $CyclesCount"; 336 for $CyclicPath (@{$This->{IndependentCyclicPaths}}) { 337 $CyclesDetectionString .= "; Cycle: " . join('-', $CyclicPath->GetVertices()); 338 } 339 340 return $CyclesDetectionString; 341 } 342 343 # Return a reference to new cycle detection object... 344 sub Copy { 345 my($This) = @_; 346 my($NewCyclesDetection); 347 348 $NewCyclesDetection = Storable::dclone($This); 349 350 return $NewCyclesDetection; 351 } 352