Combinatorica: Girth[] and FindCycle[] disagreement

Question

Warning: run the following code in a fresh Mma session, as some symbols could be shadowed (depending on your Mma version)

While trying to answer this question, I fell into the following:

 (* Let's load a large Directed Graph and convert it to Combinatorica *)
g = Graph@Union@Flatten[
     Thread[DirectedEdge @@ ##] & /@ Select[{#, IsotopeData[#, "DaughterNuclides"]} & /@ 
        IsotopeData[], #[[2]] != {} &]];
Needs["GraphUtilities`"]
<< Combinatorica`
cg = ToCombinatoricaGraph[g];

Girth[cg] gives the length of a shortest cycle in a simple graph g.

So, let's check if cg is Simple and calculate its Girth:

{SimpleQ@cg, Girth@cg}
(*
-> {True, 3}
*)

So there is at least one Cycle in cg of length 3.

But look what happens when we try to find it by the two available methods in Combinatorica:

{ExtractCycles@cg, FindCycle@cg}
(*
-> {{},{}}
*)

So, two questions:

Is this a bug?
What is the easiest way to find all cycles in g without using Combinatorica?

Edit

BTW, the (now) standard Graph functionality also detects cycles:

AcyclicGraphQ[g]
(* -> False *)

@J.M. Nope, thanks. If the answers there work for a directed graph I'll delete this q — Dr. belisarius, Aug 24 '12 at 17:37
@J.M. I don't think those solutions scale well for this graph (VertexCount > 3000) — Dr. belisarius, Aug 24 '12 at 17:46
Apparently so, since I just tried them out. An efficient (in space, time, or preferably both) cycle finder would be a very nice thing. — J. M.'s missing motivation, Aug 24 '12 at 17:48
I don't have access to Mathematica 8, but I suspect that the problem may be connnected with ToCombinatoricaGraph (see my old comment here). What happens if you use FromOrderedPairs@EdgeList@g to convert to a Combinatorica Graph? — user1066, Aug 24 '12 at 18:29
@TomD I used g1 = FromOrderedPairs[ EdgeList@g /. Thread[Rule[VertexList[g], Range@VertexCount@g]] /. DirectedEdge -> List, Type -> Directed]; Print@{Girth@g1,FindCycle@g1}; and the result is the same. Thanks for the pointer! — Dr. belisarius, Aug 24 '12 at 18:46

score 7 · Answer 1 · edited Apr 13 '17 at 12:55

7

Edited for correctness:

I use the variant for directed graphs from here. I take your graph as above, extract edges, rename so vertices are integers from 1 to #vertices. After finishing we revert to the original names.

ee = EdgeList[g];
vv = VertexList[g];
reprule = Thread[vv -> Range[Length[vv]]];
revrule = Map[Reverse, reprule];
pairs = ee /. reprule /. DirectedEdge -> List;

extendCycle[cyc_List, edges_List] := 
 Map[If[# > First[cyc] && ! MemberQ[cyc, #], Append[cyc, #], 
    Null ] &, edges[[Last[cyc]]]] /. Null :> Sequence[]

cycles[omat_, k_] := Module[
  {n = Length[Union[Flatten@omat]], m2, cyc, cyclist, mat},
  mat = Join[omat, Thread[{Range[n], 0}]];
  m2 = Map[Last, SplitBy[Sort[mat], First], {2}];
  m2 = m2 /. 0 :> Sequence[];
  cyclist = 
   Flatten[Drop[MapIndexed[{#2[[1]], #1} &, m2, {2}], -k + 1], 1];
  cyclist = Select[cyclist, #[[2]] > #[[1]] &];
  Do[cyclist = 
    Flatten[Map[extendCycle[#, m2] &, cyclist], 1], {k - 2}];
  Map[If[MemberQ[m2[[Last[#]]], First[#]], Append[#, First[#]], 
      Null] &, cyclist] /. Null :> Sequence[]]

I get no cycles in any length between 3 and 20. Have not tried further. I now believe the code will behave correctly. Let me know if not.

edited Apr 13 '17 at 12:55

Community

1

answered Aug 24 '12 at 19:17

Daniel Lichtblau

58,970
2
101
199

Those cyc3 don't look good (from a physical viewpoint). For example, in the first one check IsotopeData["Iron62", "DaughterNuclides"] ... its only daughter is {"Cobalt62"} – Dr. belisarius Aug 24 '12 at 19:53
@Verde Thanks. Probably programmer error. Though I suspect the underlying method is sound. – Daniel Lichtblau Aug 24 '12 at 20:00
I corrected the code. It was not taking into account vertices with no outgoing edges. Should behave better now. – Daniel Lichtblau Aug 24 '12 at 23:22
cycles[pairs, 3] gives me {{1025, 2035, 1823, 1025}, {1027, 2036, 1825, 1027}, {1027, 2037, 1827, 1027}}. Am I missing something? – Dr. belisarius Aug 25 '12 at 00:13
@Verde Sorry, I left out an important part of the preprocessing. It now needs: pairs = ee /. reprule /. DirectedEdge -> List; Should be better with that. – Daniel Lichtblau Aug 25 '12 at 15:59
Now works like a charm! Thanks! – Dr. belisarius Aug 25 '12 at 18:00

score 4 · Accepted Answer · answered Aug 24 '12 at 20:40

4

I don't think it's a bug. Directed radioactive decay graphs shouldn't have cycles by definition and AcyclicGraphQ doesn't see them either:

AcyclicGraphQ@g

True

You don't seem to have heeded your own warning to start with a fresh kernel, as playing around with the code for a while gave me False too.

The output of ExtractCycles and FindCycle is therefore correct.

Girth doesn't seem to take directionality into account when determining cycles:

Cyclic directed graph:

Girth@ ToCombinatoricaGraph@
     System`Graph[{1 \[DirectedEdge] 2, 2 \[DirectedEdge] 3, 3 \[DirectedEdge] 1}]

3

Acyclic directed graph (note the reversed direction of the last edge):

Girth@ToCombinatoricaGraph@
     System`Graph[{1 \[DirectedEdge] 2, 2 \[DirectedEdge] 3, 1 \[DirectedEdge] 3}]

3

answered Aug 24 '12 at 20:40

Sjoerd C. de Vries

65,815
14
188
323

:) Of course I expected no cycles, that was what bothered me. I think WRI has a great opportunity to improve in the graph's features robustness – Dr. belisarius Aug 24 '12 at 20:45
1

@verde The combination of the built-in Graph and Combinatorica is just poison. – Sjoerd C. de Vries Aug 24 '12 at 20:47
I only hope v9 will bring all Combinatorica functionality into System. – Dr. belisarius Aug 25 '12 at 05:46
@verde Same here. BTW I just loph the graph. It's strange, monstrous, caterpillar-like; almost like it conveys a message. – Sjoerd C. de Vries Aug 25 '12 at 06:37

Combinatorica: Girth[] and FindCycle[] disagreement

2 Answers2

Linked