The CyMeAlg Software (Release 0.01)
various.h
Go to the documentation of this file.
1 /// @addtogroup cymealg
2 /// @{
3 
4 /////////////////////////////////////////////////////////////////////////////
5 ///
6 /// @file
7 ///
8 /// This header file contains implementation of a selection
9 /// of various graph algorithms.
10 ///
11 /// @author Pawel Pilarczyk
12 ///
13 /////////////////////////////////////////////////////////////////////////////
14 
15 // Copyright (C) 1997-2020 by Pawel Pilarczyk.
16 //
17 // This file is part of my research software. This is free software:
18 // you can redistribute it and/or modify it under the terms of the GNU
19 // General Public License as published by the Free Software Foundation,
20 // either version 3 of the License, or (at your option) any later version.
21 //
22 // This software is distributed in the hope that it will be useful,
23 // but WITHOUT ANY WARRANTY; without even the implied warranty of
24 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 // GNU General Public License for more details.
26 //
27 // You should have received a copy of the GNU General Public License
28 // along with this software; see the file "license.txt". If not,
29 // please, see <http://www.gnu.org/licenses/>.
30 
31 // Started in January 2006. Last revision: August 23, 2019.
32 
33 
34 #ifndef _CYMEALG_VARIOUS_H_
35 #define _CYMEALG_VARIOUS_H_
36 
37 // include selected header files from the CHomP library
38 #include "chomp/system/config.h"
39 
40 namespace cymealg {
41 
42 /// A helper function for "collapseVertices".
43 template <class wType>
44 inline int_t addRenumEdges (const diGraph<wType> &g, int_t vertex,
45  const int_t *newNum, int_t cur, int_t *srcVert,
46  diGraph<wType> &result)
47 {
48  int_t nEdges = g. countEdges (vertex);
49  // add all the edges that start at the component
50  for (int_t edge = 0; edge < nEdges; ++ edge)
51  {
52  // determine the dest. vertex of the edge
53  int_t dest = newNum [g. getEdge (vertex, edge)];
54  // if this is an edge to self, then ignore it
55  if (dest == cur)
56  continue;
57  // if the edge has already been added,
58  // then skip it
59  if (srcVert [dest] == cur)
60  continue;
61  // remember that the dest. vertex has an edge
62  // pointing out from the current vertex
63  srcVert [dest] = cur;
64  // add the edge to the result graph
65  result. addEdge (dest);
66  }
67  return 0;
68 } /* addRenumEdges */
69 
70 /// Collapses disjoint subsets of vertices to single vertices and creates
71 /// a corresponding graph in which the first vertices come from
72 /// the collapsed ones. The result graph must be initially empty.
73 /// Saves translation table from old to new vertex numbers.
74 /// The table must have sufficient length.
75 template <class wType, class Table1, class Table2>
76 inline int_t collapseVertices (const diGraph<wType> &g,
77  int_t compNum, Table1 &compVertices, Table2 &compEnds,
78  diGraph<wType> &result, int_t *newNumbers = 0)
79 {
80  // do nothing if the input graph is empty
81  int_t nVert = g. countVertices ();
82  if (!nVert)
83  return 0;
84 
85  // compute the new numbers of vertices: newNum [i] is the
86  // number of vertex 'i' from 'g' in the resulting graph
87  int_t *newNum = newNumbers ? newNumbers : new int_t [nVert];
88  for (int_t i = 0; i < nVert; ++ i)
89  newNum [i] = -1;
90  int_t cur = 0; // current vertex number in the new graph
91  int_t pos = 0; // position in compVertices
92  while (cur < compNum)
93  {
94  int_t posEnd = compEnds [cur];
95  while (pos < posEnd)
96  newNum [compVertices [pos ++]] = cur;
97  ++ cur;
98  }
99  for (int_t i = 0; i < nVert; ++ i)
100  {
101  if (newNum [i] < 0)
102  newNum [i] = cur ++;
103  }
104 
105  // prepare a table to mark the most recent source vertex for edges:
106  // srcVert [j] contains i if the edge i -> j has just been added
107  int_t newVert = nVert - (compNum ? compEnds [compNum - 1] :
108  static_cast<int_t> (0)) + compNum;
109  // debug:
110  if (cur != newVert)
111  throw "DEBUG: Wrong numbers.";
112  int_t *srcVert = new int_t [newVert];
113  for (int_t i = 0; i < newVert; ++ i)
114  srcVert [i] = -1;
115 
116  // scan the input graph and create the resulting graph: begin with
117  // the vertices in the collapsed groups
118  cur = 0, pos = 0;
119  while (cur < compNum)
120  {
121  // add a new vertex to the result graph
122  result. addVertex ();
123  // for all the vertices in this component...
124  int_t posEnd = compEnds [cur];
125  while (pos < posEnd)
126  {
127  // take the next vertex from the component
128  int_t vertex = compVertices [pos ++];
129  // add the appropriate edges to the result graph
130  addRenumEdges (g, vertex, newNum, cur,
131  srcVert, result);
132  }
133  ++ cur;
134  }
135  // process the remaining vertices of the graph
136  for (int_t vertex = 0; vertex < nVert; ++ vertex)
137  {
138  // debug
139  if (newNum [vertex] > cur)
140  throw "DEBUG: Wrong order.";
141  // if the vertex has already been processed, skip it
142  if (newNum [vertex] != cur)
143  continue;
144  // add the appropriate edges to the result graph
145  result. addVertex ();
146  addRenumEdges (g, vertex, newNum, cur, srcVert, result);
147  ++ cur;
148  }
149 
150  // free memory and exit
151  delete [] srcVert;
152  if (!newNumbers)
153  delete [] newNum;
154  return 0;
155 } /* collapseVertices */
156 
157 /// A helper function for "collapseVertices2".
158 template <class wType, class TabSets>
159 inline int_t addRenumEdges2 (const diGraph<wType> &g, int_t vertex,
160  const int_t *newNum, TabSets &numSets,
161  int_t cur, int_t *srcVert, diGraph<wType> &result)
162 {
163  int_t nEdges = g. countEdges (vertex);
164 
165  // add all the edges that start at this vertex
166  for (int_t edge = 0; edge < nEdges; ++ edge)
167  {
168  // determine the dest. vertex of the edge
169  int_t destNumber = newNum [g. getEdge (vertex, edge)];
170 
171  // consider all the destination vertices
172  int_t destSize = (destNumber < 0) ?
173  numSets [-destNumber]. size () :
174  static_cast<int_t> (1);
175  for (int_t i = 0; i < destSize; ++ i)
176  {
177  // determine the consecutive destination vertex
178  int_t dest = (destNumber < 0) ?
179  numSets [-destNumber] [i] : destNumber;
180 
181  // if this is an edge to self, then ignore it
182  if (dest == cur)
183  continue;
184 
185  // if the edge has already been added,
186  // then skip it
187  if (srcVert [dest] == cur)
188  continue;
189 
190  // add the edge to the result graph
191  result. addEdge (dest);
192  }
193  }
194  return 0;
195 } /* addRenumEdges2 */
196 
197 /// Collapses (possibly non-disjoint) subsets of vertices
198 /// to single vertices and creates a corresponding graph
199 /// in which the first vertices come from the collapsed ones.
200 /// The result graph must be initially empty.
201 template <class wType, class Table1, class Table2>
202 inline int_t collapseVertices2 (const diGraph<wType> &g,
203  int_t compNum, Table1 &compVertices, Table2 &compEnds,
204  diGraph<wType> &result)
205 {
206  // do nothing if the input graph is empty
207  int_t nVert = g. countVertices ();
208  if (!nVert)
209  return 0;
210 
211  // compute the new numbers of vertices: newNum [i] is the
212  // number of vertex 'i' from 'g' in the resulting graph,
213  // unless negative; then it points to a set of numbers,
214  // with -1 meaning "not defined, yet"
215  int_t *newNum = new int_t [nVert];
216  for (int_t i = 0; i < nVert; ++ i)
217  newNum [i] = -1;
218 
219  // sets of numbers of vertices; the numbers refering to the sets
220  // begin with 2, thus the first two sets are skipped
221  chomp::homology::multitable<chomp::homology::hashedset<int_t> >
222  numSets;
223  // the number of the set waiting to be used next
224  int_t numSetCur = 2;
225 
226  int_t cur = 0; // current vertex number in the new graph
227  int_t pos = 0; // position in compVertices
228  while (cur < compNum)
229  {
230  int_t posEnd = compEnds [cur];
231  while (pos < posEnd)
232  {
233  int_t number = compVertices [pos ++];
234  if (newNum [number] == -1)
235  newNum [number] = cur;
236  else if (newNum [number] < 0)
237  numSets [-newNum [number]]. add (cur);
238  else
239  {
240  numSets [numSetCur]. add (newNum [number]);
241  numSets [numSetCur]. add (cur);
242  newNum [number] = -(numSetCur ++);
243  }
244  }
245  ++ cur;
246  }
247  for (int_t i = 0; i < nVert; ++ i)
248  {
249  if (newNum [i] == -1)
250  newNum [i] = cur ++;
251  }
252 
253  // prepare a table to mark the most recent source vertex for edges:
254  // srcVert [j] contains i if the edge i -> j has just been added
255  int_t newVert = cur;
256  int_t *srcVert = new int_t [newVert];
257  for (int_t i = 0; i < newVert; ++ i)
258  srcVert [i] = -1;
259 
260  // scan the input graph and create the resulting graph: begin with
261  // the vertices in the collapsed groups
262  cur = 0;
263  pos = 0;
264  while (cur < compNum)
265  {
266  // add a new vertex to the result graph
267  result. addVertex ();
268 
269  // for all the vertices in this component...
270  int_t posEnd = compEnds [cur];
271  while (pos < posEnd)
272  {
273  // take the next vertex from the component
274  int_t vertex = compVertices [pos ++];
275 
276  // add the appropriate edges to the result graph
277  addRenumEdges2 (g, vertex, newNum, numSets, cur,
278  srcVert, result);
279  }
280  ++ cur;
281  }
282 
283  // process the remaining vertices of the graph
284  for (int_t vertex = 0; vertex < nVert; ++ vertex)
285  {
286  // debug
287  if (newNum [vertex] > cur)
288  throw "DEBUG: Wrong order 2.";
289 
290  // if the vertex has already been processed, skip it
291  if (newNum [vertex] != cur)
292  continue;
293 
294  // add the appropriate edges to the result graph
295  result. addVertex ();
296  addRenumEdges2 (g, vertex, newNum, numSets, cur,
297  srcVert, result);
298  ++ cur;
299  }
300 
301  // free memory and exit
302  delete [] srcVert;
303  delete [] newNum;
304  return 0;
305 } /* collapseVertices2 */
306 
307 /// Computes the graph that represents connections between a number
308 /// of the first vertices in the given graph.
309 /// The vertices of the result graph are the first "nVert" vertices
310 /// from the source graph. An edge is added to the new graph
311 /// iff there exists a path from one vertex to another,
312 /// without going through any other vertex in that set.
313 /// Runs DFS starting from each of the first "nVert" vertices,
314 /// and thus may be a little inefficient in some cases.
315 template <class wType>
316 inline int_t connectionGraph (const diGraph<wType> &g, int_t nVert,
317  diGraph<wType> &result)
318 {
319  // remember the number of vertices in the input graph
320  int_t nVertG = g. countVertices ();
321  if (!nVertG)
322  return 0;
323 
324  // prepare a bitfield in which visited vertices will be marked
325  chomp::homology::BitField visited;
326  visited. allocate (nVertG);
327  visited. clearall (nVertG);
328 
329  // run DFS for each of the starting vertices
330  for (int_t startVertex = 0; startVertex < nVert; ++ startVertex)
331  {
332  // add this vertex to the resulting graph
333  result. addVertex ();
334 
335  // prepare a counter and a stack of visited vertices
336  int_t nVisited = 0;
337  std::stack<int_t> s_visited;
338 
339  // prepare stacks for the DFS algorithm
340  std::stack<int_t> s_vertex;
341  std::stack<int_t> s_edge;
342 
343  // mark the starting vertex as visited
344  visited. set (startVertex);
345  s_visited. push (startVertex);
346  ++ nVisited;
347 
348  // use DFS to visit vertices reachable from that vertex
349  int_t vertex = startVertex;
350  int_t edge = 0;
351  while (1)
352  {
353  // go back with the recursion
354  // if all the edges have been processed
355  if (edge >= g. countEdges (vertex))
356  {
357  // if this was the top recursion level
358  // then quit
359  if (s_vertex. empty ())
360  break;
361 
362  // return from the recursion
363  vertex = s_vertex. top ();
364  s_vertex. pop ();
365  edge = s_edge. top ();
366  s_edge. pop ();
367  continue;
368  }
369 
370  // take the next vertex
371  int_t next = g. getEdge (vertex, edge ++);
372 
373  // go to the deeper recursion level if not visited
374  if (!visited. test (next))
375  {
376  // add an edge to the result if necessary
377  if (next < nVert)
378  result. addEdge (next);
379 
380  // store the previous variables at the stacks
381  s_vertex. push (vertex);
382  s_edge. push (edge);
383 
384  // take the new vertex
385  vertex = next;
386  edge = 0;
387 
388  // mark the new vertex as visited
389  visited. set (vertex);
390  s_visited. push (vertex);
391  ++ nVisited;
392  }
393  }
394 
395  // if this was the last vertex then skip the rest
396  if (startVertex == nVert - 1)
397  break;
398 
399  // mark each vertex as non-visited
400  if (nVisited > (nVertG >> 6))
401  {
402  visited. clearall (nVertG);
403  }
404  else
405  {
406  while (!s_visited. empty ())
407  {
408  int_t vertex = s_visited. top ();
409  s_visited. pop ();
410  visited. clear (vertex);
411  }
412  }
413  }
414 
415  // free memory and exit
416  visited. free ();
417  return 0;
418 } /* connectionGraph */
419 
420 /// Computes the period of a strongly connected graph.
421 /// The period of a graph is defined as the greatest common divisor
422 /// of the lengths of all the cycles in the graph.
423 /// Period 1 means that the graph is aperiodic.
424 /// The complexity of this operation is linear (one DFS).
425 template <class GraphType>
426 inline int_t computePeriod (const GraphType &g)
427 {
428  // remember the number of vertices in the input graph
429  int_t nVertG = g. countVertices ();
430  if (!nVertG)
431  return 0;
432 
433  // prepare an array to record the depth of each visited vertex
434  std::vector<int_t> visited (nVertG, 0);
435 
436  // run DFS starting at the first vertex
437  // prepare stacks for the DFS algorithm
438  std::stack<int_t> s_vertex;
439  std::stack<int_t> s_edge;
440 
441  // mark the starting vertex as visited
442  visited [0] = 1;
443 
444  // use DFS to visit vertices reachable from that vertex
445  int_t vertex = 0;
446  int_t edge = 0;
447  int_t level = 1;
448  int_t period = 0;
449  while (1)
450  {
451  // go back with the recursion
452  // if all the edges have been processed
453  if (edge >= g. countEdges (vertex))
454  {
455  // if this was the top recursion level then quit
456  if (s_vertex. empty ())
457  break;
458 
459  // return from the recursion
460  vertex = s_vertex. top ();
461  s_vertex. pop ();
462  edge = s_edge. top ();
463  s_edge. pop ();
464  -- level;
465  continue;
466  }
467 
468  // take the next vertex
469  int_t next = g. getEdge (vertex, edge ++);
470 
471  // update the GCD of the cycle lengths
472  if (visited [next])
473  {
474  // determine the period provided by the cycle
475  int_t newPeriod = visited [next] - level - 1;
476  if (newPeriod < 0)
477  newPeriod = -newPeriod;
478 
479  // compute the GCD of the old and new periods
480  int_t a = newPeriod;
481  int_t b = period;
482  while (b)
483  {
484  int_t t = b;
485  b = a % b;
486  a = t;
487  }
488  period = a;
489 
490  // if the graph turns out to be aperiodic
491  // then immediately return the result
492  if (period == 1)
493  return period;
494  }
495 
496  // go to the deeper recursion level if not visited
497  else
498  {
499  // store the previous variables at the stacks
500  s_vertex. push (vertex);
501  s_edge. push (edge);
502 
503  // take the new vertex
504  vertex = next;
505  edge = 0;
506  ++ level;
507 
508  // mark the new vertex as visited
509  visited [vertex] = level;
510  }
511  }
512 
513  // return the computed peirod and exit
514  return period;
515 } /* computePeriod */
516 
517 /// Computes the graph that represents flow-induced relations on Morse sets.
518 /// The vertices of the result graph are the first "nVert" vertices
519 /// from the source graph. An edge is added to the new graph
520 /// iff there exists a path from one vertex to another.
521 /// Edges that come from the transitive closure are not added.
522 template <class wType>
523 inline int_t inclusionGraph (const diGraph<wType> &g, int_t nVert,
524  diGraph<wType> &result)
525 {
526  // remember the number of vertices in the input graph
527  int_t nVertG = g. countVertices ();
528  if (!nVertG)
529  return 0;
530 
531  // mark each vertex as non-visited
532  chomp::homology::BitField visited;
533  visited. allocate (nVertG);
534  visited. clearall (nVertG);
535 
536  // create the sets of vertices reachable from each vertex
537  chomp::homology::BitSets lists (nVertG, nVert);
538 
539  // prepare stacks for the DFS algorithm
540  std::stack<int_t> s_vertex;
541  std::stack<int_t> s_edge;
542 
543  // use DFS to propagate the reachable sets
544  int_t startVertex = 0;
545  int_t vertex = 0;
546  int_t edge = 0;
547  visited. set (vertex);
548  while (1)
549  {
550  // go back with the recursion
551  // if all the edges have been processed
552  if (edge >= g. countEdges (vertex))
553  {
554  // if this was the top recursion level,
555  // then find another starting point or quit
556  if (s_vertex. empty ())
557  {
558  while ((startVertex < nVertG) &&
559  (visited. test (startVertex)))
560  ++ startVertex;
561  if (startVertex >= nVertG)
562  break;
563  vertex = startVertex;
564  edge = 0;
565  visited. set (vertex);
566  continue;
567  }
568 
569  // determine the previous vertex (to trace back to)
570  int_t prev = s_vertex. top ();
571 
572  // add the list of the current vertex
573  // to the list of the previous vertex
574  // unless that was one of the first 'nVert' vertices
575  if (vertex >= nVert)
576  lists. addset (prev, vertex);
577 
578  // return from the recursion
579  vertex = prev;
580  s_vertex. pop ();
581  edge = s_edge. top ();
582  s_edge. pop ();
583  continue;
584  }
585 
586  // take the next vertex
587  int_t next = g. getEdge (vertex, edge ++);
588 
589  // add the number to the list if appropriate
590  if (next < nVert)
591  lists. add (vertex, next);
592 
593  // go to the deeper recursion level if vertex not visited
594  if (!visited. test (next))
595  {
596  // store the previous variables at the stacks
597  s_vertex. push (vertex);
598  s_edge. push (edge);
599 
600  // take the new vertex and mark as visited
601  vertex = next;
602  edge = 0;
603  visited. set (vertex);
604  }
605 
606  // add the list of the next vertex to the current one
607  // if the vertex has already been visited before
608  // and the vertex is not one of the first 'nVert' ones
609  else if (next >= nVert)
610  {
611  lists. addset (vertex, next);
612  }
613  }
614 
615  // create the result graph based on the lists of edges
616  for (int_t vertex = 0; vertex < nVert; ++ vertex)
617  {
618  result. addVertex ();
619  int_t edge = lists. get (vertex);
620  while (edge >= 0)
621  {
622  result. addEdge (edge);
623  edge = lists. get (vertex, edge + 1);
624  }
625  }
626 
627  // free memory and exit
628  visited. free ();
629  return 0;
630 } /* inclusionGraph */
631 
632 // --------------------------------------------------
633 
634 /// A more complicated version of the 'inclusionGraph' function.
635 /// Computes the graph that represents flow-induced relations on Morse sets.
636 /// The vertices of the result graph are the first "nVert" vertices
637 /// from the source graph. An edge is added to the new graph
638 /// iff there exists a path from one vertex to another.
639 /// Edges that come from the transitive closure are not added.
640 /// Records vertices along connecting orbits using operator () with the
641 /// following arguments: source, target, vertex.
642 template<class wType, class TConn>
643 inline int_t inclusionGraph (const diGraph<wType> &g, int_t nVert,
644  diGraph<wType> &result, TConn &connections)
645 {
646  // remember the number of vertices in the input graph
647  int_t nVertG = g. countVertices ();
648  if (!nVertG)
649  return 0;
650 
651  // mark each vertex as non-visited
652  chomp::homology::BitField visited;
653  visited. allocate (nVertG);
654  visited. clearall (nVertG);
655 
656  // create the sets of vertices reachable from each vertex
657  chomp::homology::BitSets forwardlists (nVertG, nVert);
658 
659  // prepare stacks for the DFS algorithm
660  std::stack<int_t> s_vertex;
661  std::stack<int_t> s_edge;
662 
663  // use DFS to propagate the reachable sets
664  int_t startVertex = 0;
665  int_t vertex = 0;
666  int_t edge = 0;
667  visited. set (vertex);
668  while (1)
669  {
670  // go back with the recursion
671  // if all the edges have been processed
672  if (edge >= g. countEdges (vertex))
673  {
674  // if this was the top recursion level,
675  // then find another starting point or quit
676  if (s_vertex. empty ())
677  {
678  while ((startVertex < nVertG) &&
679  (visited. test (startVertex)))
680  ++ startVertex;
681  if (startVertex >= nVertG)
682  break;
683  vertex = startVertex;
684  edge = 0;
685  visited. set (vertex);
686  continue;
687  }
688 
689  // determine the previous vertex (to trace back to)
690  int_t prev = s_vertex. top ();
691 
692  // add the list of the current vertex
693  // to the list of the previous vertex
694  // unless that was one of the first 'nVert' vertices
695  if (vertex >= nVert)
696  forwardlists. addset (prev, vertex);
697 
698  // return from the recursion
699  vertex = prev;
700  s_vertex. pop ();
701  edge = s_edge. top ();
702  s_edge. pop ();
703  continue;
704  }
705 
706  // take the next vertex
707  int_t next = g. getEdge (vertex, edge ++);
708 
709  // add the number to the list if appropriate
710  if (next < nVert)
711  forwardlists. add (vertex, next);
712 
713  // go to the deeper recursion level if vertex not visited
714  if (!visited. test (next))
715  {
716  // store the previous variables at the stacks
717  s_vertex. push (vertex);
718  s_edge. push (edge);
719 
720  // take the new vertex and mark as visited
721  vertex = next;
722  edge = 0;
723  visited. set (vertex);
724  }
725 
726  // add the list of the next vertex to the current one
727  // if the vertex has already been visited before
728  // and the vertex is not one of the first 'nVert' ones
729  else if (next >= nVert)
730  {
731  forwardlists. addset (vertex, next);
732  }
733  }
734 
735  // ------------------------------------------
736 
737  // create the sets of vertices that can reach each vertex
738  chomp::homology::BitSets backlists (nVertG, nVert);
739 
740  diGraph<wType> gT;
741  g. transpose (gT);
742 
743  // do a simple test for debugging purposes
744  if (!s_vertex. empty () || !s_edge. empty ())
745  throw "DEBUG: Nonempty stacks in 'inclusionGraph'.";
746 
747  // mark all the vertices as unvisited for the backward run
748  visited. clearall (nVertG);
749 
750  // use DFS to propagate the reachable sets
751  startVertex = 0;
752  vertex = 0;
753  edge = 0;
754  visited. set (vertex);
755  while (1)
756  {
757  // go back with the recursion
758  // if all the edges have been processed
759  if (edge >= gT. countEdges (vertex))
760  {
761  // if this was the top recursion level,
762  // then find another starting point or quit
763  if (s_vertex. empty ())
764  {
765  while ((startVertex < nVertG) &&
766  (visited. test (startVertex)))
767  ++ startVertex;
768  if (startVertex >= nVertG)
769  break;
770  vertex = startVertex;
771  edge = 0;
772  visited. set (vertex);
773  continue;
774  }
775 
776  // determine the previous vertex (to trace back to)
777  int_t prev = s_vertex. top ();
778 
779  // add the list of the current vertex
780  // to the list of the previous vertex
781  // unless that was one of the first 'nVert' vertices
782  if (vertex >= nVert)
783  backlists. addset (prev, vertex);
784 
785  // return from the recursion
786  vertex = prev;
787  s_vertex. pop ();
788  edge = s_edge. top ();
789  s_edge. pop ();
790  continue;
791  }
792 
793  // take the next vertex
794  int_t next = gT. getEdge (vertex, edge ++);
795 
796  // add the number to the list if appropriate
797  if (next < nVert)
798  backlists. add (vertex, next);
799 
800  // go to the deeper recursion level if vertex not visited
801  if (!visited. test (next))
802  {
803  // store the previous variables at the stacks
804  s_vertex. push (vertex);
805  s_edge. push (edge);
806 
807  // take the new vertex and mark as visited
808  vertex = next;
809  edge = 0;
810  visited. set (vertex);
811  }
812 
813  // add the list of the next vertex to the current one
814  // if the vertex has already been visited before
815  // and the vertex is not one of the first 'nVert' ones
816  else if (next >= nVert)
817  {
818  backlists. addset (vertex, next);
819  }
820  }
821 
822  // ------------------------------------------
823 
824  // record the connections to the obtained data structure
825  for (int_t v = nVert; v < nVertG; ++ v)
826  {
827  int_t bvertex = backlists. get (v);
828  while (bvertex >= 0)
829  {
830  int_t fvertex = forwardlists. get (v);
831  while (fvertex >= 0)
832  {
833  connections (bvertex, fvertex, v);
834  fvertex = forwardlists. get (v, fvertex + 1);
835  }
836  bvertex = backlists. get (v, bvertex + 1);
837  }
838  }
839 
840  // ------------------------------------------
841 
842  // create the result graph based on the lists of edges
843  for (int_t vertex = 0; vertex < nVert; ++ vertex)
844  {
845  result. addVertex ();
846  int_t edge = forwardlists. get (vertex);
847  while (edge >= 0)
848  {
849  result. addEdge (edge);
850  edge = forwardlists. get (vertex, edge + 1);
851  }
852  }
853 
854  // free memory and exit
855  visited. free ();
856  return 0;
857 } /* inclusionGraph */
858 
859 // --------------------------------------------------
860 
861 /// Creates the adjacency matrix of the given graph.
862 /// m [i] [j] is set to 1 if the graph g contains the edge i -> j,
863 /// using the assignment operator.
864 template <class wType, class matrix>
865 inline void graph2matrix (const diGraph<wType> &g, matrix &m)
866 {
867  int_t nVert = g. countVertices ();
868  for (int_t v = 0; v < nVert; ++ v)
869  {
870  int_t nEdges = g. countEdges (v);
871  for (int_t e = 0; e < nEdges; ++ e)
872  {
873  int_t w = g. getEdge (v, e);
874  m [v] [w] = 1;
875  }
876  }
877  return;
878 } /* graph2matrix */
879 
880 /// Creates a graph based on the given adjacency matrix.
881 /// If m [i] [j] is nonzero, then the edge i -> j is added to the graph.
882 /// It is assumed that the graph g is initially empty.
883 /// The size of the matrix (the number of vertices), n, must be given.
884 template <class wType, class matrix>
885 inline void matrix2graph (const matrix &m, int_t n, diGraph<wType> &g)
886 {
887  for (int_t v = 0; v < n; ++ v)
888  {
889  g. addVertex ();
890  for (int_t w = 0; w < n; ++ w)
891  if (m [v] [w])
892  g. addEdge (w);
893  }
894  return;
895 } /* matrix2graph */
896 
897 // --------------------------------------------------
898 
899 /// Computes the transitive closure of an acyclic graph
900 /// defined by its adjacency matrix, using the Warshall's algorithm:
901 /// S. Warshall, A theorem on Boolean matrices, J. ACM 9 (1962) 11-12.
902 template <class matrix>
903 inline void transitiveClosure (matrix &m, int_t n)
904 {
905  for (int_t k = 0; k < n; ++ k)
906  {
907  for (int_t i = 0; i < n; ++ i)
908  {
909  for (int_t j = 0; j < n; ++ j)
910  {
911  if (m [i] [k] && m [k] [j])
912  m [i] [j] = 1;
913  }
914  }
915  }
916  return;
917 } /* transitiveClosure */
918 
919 /// Computes the transitive reduction of a CLOSED acyclic graph
920 /// defined by its adjacency matrix, using the algorithm by D. Gries,
921 /// A.J. Martin, J.L.A. van de Snepscheut and J.T. Udding, An algorithm
922 /// for transitive reduction of an acyclic graph, Science of Computer
923 /// Programming 12 (1989), 151-155.
924 /// WARNING: The input graph MUST BE CLOSED, use the "transitiveClosure"
925 /// algorithm first if this is not the case.
926 template <class matrix>
927 inline void transitiveReduction (matrix &m, int_t n)
928 {
929  for (int_t k = n - 1; k >= 0; -- k)
930  {
931  for (int_t i = 0; i < n; ++ i)
932  {
933  for (int_t j = 0; j < n; ++ j)
934  {
935  if (m [i] [k] && m [k] [j])
936  m [i] [j] = 0;
937  }
938  }
939  }
940  return;
941 } /* transitiveReduction */
942 
943 /// Computes the transitive reduction of an arbitrary acyclic graph.
944 /// The output graph must be initially empty.
945 template <class wType>
946 inline void transitiveReduction (const diGraph<wType> &g,
947  diGraph<wType> &gRed)
948 {
949  int_t nVert = g. countVertices ();
950  if (nVert <= 0)
951  return;
952  chomp::homology::flatMatrix<char> m (nVert);
953  m. clear (0);
954  graph2matrix (g, m);
955  cymealg::transitiveClosure (m, nVert);
956  cymealg::transitiveReduction (m, nVert);
957  matrix2graph (m, nVert, gRed);
958  return;
959 } /* transitiveReduction */
960 
961 
962 } // namespace cymealg
963 
964 #endif // _CYMEALG_VARIOUS_H_
965 
966 /// @}
967 
int_t collapseVertices(const diGraph< wType > &g, int_t compNum, Table1 &compVertices, Table2 &compEnds, diGraph< wType > &result, int_t *newNumbers=0)
Collapses disjoint subsets of vertices to single vertices and creates a corresponding graph in which ...
Definition: various.h:76
void graph2matrix(const diGraph< wType > &g, matrix &m)
Creates the adjacency matrix of the given graph.
Definition: various.h:865
void transpose(const Graph1 &g, Graph2 &result, bool copyweights=false)
Creates a transposed graph.
Definition: transpose.h:43
int_t connectionGraph(const diGraph< wType > &g, int_t nVert, diGraph< wType > &result)
Computes the graph that represents connections between a number of the first vertices in the given gr...
Definition: various.h:316
This class defines a weighted directed graph with very limited number of operations.
Definition: digraph.h:75
int_t addRenumEdges2(const diGraph< wType > &g, int_t vertex, const int_t *newNum, TabSets &numSets, int_t cur, int_t *srcVert, diGraph< wType > &result)
A helper function for "collapseVertices2".
Definition: various.h:159
int_t computePeriod(const GraphType &g)
Computes the period of a strongly connected graph.
Definition: various.h:426
int_t inclusionGraph(const diGraph< wType > &g, int_t nVert, diGraph< wType > &result)
Computes the graph that represents flow-induced relations on Morse sets.
Definition: various.h:523
int_t collapseVertices2(const diGraph< wType > &g, int_t compNum, Table1 &compVertices, Table2 &compEnds, diGraph< wType > &result)
Collapses (possibly non-disjoint) subsets of vertices to single vertices and creates a corresponding ...
Definition: various.h:202
int_t addRenumEdges(const diGraph< wType > &g, int_t vertex, const int_t *newNum, int_t cur, int_t *srcVert, diGraph< wType > &result)
A helper function for "collapseVertices".
Definition: various.h:44
void transitiveReduction(matrix &m, int_t n)
Computes the transitive reduction of a CLOSED acyclic graph defined by its adjacency matrix...
Definition: various.h:927
void transitiveClosure(matrix &m, int_t n)
Computes the transitive closure of an acyclic graph defined by its adjacency matrix, using the Warshall&#39;s algorithm: S.
Definition: various.h:903
void matrix2graph(const matrix &m, int_t n, diGraph< wType > &g)
Creates a graph based on the given adjacency matrix.
Definition: various.h:885