# Combinatorica: Girth[] and FindCycle[] disagreement

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:

1. Is this a bug?
2. 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 *)

-
You've seen this, no? – J. M. Aug 24 '12 at 17:35
@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. Aug 24 '12 at 17:48
@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

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@
SystemGraph[{1 \[DirectedEdge] 2, 2 \[DirectedEdge] 3, 3 \[DirectedEdge] 1}]


3

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

Girth@ToCombinatoricaGraph@
SystemGraph[{1 \[DirectedEdge] 2, 2 \[DirectedEdge] 3, 1 \[DirectedEdge] 3}]


3

-
:) 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
@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 5: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 6:37

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];
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},
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.

-
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 0: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